Tempus Version of the Day
Time Integration
Loading...
Searching...
No Matches
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
15namespace Tempus {
16
17
18template<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
33template<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
44template<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
58template<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
72template<class Scalar>
73Teuchos::RCP<Thyra::VectorBase<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
83template<class Scalar>
84Teuchos::RCP<Thyra::VectorBase<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
94template<class Scalar>
95Teuchos::RCP<Thyra::VectorBase<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
106template<class Scalar>
107Teuchos::RCP<Thyra::VectorBase<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
122template<class Scalar>
123void Stepper<Scalar>::describe(Teuchos::FancyOStream & out,
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
142template<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
163template<class Scalar>
165setStepperValues(Teuchos::RCP<Teuchos::ParameterList> pl)
166{
167 if (pl != Teuchos::null) {
168 this->setStepperName(pl->name());
169 auto stepperType =
170 pl->get<std::string>("Stepper Type", this->getStepperType());
171 TEUCHOS_TEST_FOR_EXCEPTION(
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
189template<class Scalar>
190Teuchos::RCP<const Teuchos::ParameterList>
192{
193 return this->getValidParametersBasic();
194}
195
196
197template<class Scalar>
198Teuchos::RCP<Teuchos::ParameterList>
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
265template<class Scalar>
267 const Teuchos::RCP<const Thyra::ModelEvaluator<Scalar> >& model)
268{
269 TEUCHOS_TEST_FOR_EXCEPT( is_null(model) );
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
292template<class Scalar>
294 const Teuchos::RCP<const Thyra::ModelEvaluator<Scalar> >& model)
295{
296 TEUCHOS_TEST_FOR_EXCEPT( is_null(model) );
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
323template<class Scalar>
325 const Teuchos::RCP<const Thyra::ModelEvaluator<Scalar> >& model)
326{
327 TEUCHOS_TEST_FOR_EXCEPT( is_null(model) );
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
368template<class Scalar>
370 const Teuchos::RCP<const Thyra::ModelEvaluator<Scalar> >& model)
371{
372 TEUCHOS_TEST_FOR_EXCEPT( is_null(model) );
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
417Teuchos::RCP<Teuchos::ParameterList> defaultSolverParameters()
418{
419 using Teuchos::RCP;
420 using Teuchos::ParameterList;
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
Solution state for integrators and steppers. SolutionState contains the metadata for solutions and th...
virtual Teuchos::RCP< Thyra::VectorBase< Scalar > > getStepperXDotDot()
Get Stepper xDotDot.
void setUseFSALFalseOnly(bool a)
virtual Teuchos::RCP< Thyra::VectorBase< Scalar > > getStepperXDot()
Get Stepper xDot.
void setStepperValues(const Teuchos::RCP< Teuchos::ParameterList > pl)
Set Stepper member data from ParameterList.
virtual void initialize()
Initialize after construction and changing input parameters.
Teuchos::RCP< Teuchos::ParameterList > getValidParametersBasic() const
Add basic parameters to Steppers ParameterList.
virtual bool isValidSetup(Teuchos::FancyOStream &out) const
virtual void checkInitialized()
Check initialization, and error out on failure.
virtual Teuchos::RCP< const Teuchos::ParameterList > getValidParameters() const
virtual Teuchos::RCP< Thyra::VectorBase< Scalar > > getStepperX()
Get Stepper x.
void setUseFSALTrueOnly(bool a)
virtual void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel) const
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].
void validImplicitODE_DAE(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &model)
Validate ME supports implicit ODE/DAE evaluation, f(xdot,x,t) [= 0].
void validExplicitODE(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &model)
Validate that the model supports explicit ODE evaluation, f(x,t) [=xdot].
Teuchos::RCP< Teuchos::ParameterList > defaultSolverParameters()
Returns the default solver ParameterList for implicit Steppers.