Tempus  Version of the Day
Time Integration
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Tempus_Stepper_impl.hpp
Go to the documentation of this file.
1 // @HEADER
2 // ****************************************************************************
3 // Tempus: Copyright (2017) Sandia Corporation
4 //
5 // Distributed under BSD 3-clause license (See accompanying file Copyright.txt)
6 // ****************************************************************************
7 // @HEADER
8 
9 #ifndef Tempus_Stepper_impl_hpp
10 #define Tempus_Stepper_impl_hpp
11 
12 #include "NOX_Thyra.H"
13 
14 
15 namespace Tempus {
16 
17 
18 template<class Scalar>
20 {
21  Teuchos::RCP<Teuchos::FancyOStream> out = this->getOStream();
22  out->setOutputToRootOnly(0);
23 
24  bool isValidSetup = this->isValidSetup(*out);
25 
26  if (isValidSetup)
27  this->isInitialized_ = true; // Only place it is set to true.
28  else
29  this->describe(*out, Teuchos::VERB_MEDIUM);
30 }
31 
32 
33 template<class Scalar>
35 {
36  if ( !this->isInitialized() ) {
37  this->describe( *(this->getOStream()), Teuchos::VERB_MEDIUM);
38  TEUCHOS_TEST_FOR_EXCEPTION( !this->isInitialized(), std::logic_error,
39  "Error - " << this->description() << " is not initialized!");
40  }
41 }
42 
43 
44 template<class Scalar>
46 {
47  if (a == false) {
48  Teuchos::RCP<Teuchos::FancyOStream> out = this->getOStream();
49  out->setOutputToRootOnly(0);
50  Teuchos::OSTab ostab(out,1,"Stepper::setUseFSALTrueOnly()");
51  *out << "Warning -- useFSAL for '" << this->getStepperType() << "'\n"
52  << "can only be set to true. Leaving set to true." << std::endl;
53  }
54  useFSAL_ = true;
55 }
56 
57 
58 template<class Scalar>
60 {
61  if (a == true) {
62  Teuchos::RCP<Teuchos::FancyOStream> out = this->getOStream();
63  out->setOutputToRootOnly(0);
64  Teuchos::OSTab ostab(out,1,"Stepper::setUseFSALFalseOnly()");
65  *out << "Warning -- useFSAL for '" << this->getStepperType() << "'\n"
66  << "can only be set to false. Leaving set to false." << std::endl;
67  }
68  useFSAL_ = false;
69 }
70 
71 
72 template<class Scalar>
75 {
76  TEUCHOS_TEST_FOR_EXCEPTION( stepperX_ == Teuchos::null, std::logic_error,
77  "Error - stepperX_ has not been set in setInitialConditions() or\n"
78  " can not be set from the state!\n");
79 
80  return stepperX_;
81 }
82 
83 template<class Scalar>
86 {
87  TEUCHOS_TEST_FOR_EXCEPTION( stepperXDot_ == Teuchos::null, std::logic_error,
88  "Error - stepperXDot_ has not been set in setInitialConditions() or\n"
89  " can not be set from the state!\n");
90 
91  return stepperXDot_;
92 }
93 
94 template<class Scalar>
97 {
98  TEUCHOS_TEST_FOR_EXCEPTION( stepperXDotDot_==Teuchos::null, std::logic_error,
99  "Error - stepperXDotDot_ has not been set in setInitialConditions() or\n"
100  " can not be set from the state!\n");
101 
102  return stepperXDotDot_;
103 }
104 
105 
106 template<class Scalar>
109 {
110  if (state->getXDotDot() != Teuchos::null) stepperXDotDot_=state->getXDotDot();
111  // Else use temporary storage stepperXDotDot_ which should have been set in
112  // setInitialConditions().
113 
114  TEUCHOS_TEST_FOR_EXCEPTION( stepperXDotDot_==Teuchos::null, std::logic_error,
115  "Error - stepperXDotDot_ has not been set in setInitialConditions() or\n"
116  " can not be set from the state!\n");
117 
118  return stepperXDotDot_;
119 }
120 
121 
122 template<class Scalar>
124  const Teuchos::EVerbosityLevel verbLevel) const
125 {
126  auto l_out = Teuchos::fancyOStream( out.getOStream() );
127  Teuchos::OSTab ostab(*l_out, 2, this->description());
128  l_out->setOutputToRootOnly(0);
129 
130  *l_out << "--- Stepper ---\n"
131  << " isInitialized_ = " << Teuchos::toString(isInitialized_) << std::endl
132  << " stepperType_ = " << stepperType_ << std::endl
133  << " useFSAL_ = " << Teuchos::toString(useFSAL_) << std::endl
134  << " ICConsistency_ = " << ICConsistency_ << std::endl
135  << " ICConsistencyCheck_ = " << Teuchos::toString(ICConsistencyCheck_) << std::endl
136  << " stepperX_ = " << stepperX_ << std::endl
137  << " stepperXDot_ = " << stepperXDot_ << std::endl
138  << " stepperXDotDot_ = " << stepperXDotDot_ << std::endl;
139 }
140 
141 
142 template<class Scalar>
144  Teuchos::FancyOStream & out) const
145 {
146  out.setOutputToRootOnly(0);
147  bool isValidSetup = true;
148 
149  if ( !(ICConsistency_ == "None" || ICConsistency_ == "Zero" ||
150  ICConsistency_ == "App" || ICConsistency_ == "Consistent") ) {
151  isValidSetup = false;
152  auto l_out = Teuchos::fancyOStream( out.getOStream() );
153  l_out->setOutputToRootOnly(0);
154  *l_out << "The IC consistency does not have a valid value!\n"
155  << "('None', 'Zero', 'App' or 'Consistent')\n"
156  << " ICConsistency = " << ICConsistency_ << "\n";
157  }
158 
159  return isValidSetup;
160 }
161 
162 
163 template<class Scalar>
166 {
167  if (pl != Teuchos::null) {
168  this->setStepperName(pl->name());
169  auto stepperType =
170  pl->get<std::string>("Stepper Type", this->getStepperType());
172  stepperType != this->getStepperType() ,std::runtime_error,
173  " ParameterList 'Stepper Type' (='" + stepperType +"')\n"
174  " does not match type for this Stepper (='"
175  + this->getStepperType() + "').");
176  this->setStepperType(stepperType);
177 
178  this->setUseFSAL(pl->get<bool>("Use FSAL", this->getUseFSAL()));
179  this->setICConsistency(
180  pl->get<std::string>("Initial Condition Consistency",
181  this->getICConsistency()));
182  this->setICConsistencyCheck(
183  pl->get<bool>("Initial Condition Consistency Check",
184  this->getICConsistencyCheck()));
185  }
186 }
187 
188 
189 template<class Scalar>
192 {
193  return this->getValidParametersBasic();
194 }
195 
196 
197 template<class Scalar>
200 {
201  auto pl = Teuchos::parameterList(this->getStepperName());
202  pl->template set<std::string>("Stepper Type", this->getStepperType());
203 
204  pl->template set<bool>("Use FSAL", this->getUseFSAL(),
205  "The First-Same-As-Last (FSAL) principle is the situation where the\n"
206  "last function evaluation, f(x^{n-1},t^{n-1}) [a.k.a. xDot^{n-1}],\n"
207  "can be used for the first function evaluation, f(x^n,t^n)\n"
208  "[a.k.a. xDot^n]. For RK methods, this applies to the stages.\n"
209  "\n"
210  "Often the FSAL priniciple can be used to save an evaluation.\n"
211  "However there are cases when it cannot be used, e.g., operator\n"
212  "splitting where other steppers/operators have modified the solution,\n"
213  "x^*, and thus require the function evaluation, f(x^*, t^{n-1}).\n"
214  "\n"
215  "It should be noted that when the FSAL priniciple can be used\n"
216  "(can set useFSAL=true), setting useFSAL=false will give the\n"
217  "same solution but at additional expense. However, the reverse\n"
218  "is not true. When the FSAL priniciple can not be used\n"
219  "(need to set useFSAL=false), setting useFSAL=true will produce\n"
220  "incorrect solutions.\n"
221  "\n"
222  "Default in general for explicit and implicit steppers is false,\n"
223  "but individual steppers can override this default.");
224 
225  pl->template set<std::string>("Initial Condition Consistency",this->getICConsistency(),
226  "This indicates which type of consistency should be applied to\n"
227  "the initial conditions (ICs):\n"
228  "\n"
229  " 'None' - Do nothing to the ICs provided in the SolutionHistory.\n"
230  " 'Zero' - Set the derivative of the SolutionState to zero in the\n"
231  " SolutionHistory provided, e.g., xDot^0 = 0, or \n"
232  " xDotDot^0 = 0.\n"
233  " 'App' - Use the application's ICs, e.g., getNominalValues().\n"
234  " 'Consistent' - Make the initial conditions for x and xDot\n"
235  " consistent with the governing equations, e.g.,\n"
236  " xDot = f(x,t), and f(x, xDot, t) = 0. For implicit\n"
237  " ODEs, this requires a solve of f(x, xDot, t) = 0 for\n"
238  " xDot, and another Jacobian and residual may be\n"
239  " needed, e.g., boundary conditions on xDot may need\n"
240  " to replace boundary conditions on x.\n"
241  "\n"
242  "In general for explicit steppers, the default is 'Consistent',\n"
243  "because it is fairly cheap with just one residual evaluation.\n"
244  "In general for implicit steppers, the default is 'None', because\n"
245  "the application often knows its IC and can set it the initial\n"
246  "SolutionState. Also, as noted above, 'Consistent' may require\n"
247  "another Jacobian from the application. Individual steppers may\n"
248  "override these defaults.");
249 
250  pl->template set<bool>("Initial Condition Consistency Check", this->getICConsistencyCheck(),
251  "Check if the initial condition, x and xDot, is consistent with the\n"
252  "governing equations, xDot = f(x,t), or f(x, xDot, t) = 0.\n"
253  "\n"
254  "In general for explicit and implicit steppers, the default is true,\n"
255  "because it is fairly cheap with just one residual evaluation.\n"
256  "Individual steppers may override this default.");
257 
258  return pl;
259 }
260 
261 
262 // Nonmember Helper Functions.
263 // ------------------------------------------------------------------------
264 
265 template<class Scalar>
267  const Teuchos::RCP<const Thyra::ModelEvaluator<Scalar> >& model)
268 {
270  typedef Thyra::ModelEvaluatorBase MEB;
271  const MEB::InArgs<Scalar> inArgs = model->createInArgs();
272  const MEB::OutArgs<Scalar> outArgs = model->createOutArgs();
273  const bool supports = inArgs.supports(MEB::IN_ARG_x) &&
274  outArgs.supports(MEB::OUT_ARG_f);
275 
276  TEUCHOS_TEST_FOR_EXCEPTION( supports == false, std::logic_error,
277  model->description() << " can not support an explicit ODE with\n"
278  << " IN_ARG_x = " << inArgs.supports(MEB::IN_ARG_x) << "\n"
279  << " OUT_ARG_f = " << outArgs.supports(MEB::OUT_ARG_f) << "\n"
280  << "Explicit ODE requires:\n"
281  << " IN_ARG_x = true\n"
282  << " OUT_ARG_f = true\n"
283  << "\n"
284  << "NOTE: Currently the convention to evaluate f(x,t) is to set\n"
285  << "xdot=null! There is no InArgs support to test if xdot is null,\n"
286  << "so we set xdot=null and hope the ModelEvaluator can handle it.\n");
287 
288  return;
289 }
290 
291 
292 template<class Scalar>
294  const Teuchos::RCP<const Thyra::ModelEvaluator<Scalar> >& model)
295 {
297  typedef Thyra::ModelEvaluatorBase MEB;
298  const MEB::InArgs<Scalar> inArgs = model->createInArgs();
299  const MEB::OutArgs<Scalar> outArgs = model->createOutArgs();
300  const bool supports = inArgs.supports(MEB::IN_ARG_x) &&
301  inArgs.supports(MEB::IN_ARG_x_dot) &&
302  outArgs.supports(MEB::OUT_ARG_f);
303 
304  TEUCHOS_TEST_FOR_EXCEPTION( supports == false, std::logic_error,
305  model->description() << "can not support an explicit ODE with\n"
306  << " IN_ARG_x = " << inArgs.supports(MEB::IN_ARG_x) << "\n"
307  << " IN_ARG_x_dot = " << inArgs.supports(MEB::IN_ARG_x_dot) << "\n"
308  << " OUT_ARG_f = " << outArgs.supports(MEB::OUT_ARG_f) << "\n"
309  << "Explicit ODE requires:\n"
310  << " IN_ARG_x = true\n"
311  << " IN_ARG_x_dot = true\n"
312  << " OUT_ARG_f = true\n"
313  << "\n"
314  << "NOTE: Currently the convention to evaluate f(x, xdot, t) is to\n"
315  << "set xdotdot=null! There is no InArgs support to test if xdotdot\n"
316  << "is null, so we set xdotdot=null and hope the ModelEvaluator can\n"
317  << "handle it.\n");
318 
319  return;
320 }
321 
322 
323 template<class Scalar>
325  const Teuchos::RCP<const Thyra::ModelEvaluator<Scalar> >& model)
326 {
328  typedef Thyra::ModelEvaluatorBase MEB;
329  const MEB::InArgs<Scalar> inArgs = model->createInArgs();
330  const MEB::OutArgs<Scalar> outArgs = model->createOutArgs();
331  const bool supports = inArgs.supports(MEB::IN_ARG_x) &&
332  inArgs.supports(MEB::IN_ARG_x_dot) &&
333  inArgs.supports(MEB::IN_ARG_alpha) &&
334  inArgs.supports(MEB::IN_ARG_beta) &&
335  !inArgs.supports(MEB::IN_ARG_W_x_dot_dot_coeff) &&
336  outArgs.supports(MEB::OUT_ARG_f) &&
337  outArgs.supports(MEB::OUT_ARG_W);
338 
339  TEUCHOS_TEST_FOR_EXCEPTION( supports == false, std::logic_error,
340  model->description() << " can not support an implicit ODE with\n"
341  << " IN_ARG_x = "
342  << inArgs.supports(MEB::IN_ARG_x) << "\n"
343  << " IN_ARG_x_dot = "
344  << inArgs.supports(MEB::IN_ARG_x_dot) << "\n"
345  << " IN_ARG_alpha = "
346  << inArgs.supports(MEB::IN_ARG_alpha) << "\n"
347  << " IN_ARG_beta = "
348  << inArgs.supports(MEB::IN_ARG_beta) << "\n"
349  << " IN_ARG_W_x_dot_dot_coeff = "
350  << inArgs.supports(MEB::IN_ARG_W_x_dot_dot_coeff) << "\n"
351  << " OUT_ARG_f = "
352  << outArgs.supports(MEB::OUT_ARG_f) << "\n"
353  << " OUT_ARG_W = "
354  << outArgs.supports(MEB::OUT_ARG_W) << "\n"
355  << "Implicit ODE requires:\n"
356  << " IN_ARG_x = true\n"
357  << " IN_ARG_x_dot = true\n"
358  << " IN_ARG_alpha = true\n"
359  << " IN_ARG_beta = true\n"
360  << " IN_ARG_W_x_dot_dot_coeff = false\n"
361  << " OUT_ARG_f = true\n"
362  << " OUT_ARG_W = true\n");
363 
364  return;
365 }
366 
367 
368 template<class Scalar>
370  const Teuchos::RCP<const Thyra::ModelEvaluator<Scalar> >& model)
371 {
373  typedef Thyra::ModelEvaluatorBase MEB;
374  const MEB::InArgs<Scalar> inArgs = model->createInArgs();
375  const MEB::OutArgs<Scalar> outArgs = model->createOutArgs();
376  const bool supports = inArgs.supports(MEB::IN_ARG_x) &&
377  inArgs.supports(MEB::IN_ARG_x_dot) &&
378  inArgs.supports(MEB::IN_ARG_x_dot_dot) &&
379  inArgs.supports(MEB::IN_ARG_alpha) &&
380  inArgs.supports(MEB::IN_ARG_beta) &&
381  inArgs.supports(MEB::IN_ARG_W_x_dot_dot_coeff) &&
382  outArgs.supports(MEB::OUT_ARG_f) &&
383  outArgs.supports(MEB::OUT_ARG_W);
384 
385  TEUCHOS_TEST_FOR_EXCEPTION( supports == false, std::logic_error,
386  model->description() << " can not support an implicit ODE with\n"
387  << " IN_ARG_x = "
388  << inArgs.supports(MEB::IN_ARG_x) << "\n"
389  << " IN_ARG_x_dot = "
390  << inArgs.supports(MEB::IN_ARG_x_dot) << "\n"
391  << " IN_ARG_x_dot_dot = "
392  << inArgs.supports(MEB::IN_ARG_x_dot_dot) << "\n"
393  << " IN_ARG_alpha = "
394  << inArgs.supports(MEB::IN_ARG_alpha) << "\n"
395  << " IN_ARG_beta = "
396  << inArgs.supports(MEB::IN_ARG_beta) << "\n"
397  << " IN_ARG_W_x_dot_dot_coeff = "
398  << inArgs.supports(MEB::IN_ARG_W_x_dot_dot_coeff) << "\n"
399  << " OUT_ARG_f = "
400  << outArgs.supports(MEB::OUT_ARG_f) << "\n"
401  << " OUT_ARG_W = "
402  << outArgs.supports(MEB::OUT_ARG_W) << "\n"
403  << "Implicit Second Order ODE requires:\n"
404  << " IN_ARG_x = true\n"
405  << " IN_ARG_x_dot = true\n"
406  << " IN_ARG_x_dot_dot = true\n"
407  << " IN_ARG_alpha = true\n"
408  << " IN_ARG_beta = true\n"
409  << " IN_ARG_W_x_dot_dot_coeff = true\n"
410  << " OUT_ARG_f = true\n"
411  << " OUT_ARG_W = true\n");
412 
413  return;
414 }
415 
416 
418 {
419  using Teuchos::RCP;
421 
422  // NOX Solver ParameterList
423  RCP<ParameterList> noxPL = Teuchos::parameterList();
424 
425  // Direction ParameterList
426  RCP<ParameterList> directionPL = Teuchos::parameterList();
427  directionPL->set<std::string>("Method", "Newton");
428  RCP<ParameterList> newtonPL = Teuchos::parameterList();
429  newtonPL->set<std::string>("Forcing Term Method", "Constant");
430  newtonPL->set<bool> ("Rescue Bad Newton Solve", 1);
431  directionPL->set("Newton", *newtonPL);
432  noxPL->set("Direction", *directionPL);
433 
434  // Line Search ParameterList
435  RCP<ParameterList> lineSearchPL = Teuchos::parameterList();
436  lineSearchPL->set<std::string>("Method", "Full Step");
437  RCP<ParameterList> fullStepPL = Teuchos::parameterList();
438  fullStepPL->set<double>("Full Step", 1);
439  lineSearchPL->set("Full Step", *fullStepPL);
440  noxPL->set("Line Search", *lineSearchPL);
441 
442  noxPL->set<std::string>("Nonlinear Solver", "Line Search Based");
443 
444  // Printing ParameterList
445  RCP<ParameterList> printingPL = Teuchos::parameterList();
446  printingPL->set<int>("Output Precision", 3);
447  printingPL->set<int>("Output Processor", 0);
448  RCP<ParameterList> outputPL = Teuchos::parameterList();
449  outputPL->set<bool>("Error", 1);
450  outputPL->set<bool>("Warning", 1);
451  outputPL->set<bool>("Outer Iteration", 0);
452  outputPL->set<bool>("Parameters", 0);
453  outputPL->set<bool>("Details", 0);
454  outputPL->set<bool>("Linear Solver Details", 1);
455  outputPL->set<bool>("Stepper Iteration", 1);
456  outputPL->set<bool>("Stepper Details", 1);
457  outputPL->set<bool>("Stepper Parameters", 1);
458  printingPL->set("Output Information", *outputPL);
459  noxPL->set("Printing", *printingPL);
460 
461  // Solver Options ParameterList
462  RCP<ParameterList> solverOptionsPL = Teuchos::parameterList();
463  solverOptionsPL->set<std::string>("Status Test Check Type", "Minimal");
464  noxPL->set("Solver Options", *solverOptionsPL);
465 
466  // Status Tests ParameterList
467  RCP<ParameterList> statusTestsPL = Teuchos::parameterList();
468  statusTestsPL->set<std::string>("Test Type", "Combo");
469  statusTestsPL->set<std::string>("Combo Type", "OR");
470  statusTestsPL->set<int>("Number of Tests", 2);
471  RCP<ParameterList> test0PL = Teuchos::parameterList();
472  test0PL->set<std::string>("Test Type", "NormF");
473  test0PL->set<double>("Tolerance", 1e-08);
474  statusTestsPL->set("Test 0", *test0PL);
475  RCP<ParameterList> test1PL = Teuchos::parameterList();
476  test1PL->set<std::string>("Test Type", "MaxIters");
477  test1PL->set<int>("Maximum Iterations", 10);
478  statusTestsPL->set("Test 1", *test1PL);
479  noxPL->set("Status Tests", *statusTestsPL);
480 
481  // Solver ParameterList
482  RCP<ParameterList> solverPL = Teuchos::parameterList("Default Solver");
483  solverPL->set("NOX", *noxPL);
484 
485  return solverPL;
486 }
487 
488 
489 } // namespace Tempus
490 #endif // Tempus_Stepper_impl_hpp
const std::string & name() const
virtual bool isValidSetup(Teuchos::FancyOStream &out) const
bool is_null(const boost::shared_ptr< T > &p)
T & get(const std::string &name, T def_value)
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
Teuchos::RCP< Teuchos::ParameterList > defaultSolverParameters()
Returns the default solver ParameterList for implicit Steppers.
void validExplicitODE(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &model)
Validate that the model supports explicit ODE evaluation, f(x,t) [=xdot].
virtual void initialize()
Initialize after construction and changing input parameters.
Teuchos::RCP< Teuchos::ParameterList > getValidParametersBasic() const
Add basic parameters to Steppers ParameterList.
void setUseFSALFalseOnly(bool a)
virtual void checkInitialized()
Check initialization, and error out on failure.
virtual void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel) const
void setUseFSALTrueOnly(bool a)
virtual Teuchos::RCP< Thyra::VectorBase< Scalar > > getStepperX()
Get Stepper x.
basic_FancyOStream & setOutputToRootOnly(const int rootRank)
void setStepperValues(const Teuchos::RCP< Teuchos::ParameterList > pl)
Set Stepper member data from ParameterList.
void validSecondOrderExplicitODE(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &model)
Validate that the model supports explicit second order ODE evaluation, f(x,xdot,t) [=xdotdot]...
void validSecondOrderODE_DAE(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &model)
Validate ME supports 2nd order implicit ODE/DAE evaluation, f(xdotdot,xdot,x,t) [= 0]...
RCP< std::basic_ostream< char_type, traits_type > > getOStream()
virtual Teuchos::RCP< Thyra::VectorBase< Scalar > > getStepperXDotDot()
Get Stepper xDotDot.
void validImplicitODE_DAE(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &model)
Validate ME supports implicit ODE/DAE evaluation, f(xdot,x,t) [= 0].
#define TEUCHOS_TEST_FOR_EXCEPT(throw_exception_test)
Solution state for integrators and steppers.
std::string toString(const T &t)
virtual Teuchos::RCP< Thyra::VectorBase< Scalar > > getStepperXDot()
Get Stepper xDot.
virtual Teuchos::RCP< const Teuchos::ParameterList > getValidParameters() const