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