Tempus  Version of the Day
Time Integration
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Tempus_UnitTest_BDF2.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 
12 
13 #include "Tempus_StepperBDF2.hpp"
18 
19 
20 namespace Tempus_Unit_Test {
21 
22 using Teuchos::RCP;
23 using Teuchos::rcp;
24 using Teuchos::rcp_const_cast;
25 using Teuchos::rcp_dynamic_cast;
27 using Teuchos::sublist;
28 
29 
30 // ************************************************************
31 // ************************************************************
32 TEUCHOS_UNIT_TEST(BDF2, Default_Construction)
33 {
34  auto model = rcp(new Tempus_Test::SinCosModel<double>());
35 
36  // Default construction.
37  auto stepper = rcp(new Tempus::StepperBDF2<double>());
38  stepper->setModel(model);
39  stepper->initialize();
40  TEUCHOS_TEST_FOR_EXCEPT(!stepper->isInitialized());
41 
42 
43  // Default values for construction.
44  auto modifier = rcp(new Tempus::StepperBDF2ModifierDefault<double>());
45  auto solver = rcp(new Thyra::NOXNonlinearSolver());
46  solver->setParameterList(Tempus::defaultSolverParameters());
47 
48  auto startUpStepper = rcp(new Tempus::StepperDIRK_1StageTheta<double>());
49  startUpStepper->setModel(model); // Can use the same model since both steppers are implicit ODEs.
50  startUpStepper->initialize();
51 
52  auto defaultStepper = rcp(new Tempus::StepperBDF2<double>());
53  bool useFSAL = defaultStepper->getUseFSAL();
54  std::string ICConsistency = defaultStepper->getICConsistency();
55  bool ICConsistencyCheck = defaultStepper->getICConsistencyCheck();
56  bool zeroInitialGuess = defaultStepper->getZeroInitialGuess();
57 
58  // Test the set functions.
59  stepper->setAppAction(modifier); stepper->initialize(); TEUCHOS_TEST_FOR_EXCEPT(!stepper->isInitialized());
60  stepper->setSolver(solver); stepper->initialize(); TEUCHOS_TEST_FOR_EXCEPT(!stepper->isInitialized());
61  stepper->setStartUpStepper(startUpStepper); stepper->initialize(); TEUCHOS_TEST_FOR_EXCEPT(!stepper->isInitialized());
62  stepper->setUseFSAL(useFSAL); stepper->initialize(); TEUCHOS_TEST_FOR_EXCEPT(!stepper->isInitialized());
63  stepper->setICConsistency(ICConsistency); stepper->initialize(); TEUCHOS_TEST_FOR_EXCEPT(!stepper->isInitialized());
64  stepper->setICConsistencyCheck(ICConsistencyCheck); stepper->initialize(); TEUCHOS_TEST_FOR_EXCEPT(!stepper->isInitialized());
65  stepper->setZeroInitialGuess(zeroInitialGuess); stepper->initialize(); TEUCHOS_TEST_FOR_EXCEPT(!stepper->isInitialized());
66 
67  stepper = rcp(new Tempus::StepperBDF2<double>(model, solver, startUpStepper, useFSAL,
68  ICConsistency, ICConsistencyCheck, zeroInitialGuess,modifier));
69  TEUCHOS_TEST_FOR_EXCEPT(!stepper->isInitialized());
70  // Test stepper properties.
71  TEUCHOS_ASSERT(stepper->getOrder() == 2);
72 }
73 
74 
75 // ************************************************************
76 // ************************************************************
77 TEUCHOS_UNIT_TEST(BDF2, StepperFactory_Construction)
78 {
79  auto model = rcp(new Tempus_Test::SinCosModel<double>());
80  testFactoryConstruction("BDF2", model);
81 }
82 
83 
84 // ************************************************************
85 // ************************************************************
86 class StepperBDF2ModifierTest
87  : virtual public Tempus::StepperBDF2ModifierBase<double>
88 {
89 public:
90 
92  StepperBDF2ModifierTest()
93  : testBEGIN_STEP(false),testBEFORE_SOLVE(false),
94  testAFTER_SOLVE(false),testEND_STEP(false),
95  testCurrentValue(-0.99), testWorkingValue(-0.99),
96  testDt(.99), testType("")
97  {}
98 
100  virtual ~StepperBDF2ModifierTest(){}
101 
103  virtual void modify(Teuchos::RCP<Tempus::SolutionHistory<double> > sh,
106  {
107  switch(actLoc) {
109  {
110  testBEGIN_STEP = true;
111  auto x = sh->getWorkingState()->getX();
112  testCurrentValue = get_ele(*(x), 0);
113  break;
114  }
116  {
117  testBEFORE_SOLVE = true;
118  testType = "BDF2 - Modifier";
119  break;
120  }
122  {
123  testAFTER_SOLVE = true;
124  testDt = sh->getCurrentState()->getTimeStep()/10.0;
125  sh->getCurrentState()->setTimeStep(testDt);
126  break;
127  }
129  {
130  testEND_STEP = true;
131  auto x = sh->getWorkingState()->getX();
132  testWorkingValue = get_ele(*(x), 0);
133  break;
134  }
135  default:
136  TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error,
137  "Error - unknown action location.\n");
138  }
139  }
140  bool testBEGIN_STEP;
141  bool testBEFORE_SOLVE;
142  bool testAFTER_SOLVE;
143  bool testEND_STEP;
144  double testCurrentValue;
145  double testWorkingValue;
146  double testDt;
147  std::string testType;
148 };
149 
150 TEUCHOS_UNIT_TEST(BDF2, AppAction_Modifier)
151 {
152  auto model = rcp(new Tempus_Test::SinCosModel<double>());
153 
154  // Setup Stepper for field solve ----------------------------
155  auto stepper = rcp(new Tempus::StepperBDF2<double>());
156  stepper->setModel(model);
157  auto modifier = rcp(new StepperBDF2ModifierTest());
158  stepper->setAppAction(modifier);
159  stepper->initialize();
160 
161  // Setup initial condition SolutionState --------------------
162  auto inArgsIC = model->getNominalValues();
163  auto icSolution =
164  rcp_const_cast<Thyra::VectorBase<double> > (inArgsIC.get_x());
165  auto icState = Tempus::createSolutionStateX(icSolution);
166  icState->setTime (0.0);
167  icState->setIndex (0);
168  icState->setTimeStep(1.0);
169  icState->setOrder (stepper->getOrder());
170  icState->setSolutionStatus(Tempus::Status::PASSED); // ICs are passing.
171 
172  // Setup TimeStepControl ------------------------------------
173  auto timeStepControl = rcp(new Tempus::TimeStepControl<double>());
174  timeStepControl->setInitIndex(0);
175  timeStepControl->setInitTime (0.0);
176  timeStepControl->setFinalTime(2.0);
177  timeStepControl->setInitTimeStep(1.0);
178  timeStepControl->initialize();
179 
180  // Setup SolutionHistory ------------------------------------
181  auto solutionHistory = rcp(new Tempus::SolutionHistory<double>());
182  solutionHistory->setName("Forward States");
183  solutionHistory->setStorageType(Tempus::STORAGE_TYPE_STATIC);
184  solutionHistory->setStorageLimit(3);
185  solutionHistory->addState(icState);
186 
187  // Take two time steps (the first step will not test BDF2's modifier)
188  stepper->setInitialConditions(solutionHistory);
189  solutionHistory->initWorkingState();
190  double dt = 1.0;
191  solutionHistory->getWorkingState()->setTimeStep(dt);
192  stepper->takeStep(solutionHistory);
193  solutionHistory->promoteWorkingState();
194  solutionHistory->initWorkingState();
195  stepper->takeStep(solutionHistory);
196 
197  TEST_COMPARE(modifier->testBEGIN_STEP, ==, true);
198  TEST_COMPARE(modifier->testBEFORE_SOLVE, ==, true);
199  TEST_COMPARE(modifier->testAFTER_SOLVE, ==, true);
200  TEST_COMPARE(modifier->testEND_STEP, ==, true);
201 
202  auto Dt = solutionHistory->getCurrentState()->getTimeStep();
203  TEST_FLOATING_EQUALITY(modifier->testDt, Dt, 1.0e-15);
204  auto x = solutionHistory->getCurrentState()->getX();
205  TEST_FLOATING_EQUALITY(modifier->testCurrentValue, get_ele(*(x), 0), 1.0e-15);
206  x = solutionHistory->getWorkingState()->getX();
207  TEST_FLOATING_EQUALITY(modifier->testWorkingValue, get_ele(*(x), 0), 1.0e-15);
208  TEST_COMPARE(modifier->testType, ==, "BDF2 - Modifier");
209 }
210 
211 
212 // ************************************************************
213 // ************************************************************
214 class StepperBDF2ObserverTest
215  : virtual public Tempus::StepperBDF2ObserverBase<double>
216 {
217 public:
218 
220  StepperBDF2ObserverTest()
221  : testBEGIN_STEP(false),testBEFORE_SOLVE(false),
222  testAFTER_SOLVE(false),testEND_STEP(false),
223  testCurrentValue(0.99), testWorkingValue(0.99),
224  testDt(.99), testType("")
225  {}
226 
228  virtual ~StepperBDF2ObserverTest(){}
229 
231  virtual void modify(Teuchos::RCP<Tempus::SolutionHistory<double> > sh,
234  {
235  switch(actLoc) {
237  {
238  testBEGIN_STEP = true;
239  auto x = sh->getCurrentState()->getX();
240  testCurrentValue = get_ele(*(x), 0);
241  break;
242  }
244  {
245  testBEFORE_SOLVE = true;
246  testType = stepper->getStepperType();
247  break;
248  }
250  {
251  testAFTER_SOLVE = true;
252  testDt = sh->getCurrentState()->getTimeStep();
253  break;
254  }
256  {
257  testEND_STEP = true;
258  auto x = sh->getWorkingState()->getX();
259  testWorkingValue = get_ele(*(x), 0);
260  break;
261  }
262 
263  default:
264  TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error,
265  "Error - unknown action location.\n");
266  }
267  }
268  bool testBEGIN_STEP;
269  bool testBEFORE_SOLVE;
270  bool testAFTER_SOLVE;
271  bool testEND_STEP;
272  double testCurrentValue;
273  double testWorkingValue;
274  double testDt;
275  std::string testType;
276 };
277 
278 
279 // ************************************************************
280 // ************************************************************
281 TEUCHOS_UNIT_TEST(BDF2, AppAction_Observer)
282 {
283  auto model = rcp(new Tempus_Test::SinCosModel<double>());
284 
285  // Setup Stepper for field solve ----------------------------
286  auto stepper = rcp(new Tempus::StepperBDF2<double>());
287  stepper->setModel(model);
288  auto observer = rcp(new StepperBDF2ModifierTest());
289  stepper->setAppAction(observer);
290  stepper->initialize();
291 
292  // Setup initial condition SolutionState --------------------
293  auto inArgsIC = model->getNominalValues();
294  auto icSolution =
295  rcp_const_cast<Thyra::VectorBase<double> > (inArgsIC.get_x());
296  auto icState = Tempus::createSolutionStateX(icSolution);
297  icState->setTime (0.0);
298  icState->setIndex (0);
299  icState->setTimeStep(1.0);
300  icState->setOrder (stepper->getOrder());
301  icState->setSolutionStatus(Tempus::Status::PASSED); // ICs are passing.
302 
303  // Setup TimeStepControl ------------------------------------
304  auto timeStepControl = rcp(new Tempus::TimeStepControl<double>());
305  timeStepControl->setInitIndex(0);
306  timeStepControl->setInitTime (0.0);
307  timeStepControl->setFinalTime(2.0);
308  timeStepControl->setInitTimeStep(1.0);
309  timeStepControl->initialize();
310 
311  // Setup SolutionHistory ------------------------------------
312  auto solutionHistory = rcp(new Tempus::SolutionHistory<double>());
313  solutionHistory->setName("Forward States");
314  solutionHistory->setStorageType(Tempus::STORAGE_TYPE_STATIC);
315  solutionHistory->setStorageLimit(3);
316  solutionHistory->addState(icState);
317 
318  // Take two time steps (the first step will not test BDF2's modifier)
319  stepper->setInitialConditions(solutionHistory);
320  solutionHistory->initWorkingState();
321  double dt = 1.0;
322  solutionHistory->getWorkingState()->setTimeStep(dt);
323  stepper->takeStep(solutionHistory);
324  solutionHistory->promoteWorkingState();
325  solutionHistory->initWorkingState();
326  stepper->takeStep(solutionHistory);
327 
328  TEST_COMPARE(observer->testBEGIN_STEP, ==, true);
329  TEST_COMPARE(observer->testBEFORE_SOLVE, ==, true);
330  TEST_COMPARE(observer->testAFTER_SOLVE, ==, true);
331  TEST_COMPARE(observer->testEND_STEP, ==, true);
332 
333  auto Dt = solutionHistory->getCurrentState()->getTimeStep();
334  TEST_FLOATING_EQUALITY(observer->testDt, Dt, 1.0e-15);
335 
336  auto x = solutionHistory->getCurrentState()->getX();
337  TEST_FLOATING_EQUALITY(observer->testCurrentValue, get_ele(*(x), 0), 1.0e-15);
338  x = solutionHistory->getWorkingState()->getX();
339  TEST_FLOATING_EQUALITY(observer->testWorkingValue, get_ele(*(x), 0), 1.0e-15);
340  TEST_COMPARE(observer->testType, ==, "BDF2 - Modifier");
341 }
342 
343 
344 // ************************************************************
345 // ************************************************************
346 class StepperBDF2ModifierXTest
347  : virtual public Tempus::StepperBDF2ModifierXBase<double>
348 {
349 public:
350 
352  StepperBDF2ModifierXTest()
353  : testX_BEGIN_STEP(false),testX_BEFORE_SOLVE(false),
354  testX_AFTER_SOLVE(false),testX_END_STEP(false),
355  testXbegin(-.99),testXend(-.99),testTime(0.0),testDt(0.0)
356  {}
357 
359  virtual ~StepperBDF2ModifierXTest(){}
360 
362  virtual void modify(
364  const double time, const double dt,
366  {
367  switch(modType) {
369  {
370  testX_BEGIN_STEP = true;
371  testXbegin = get_ele(*(x), 0);
372  break;
373  }
375  {
376  testX_BEFORE_SOLVE = true;
377  testDt = dt;
378  break;
379  }
381  {
382  testX_AFTER_SOLVE = true;
383  testTime = time;
384  break;
385  }
387  {
388  testX_END_STEP = true;
389  testXend = get_ele(*(x), 0);
390  break;
391  }
392  default:
393  TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error,
394  "Error - unknown action location.\n");
395  }
396  }
397  bool testX_BEGIN_STEP;
398  bool testX_BEFORE_SOLVE;
399  bool testX_AFTER_SOLVE;
400  bool testX_END_STEP;
401  double testXbegin;
402  double testXend;
403  double testTime;
404  double testDt;
405 };
406 
407 
408 // ************************************************************
409 // ************************************************************
410 TEUCHOS_UNIT_TEST(BDF2, AppAction_ModifierX)
411 {
412  auto model = rcp(new Tempus_Test::SinCosModel<double>());
413  // Setup Stepper for field solve ----------------------------
414  auto stepper = rcp(new Tempus::StepperBDF2<double>());
415  stepper->setModel(model);
416  auto modifierX = rcp(new StepperBDF2ModifierXTest());
417  stepper->setAppAction(modifierX);
418  stepper->initialize();
419 
420  // Setup initial condition SolutionState --------------------
421  auto inArgsIC = model->getNominalValues();
422  auto icSolution =
423  rcp_const_cast<Thyra::VectorBase<double> > (inArgsIC.get_x());
424  auto icState = Tempus::createSolutionStateX(icSolution);
425  icState->setTime (0.0);
426  icState->setIndex (0);
427  icState->setTimeStep(1.0);
428  icState->setOrder (stepper->getOrder());
429  icState->setSolutionStatus(Tempus::Status::PASSED); // ICs are passing.
430 
431  // Setup TimeStepControl ------------------------------------
432  auto timeStepControl = rcp(new Tempus::TimeStepControl<double>());
433  timeStepControl->setInitIndex(0);
434  timeStepControl->setInitTime (0.0);
435  timeStepControl->setFinalTime(2.0);
436  timeStepControl->setInitTimeStep(1.0);
437  timeStepControl->initialize();
438 
439  // Setup SolutionHistory ------------------------------------
440  auto solutionHistory = rcp(new Tempus::SolutionHistory<double>());
441  solutionHistory->setName("Forward States");
442  solutionHistory->setStorageType(Tempus::STORAGE_TYPE_STATIC);
443  solutionHistory->setStorageLimit(3);
444  solutionHistory->addState(icState);
445 
446 
447  // Take two time steps (the first step will not test BDF2's modifier)
448  stepper->setInitialConditions(solutionHistory);
449  solutionHistory->initWorkingState();
450  double dt = 1.0;
451  solutionHistory->getWorkingState()->setTimeStep(dt);
452  stepper->takeStep(solutionHistory);
453  solutionHistory->promoteWorkingState();
454  solutionHistory->initWorkingState();
455  stepper->takeStep(solutionHistory);
456 
457  TEST_COMPARE(modifierX->testX_BEGIN_STEP, ==, true);
458  TEST_COMPARE(modifierX->testX_BEFORE_SOLVE, ==, true);
459  TEST_COMPARE(modifierX->testX_AFTER_SOLVE, ==, true);
460  TEST_COMPARE(modifierX->testX_END_STEP, ==, true);
461 
462  // Testing that values can be set through the Modifier.
463  auto x = solutionHistory->getCurrentState()->getX();
464  TEST_FLOATING_EQUALITY(modifierX->testXbegin, get_ele(*(x), 0), 1.0e-15);
465  auto Dt = solutionHistory->getWorkingState()->getTimeStep();
466  x = solutionHistory->getWorkingState()->getX();
467  TEST_FLOATING_EQUALITY(modifierX->testXend, get_ele(*(x), 0), 1.0e-15);
468  TEST_FLOATING_EQUALITY(modifierX->testDt, Dt, 1.0e-15);
469  auto time = solutionHistory->getWorkingState()->getTime();
470  TEST_FLOATING_EQUALITY(modifierX->testTime, time, 1.0e-15);
471 }
472 
473 
474 } // namespace Tempus_Unit_Test
BDF2 (Backward-Difference-Formula-2) time stepper.
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.
ACTION_LOCATION
Indicates the location of application action (see algorithm).
void testFactoryConstruction(std::string stepperType, const Teuchos::RCP< const Thyra::ModelEvaluator< double > > &model)
Unit test utility for Stepper construction through StepperFactory.
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
#define TEST_COMPARE(v1, comp, v2)
#define TEST_FLOATING_EQUALITY(v1, v2, tol)
Sine-Cosine model problem from Rythmos. This is a canonical Sine-Cosine differential equation with a...
Teuchos::RCP< Teuchos::ParameterList > defaultSolverParameters()
Returns the default solver ParameterList for implicit Steppers.
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
virtual void modify(Teuchos::RCP< Thyra::VectorBase< double > >, const double, const double, const MODIFIER_TYPE modType)=0
Modify solution based on the MODIFIER_TYPE.
TEUCHOS_UNIT_TEST(BackwardEuler, Default_Construction)
TimeStepControl manages the time step size. There several mechanisms that effect the time step size a...
SolutionHistory is basically a container of SolutionStates. SolutionHistory maintains a collection of...
Keep a fix number of states.
#define TEUCHOS_ASSERT(assertion_test)
MODIFIER_TYPE
Indicates the location of application action (see algorithm).
#define TEUCHOS_TEST_FOR_EXCEPT(throw_exception_test)