Tempus Version of the Day
Time Integration
Loading...
Searching...
No Matches
Tempus_UnitTest_Subcycling.cpp
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
10
11#include "Tempus_StepperForwardEuler.hpp"
12#include "Tempus_StepperBackwardEuler.hpp"
13
14#include "Tempus_StepperSubcycling.hpp"
21
22
23namespace Tempus_Unit_Test {
24
25using Teuchos::RCP;
26using Teuchos::rcp;
27using Teuchos::rcp_const_cast;
28using Teuchos::rcp_dynamic_cast;
29using Teuchos::ParameterList;
30using Teuchos::sublist;
31
33
34
35// ************************************************************
36// ************************************************************
37TEUCHOS_UNIT_TEST(Subcycling, Default_Construction)
38{
39 auto model = rcp(new Tempus_Test::SinCosModel<double>());
40 auto modelME = rcp_dynamic_cast<const Thyra::ModelEvaluator<double> > (model);
41
42 // Setup SolutionHistory ------------------------------------
43 auto inArgsIC = model->getNominalValues();
44 auto icSolution =rcp_const_cast<Thyra::VectorBase<double> >(inArgsIC.get_x());
45 auto icState = Tempus::createSolutionStateX(icSolution);
46 auto solutionHistory = rcp(new Tempus::SolutionHistory<double>());
47 solutionHistory->addState(icState);
48 solutionHistory->initWorkingState();
49
50 // Default construction.
51 auto stepper = rcp(new Tempus::StepperSubcycling<double>());
52 auto stepperBE = Tempus::createStepperBackwardEuler(modelME, Teuchos::null);
53 stepper->setSubcyclingStepper(stepperBE);
54 stepper->setInitialConditions(solutionHistory);
55 stepper->initialize();
56 TEUCHOS_TEST_FOR_EXCEPT(!stepper->isInitialized());
57
58 // Default values for construction.
59 auto modifier = rcp(new Tempus::StepperSubcyclingModifierDefault<double>());
60 auto modifierX = rcp(new Tempus::StepperSubcyclingModifierXDefault<double>());
61 auto observer = rcp(new Tempus::StepperSubcyclingObserverDefault<double>());
62 auto solver = rcp(new Thyra::NOXNonlinearSolver());
63 solver->setParameterList(Tempus::defaultSolverParameters());
64
65 bool useFSAL = stepper->getUseFSAL();
66 std::string ICConsistency = stepper->getICConsistency();
67 bool ICConsistencyCheck = stepper->getICConsistencyCheck();
68
69 // Test the set functions
70 stepper->setSolver(solver); stepper->initialize(); TEUCHOS_TEST_FOR_EXCEPT(!stepper->isInitialized());
71 stepper->setAppAction(modifier); stepper->initialize(); TEUCHOS_TEST_FOR_EXCEPT(!stepper->isInitialized());
72 stepper->setAppAction(modifierX); stepper->initialize(); TEUCHOS_TEST_FOR_EXCEPT(!stepper->isInitialized());
73 stepper->setAppAction(observer); stepper->initialize(); TEUCHOS_TEST_FOR_EXCEPT(!stepper->isInitialized());
74 stepper->setUseFSAL(useFSAL); stepper->initialize(); TEUCHOS_TEST_FOR_EXCEPT(!stepper->isInitialized());
75 stepper->setICConsistency(ICConsistency); stepper->initialize(); TEUCHOS_TEST_FOR_EXCEPT(!stepper->isInitialized());
76 stepper->setICConsistencyCheck(ICConsistencyCheck); stepper->initialize(); TEUCHOS_TEST_FOR_EXCEPT(!stepper->isInitialized());
77
78
79 // Full argument list construction.
80 auto scIntegrator = Teuchos::rcp(new Tempus::IntegratorBasic<double>());
81 auto stepperFE = Tempus::createStepperForwardEuler(modelME, Teuchos::null);
82 scIntegrator->setStepper(stepperFE);
83 scIntegrator->setSolutionHistory(solutionHistory);
84 scIntegrator->initialize();
85
86 stepper = rcp(new Tempus::StepperSubcycling<double>(
87 model,scIntegrator, useFSAL, ICConsistency, ICConsistencyCheck,modifier));
88 TEUCHOS_TEST_FOR_EXCEPT(!stepper->isInitialized());
89
90 // Test stepper properties.
91 TEUCHOS_ASSERT(stepper->getOrder() == 1);
92}
93
94
95// ************************************************************
96// ************************************************************
97TEUCHOS_UNIT_TEST(Subcycling, MaxTimeStepDoesNotChangeDuring_takeStep)
98{
99 // Setup the stepper ----------------------------------------
100 auto model = rcp(new Tempus_Test::SinCosModel<double>());
101 auto modelME = rcp_dynamic_cast<const Thyra::ModelEvaluator<double> > (model);
102 auto stepper = rcp(new Tempus::StepperSubcycling<double>());
103 auto stepperBE = Tempus::createStepperBackwardEuler(modelME, Teuchos::null);
104 stepper->setSubcyclingStepper(stepperBE);
105
106 // Setup SolutionHistory ------------------------------------
107 auto inArgsIC = model->getNominalValues();
108 auto icSolution =rcp_const_cast<Thyra::VectorBase<double> >(inArgsIC.get_x());
109 auto icState = Tempus::createSolutionStateX(icSolution);
110 auto solutionHistory = rcp(new Tempus::SolutionHistory<double>());
111 solutionHistory->addState(icState);
112 solutionHistory->initWorkingState();
113
114 // Ensure ICs are consistent and stepper memory is set (e.g., xDot is set).
115 stepper->setInitialConditions(solutionHistory);
116 stepper->initialize();
117
118 // Test
119 stepper->setSubcyclingInitTimeStep(0.25);
120 stepper->setSubcyclingMaxTimeStep(0.5);
121 double maxTimeStep_Set = stepper->getSubcyclingMaxTimeStep();
122 stepper->takeStep(solutionHistory);
123 double maxTimeStep_After = stepper->getSubcyclingMaxTimeStep();
124
125 TEST_FLOATING_EQUALITY(maxTimeStep_Set, maxTimeStep_After, 1.0e-14 );
126}
127
128
129// ************************************************************
130// ************************************************************
131class StepperSubcyclingModifierTest
132 : virtual public Tempus::StepperSubcyclingModifierBase<double>
133{
134public:
135
137 StepperSubcyclingModifierTest()
138 : testBEGIN_STEP(false), testEND_STEP(false),
139 testCurrentValue(-0.99), testWorkingValue(-0.99),
140 testDt(1.5), testName("")
141 {}
143 virtual ~StepperSubcyclingModifierTest(){}
144
146 virtual void modify(
147 Teuchos::RCP<Tempus::SolutionHistory<double> > sh,
148 Teuchos::RCP<Tempus::StepperSubcycling<double> > stepper,
150 {
151 switch(actLoc) {
152 case StepperSubcyclingAppAction<double>::BEGIN_STEP:
153 {
154 testBEGIN_STEP = true;
155 auto x = sh->getCurrentState()->getX();
156 testCurrentValue = get_ele(*(x), 0);
157 testName = "Subcycling - Modifier";
158 stepper->setStepperName(testName);
159 break;
160 }
161 case StepperSubcyclingAppAction<double>::END_STEP:
162 {
163 testEND_STEP = true;
164 auto x = sh->getWorkingState()->getX();
165 testWorkingValue = get_ele(*(x), 0);
166 testDt = sh->getWorkingState()->getTimeStep()/10.0;
167 sh->getWorkingState()->setTimeStep(testDt);
168 break;
169 }
170 default:
171 TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error,
172 "Error - unknown action location.\n");
173 }
174 }
175
176 bool testBEGIN_STEP;
177 bool testEND_STEP;
178 double testCurrentValue;
179 double testWorkingValue;
180 double testDt;
181 std::string testName;
182};
183
184TEUCHOS_UNIT_TEST(Subcycling, AppAction_Modifier)
185{
186 // Setup the SinCosModel ------------------------------------
187 auto model = rcp(new Tempus_Test::SinCosModel<double>());
188 auto modelME = rcp_dynamic_cast<const Thyra::ModelEvaluator<double> > (model);
189
190 // Setup Stepper for field solve ----------------------------
191 auto stepper = rcp(new Tempus::StepperSubcycling<double>());
192 auto stepperFE = Tempus::createStepperForwardEuler(modelME, Teuchos::null);
193 auto modifier = rcp(new StepperSubcyclingModifierTest());
194 stepper->setAppAction(modifier);
195 stepper->setSubcyclingStepper(stepperFE);
196
197 stepper->setSubcyclingMinTimeStep (15);
198 stepper->setSubcyclingInitTimeStep (15.0);
199 stepper->setSubcyclingMaxTimeStep (15.0);
200 stepper->setSubcyclingMaxFailures (10);
201 stepper->setSubcyclingMaxConsecFailures(5);
202 stepper->setSubcyclingScreenOutputIndexInterval(1);
203 stepper->setSubcyclingPrintDtChanges(true);
204
205 // Setup TimeStepControl ------------------------------------
206 auto timeStepControl = rcp(new Tempus::TimeStepControl<double>());
207 timeStepControl->setInitIndex(0);
208 timeStepControl->setInitTime (0.0);
209 timeStepControl->setFinalTime(1.0);
210 timeStepControl->setInitTimeStep(15.0);
211 timeStepControl->initialize();
212
213 // Setup initial condition SolutionState --------------------
214 auto inArgsIC = model->getNominalValues();
215 auto icSolution = rcp_const_cast<Thyra::VectorBase<double> > (inArgsIC.get_x());
216 auto icState = Tempus::createSolutionStateX(icSolution);
217 icState->setTime (timeStepControl->getInitTime());;
218 icState->setIndex (timeStepControl->getInitIndex());
219 icState->setTimeStep(0.0); // dt for ICs are indicated by zero.
220 icState->setSolutionStatus(Tempus::Status::PASSED); // ICs are passing.
221
222 // Setup SolutionHistory ------------------------------------
223 auto solutionHistory = rcp(new Tempus::SolutionHistory<double>());
224 solutionHistory->setName("Forward States");
225 solutionHistory->setStorageType(Tempus::STORAGE_TYPE_STATIC);
226 solutionHistory->setStorageLimit(2);
227 solutionHistory->addState(icState);
228
229 // Ensure ICs are consistent and stepper memory is set (e.g., xDot is set).
230 stepper->setInitialConditions(solutionHistory);
231 stepper->initialize();
232
233 // Take one time step.
234 stepper->setInitialConditions(solutionHistory);
235 solutionHistory->initWorkingState();
236 solutionHistory->getWorkingState()->setTimeStep(15.0);
237 stepper->takeStep(solutionHistory);
238
239 // Testing that each ACTION_LOCATION has been called.
240 TEST_COMPARE(modifier->testBEGIN_STEP, ==, true);
241 TEST_COMPARE(modifier->testEND_STEP, ==, true);
242
243 // Testing that values can be set through the Modifier.
244 auto x = solutionHistory->getCurrentState()->getX();
245 TEST_FLOATING_EQUALITY(modifier->testCurrentValue, get_ele(*(x), 0), 1.0e-14);
246 x = solutionHistory->getWorkingState()->getX();
247 TEST_FLOATING_EQUALITY(modifier->testWorkingValue, get_ele(*(x), 0), 1.0e-14);
248 auto Dt = solutionHistory->getWorkingState()->getTimeStep();
249 TEST_FLOATING_EQUALITY(modifier->testDt, Dt, 1.0e-14);
250
251 TEST_COMPARE(modifier->testName, ==, "Subcycling - Modifier");
252}
253
254// ************************************************************
255// ************************************************************
256class StepperSubcyclingObserverTest
257 : virtual public Tempus::StepperSubcyclingObserverBase<double>
258{
259public:
260
262 StepperSubcyclingObserverTest()
263 : testBEGIN_STEP(false), testEND_STEP(false),
264 testCurrentValue(-0.99), testWorkingValue(-0.99),
265 testDt(15.0), testName("Subcyling")
266 {}
267
269 virtual ~StepperSubcyclingObserverTest(){}
270
272 virtual void observe(
273 Teuchos::RCP<const Tempus::SolutionHistory<double> > sh,
274 Teuchos::RCP<const Tempus::StepperSubcycling<double> > stepper,
276 {
277 switch(actLoc) {
278 case StepperSubcyclingAppAction<double>::BEGIN_STEP:
279 {
280 testBEGIN_STEP = true;
281 auto x = sh->getCurrentState()->getX();
282 testCurrentValue = get_ele(*(x), 0);
283 break;
284 }
285 case StepperSubcyclingAppAction<double>::END_STEP:
286 {
287 testEND_STEP = true;
288 auto x = sh->getWorkingState()->getX();
289 testWorkingValue = get_ele(*(x), 0);
290 break;
291 }
292 default:
293 TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error,
294 "Error - unknown action location.\n");
295 }
296 }
297
298 bool testBEGIN_STEP;
299 bool testEND_STEP;
300 double testCurrentValue;
301 double testWorkingValue;
302 double testDt;
303 std::string testName;
304};
305
306TEUCHOS_UNIT_TEST(Subcycling, AppAction_Observer)
307{
308 // Setup the SinCosModel ------------------------------------
309 auto model = rcp(new Tempus_Test::SinCosModel<double>());
310 auto modelME = rcp_dynamic_cast<const Thyra::ModelEvaluator<double> > (model);
311
312 // Setup Stepper for field solve ----------------------------
313 auto stepper = rcp(new Tempus::StepperSubcycling<double>());
314 auto stepperFE = Tempus::createStepperForwardEuler(modelME, Teuchos::null);
315 auto observer = rcp(new StepperSubcyclingObserverTest());
316 stepper->setAppAction(observer);
317 stepper->setSubcyclingStepper(stepperFE);
318
319 stepper->setSubcyclingMinTimeStep (15);
320 stepper->setSubcyclingInitTimeStep (15.0);
321 stepper->setSubcyclingMaxTimeStep (15.0);
322 stepper->setSubcyclingMaxFailures (10);
323 stepper->setSubcyclingMaxConsecFailures(5);
324 stepper->setSubcyclingScreenOutputIndexInterval(1);
325 stepper->setSubcyclingPrintDtChanges(true);
326
327 // Setup TimeStepControl ------------------------------------
328 auto timeStepControl = rcp(new Tempus::TimeStepControl<double>());
329 timeStepControl->setInitIndex(0);
330 timeStepControl->setInitTime (0.0);
331 timeStepControl->setFinalTime(1.0);
332 timeStepControl->setInitTimeStep(15.0);
333 timeStepControl->initialize();
334
335 // Setup initial condition SolutionState --------------------
336 auto inArgsIC = model->getNominalValues();
337 auto icSolution = rcp_const_cast<Thyra::VectorBase<double> > (inArgsIC.get_x());
338 auto icState = Tempus::createSolutionStateX(icSolution);
339 icState->setTime (timeStepControl->getInitTime());;
340 icState->setIndex (timeStepControl->getInitIndex());
341 icState->setTimeStep(0.0); // dt for ICs are indicated by zero.
342 icState->setSolutionStatus(Tempus::Status::PASSED); // ICs are passing.
343
344 // Setup SolutionHistory ------------------------------------
345 auto solutionHistory = rcp(new Tempus::SolutionHistory<double>());
346 solutionHistory->setName("Forward States");
347 solutionHistory->setStorageType(Tempus::STORAGE_TYPE_STATIC);
348 solutionHistory->setStorageLimit(2);
349 solutionHistory->addState(icState);
350
351 // Ensure ICs are consistent and stepper memory is set (e.g., xDot is set).
352 stepper->setInitialConditions(solutionHistory);
353 stepper->initialize();
354
355 // Take one time step.
356 stepper->setInitialConditions(solutionHistory);
357 solutionHistory->initWorkingState();
358 solutionHistory->getWorkingState()->setTimeStep(15.0);
359 stepper->takeStep(solutionHistory);
360
361 // Testing that each ACTION_LOCATION has been called.
362 TEST_COMPARE(observer->testBEGIN_STEP, ==, true);
363 TEST_COMPARE(observer->testEND_STEP, ==, true);
364
365 // Testing that values can be observed through the observer.
366 auto x = solutionHistory->getCurrentState()->getX();
367 TEST_FLOATING_EQUALITY(observer->testCurrentValue, get_ele(*(x), 0), 1.0e-14);
368 x = solutionHistory->getWorkingState()->getX();
369 TEST_FLOATING_EQUALITY(observer->testWorkingValue, get_ele(*(x), 0), 1.0e-14);
370 TEST_FLOATING_EQUALITY(observer->testDt, 15.0, 1.0e-14);
371
372 TEST_COMPARE(observer->testName, ==, "Subcyling");
373}
374
375 // ************************************************************
376 // ************************************************************
377class StepperSubcyclingModifierXTest
378 : virtual public Tempus::StepperSubcyclingModifierXBase<double>
379{
380public:
381
383 StepperSubcyclingModifierXTest()
384 : testX_BEGIN_STEP(false), testXDOT_END_STEP(false),
385 testX(-0.99), testXDot(-0.99),
386 testDt(1.5), testTime(1.5)
387 {}
388
390 virtual ~StepperSubcyclingModifierXTest(){}
391
393 virtual void modify(
394 Teuchos::RCP<Thyra::VectorBase<double> > x,
395 const double time, const double dt,
397 {
398 switch(modType) {
399 case StepperSubcyclingModifierXBase<double>::X_BEGIN_STEP:
400 {
401 testX_BEGIN_STEP = true;
402 testX = get_ele(*(x), 0);
403 testDt = dt;
404 break;
405 }
406 case StepperSubcyclingModifierXBase<double>::XDOT_END_STEP:
407 {
408 testXDOT_END_STEP = true;
409 testXDot = get_ele(*(x), 0);
410 testTime = time;
411 break;
412 }
413 default:
414 TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error,
415 "Error - unknown action location.\n");
416 }
417 }
418
419 bool testX_BEGIN_STEP;
420 bool testXDOT_END_STEP;
421 double testX;
422 double testXDot;
423 double testDt;
424 double testTime;
425};
426
427TEUCHOS_UNIT_TEST(Subcycling, AppAction_ModifierX)
428{
429 // Setup the SinCosModel ------------------------------------
430 auto model = rcp(new Tempus_Test::SinCosModel<double>());
431 auto modelME = rcp_dynamic_cast<const Thyra::ModelEvaluator<double> > (model);
432
433 // Setup Stepper for field solve ----------------------------
434 auto stepper = rcp(new Tempus::StepperSubcycling<double>());
435 auto stepperFE = Tempus::createStepperForwardEuler(modelME, Teuchos::null);
436 auto modifierX = rcp(new StepperSubcyclingModifierXTest());
437 stepper->setAppAction(modifierX);
438 stepper->setSubcyclingStepper(stepperFE);
439
440 stepper->setSubcyclingMinTimeStep (15);
441 stepper->setSubcyclingInitTimeStep (15.0);
442 stepper->setSubcyclingMaxTimeStep (15.0);
443 stepper->setSubcyclingMaxFailures (10);
444 stepper->setSubcyclingMaxConsecFailures(5);
445 stepper->setSubcyclingScreenOutputIndexInterval(1);
446 stepper->setSubcyclingPrintDtChanges(true);
447
448 // Setup TimeStepControl ------------------------------------
449 auto timeStepControl = rcp(new Tempus::TimeStepControl<double>());
450 timeStepControl->setInitIndex(0);
451 timeStepControl->setInitTime (0.0);
452 timeStepControl->setFinalTime(1.0);
453 timeStepControl->setInitTimeStep(15.0);
454 timeStepControl->initialize();
455
456 // Setup initial condition SolutionState --------------------
457 auto inArgsIC = model->getNominalValues();
458 auto icSolution = rcp_const_cast<Thyra::VectorBase<double> > (inArgsIC.get_x());
459 auto icSolutionDot = rcp_const_cast<Thyra::VectorBase<double> > (inArgsIC.get_x_dot());
460 auto icState = Tempus::createSolutionStateX(icSolution,icSolutionDot);
461 icState->setTime (timeStepControl->getInitTime());;
462 icState->setIndex (timeStepControl->getInitIndex());
463 icState->setTimeStep(0.0); // dt for ICs are indicated by zero.
464 icState->setSolutionStatus(Tempus::Status::PASSED); // ICs are passing.
465
466 // Setup SolutionHistory ------------------------------------
467 auto solutionHistory = rcp(new Tempus::SolutionHistory<double>());
468 solutionHistory->setName("Forward States");
469 solutionHistory->setStorageType(Tempus::STORAGE_TYPE_STATIC);
470 solutionHistory->setStorageLimit(2);
471 solutionHistory->addState(icState);
472
473 // Ensure ICs are consistent and stepper memory is set (e.g., xDot is set).
474 stepper->setInitialConditions(solutionHistory);
475 stepper->initialize();
476
477 // Take one time step.
478 stepper->setInitialConditions(solutionHistory);
479 solutionHistory->initWorkingState();
480 solutionHistory->getWorkingState()->setTimeStep(15.0);
481 stepper->takeStep(solutionHistory);
482
483 // Take one time step.
484 stepper->setInitialConditions(solutionHistory);
485 solutionHistory->initWorkingState();
486 solutionHistory->getWorkingState()->setTimeStep(15.0);
487 stepper->takeStep(solutionHistory);
488
489 // Testing that each ACTION_LOCATION has been called.
490 TEST_COMPARE(modifierX->testX_BEGIN_STEP, ==, true);
491 TEST_COMPARE(modifierX->testXDOT_END_STEP, ==, true);
492
493 // Testing that values can be set through the Modifier.
494 auto x = solutionHistory->getCurrentState()->getX();
495 TEST_FLOATING_EQUALITY(modifierX->testX, get_ele(*(x), 0), 1.0e-14);
496 // Temporary memory for xDot is not guarranteed to exist outside the Stepper.
497 auto xDot = solutionHistory->getWorkingState()->getXDot();
498 if (xDot == Teuchos::null) xDot = stepper->getStepperXDot();
499
500 TEST_FLOATING_EQUALITY(modifierX->testXDot, get_ele(*(xDot), 0),1.0e-14);
501 auto Dt = solutionHistory->getWorkingState()->getTimeStep();
502 TEST_FLOATING_EQUALITY(modifierX->testDt, Dt, 1.0e-14);
503
504 auto time = solutionHistory->getWorkingState()->getTime();
505 TEST_FLOATING_EQUALITY(modifierX->testTime, time, 1.0e-14);
506}
507
508} // namespace Tempus_Test
SolutionHistory is basically a container of SolutionStates. SolutionHistory maintains a collection of...
ACTION_LOCATION
Indicates the location of application action (see algorithm).
MODIFIER_TYPE
Indicates the location of application action (see algorithm).
TimeStepControl manages the time step size. There several mechanisms that effect the time step size a...
TEUCHOS_UNIT_TEST(BackwardEuler, Default_Construction)
@ STORAGE_TYPE_STATIC
Keep a fix number of states.
Teuchos::RCP< StepperBackwardEuler< Scalar > > createStepperBackwardEuler(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &model, Teuchos::RCP< Teuchos::ParameterList > pl)
Nonmember constructor - ModelEvaluator and ParameterList.
Teuchos::RCP< SolutionState< Scalar > > createSolutionStateX(const Teuchos::RCP< Thyra::VectorBase< Scalar > > &x, const Teuchos::RCP< Thyra::VectorBase< Scalar > > &xdot=Teuchos::null, const Teuchos::RCP< Thyra::VectorBase< Scalar > > &xdotdot=Teuchos::null)
Nonmember constructor from non-const solution vectors, x.
Teuchos::RCP< StepperForwardEuler< Scalar > > createStepperForwardEuler(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &model, Teuchos::RCP< Teuchos::ParameterList > pl)
Nonmember constructor - ModelEvaluator and ParameterList.
Teuchos::RCP< Teuchos::ParameterList > defaultSolverParameters()
Returns the default solver ParameterList for implicit Steppers.