Tempus  Version of the Day
Time Integration
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Tempus_UnitTest_NewmarkImplicitDForm.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 
9 #include "Teuchos_UnitTestHarness.hpp"
10 #include "Teuchos_XMLParameterListHelpers.hpp"
11 #include "Teuchos_TimeMonitor.hpp"
12 #include "Teuchos_DefaultComm.hpp"
13 
14 #include "Thyra_VectorStdOps.hpp"
15 
16 #include "Tempus_IntegratorBasic.hpp"
17 #include "Tempus_StepperFactory.hpp"
18 #include "Tempus_UnitTest_Utils.hpp"
19 #include "Tempus_TimeStepControl.hpp"
20 
21 #include "Tempus_StepperNewmarkImplicitDForm.hpp"
22 #include "Tempus_StepperNewmarkImplicitDFormModifierBase.hpp"
23 #include "Tempus_StepperNewmarkImplicitDFormModifierXBase.hpp"
24 #include "Tempus_StepperNewmarkImplicitDFormModifierDefault.hpp"
25 #include "Tempus_StepperNewmarkImplicitDFormModifierXDefault.hpp"
26 
27 #include "../TestModels/SinCosModel.hpp"
28 #include "../TestModels/VanDerPolModel.hpp"
29 #include "../TestModels/HarmonicOscillatorModel.hpp"
30 #include "../TestUtils/Tempus_ConvergenceTestUtils.hpp"
31 
32 #include <fstream>
33 #include <vector>
34 
35 namespace Tempus_Unit_Test {
36 
37 using Teuchos::RCP;
38 using Teuchos::rcp;
39 using Teuchos::rcp_const_cast;
40 using Teuchos::rcp_dynamic_cast;
41 using Teuchos::ParameterList;
42 using Teuchos::sublist;
43 using Teuchos::getParametersFromXmlFile;
44 
45 using Tempus::StepperFactory;
46 
47 // ************************************************************
48 // ************************************************************
49 class StepperNewmarkImplicitDFormModifierTest
50  : virtual public Tempus::StepperNewmarkImplicitDFormModifierBase<double>
51 {
52 public:
53 
55  StepperNewmarkImplicitDFormModifierTest()
56  : testBEGIN_STEP(false), testBEFORE_SOLVE(false),
57  testAFTER_SOLVE(false), testEND_STEP(false),
58  testCurrentValue(-0.99),
59  testDt(-1.5), testName("")
60  {}
61 
63  virtual ~StepperNewmarkImplicitDFormModifierTest(){}
64 
66  virtual void modify(
67  Teuchos::RCP<Tempus::SolutionHistory<double> > sh,
68  Teuchos::RCP<Tempus::StepperNewmarkImplicitDForm<double> > stepper,
69  const typename Tempus::StepperNewmarkImplicitDFormAppAction<double>::ACTION_LOCATION actLoc)
70  {
71  switch(actLoc) {
72  case StepperNewmarkImplicitDFormAppAction<double>::BEGIN_STEP:
73  {
74  testBEGIN_STEP = true;
75  break;
76  }
77  case StepperNewmarkImplicitDFormAppAction<double>::BEFORE_SOLVE:
78  {
79  testBEFORE_SOLVE = true;
80  testDt = sh->getWorkingState()->getTimeStep();
81  //sh->getWorkingState()->setTimeStep(testDt);
82  break;
83  }
84  case StepperNewmarkImplicitDFormAppAction<double>::AFTER_SOLVE:
85  {
86  testAFTER_SOLVE = true;
87  testName = "Newmark Implicit A Form - Modifier";
88  stepper->setStepperName(testName);
89  break;
90  }
91  case StepperNewmarkImplicitDFormAppAction<double>::END_STEP:
92  {
93  testEND_STEP = true;
94  auto x = sh->getWorkingState()->getX();
95  testCurrentValue = get_ele(*(x), 0);
96  break;
97  }
98  default:
99  TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error,
100  "Error - unknown action location.\n");
101  }
102  }
103 
104  bool testBEGIN_STEP;
105  bool testBEFORE_SOLVE;
106  bool testAFTER_SOLVE;
107  bool testEND_STEP;
108  double testCurrentValue;
109  double testDt;
110  std::string testName;
111 };
112 
113 
114 // ************************************************************
115 // ************************************************************
116 class StepperNewmarkImplicitDFormModifierXTest
117  : virtual public Tempus::StepperNewmarkImplicitDFormModifierXBase<double>
118 {
119 public:
120 
122  StepperNewmarkImplicitDFormModifierXTest()
123  : testX_BEGIN_STEP(false), testX_BEFORE_SOLVE(false),
124  testX_AFTER_SOLVE(false), testX_END_STEP(false),
125  testX(-0.99), testXDot(-0.99),
126  testDt(-1.5), testTime(-1.5)
127  {}
128 
130  virtual ~StepperNewmarkImplicitDFormModifierXTest(){}
131 
133  virtual void modify(
134  Teuchos::RCP<Thyra::VectorBase<double> > x,
135  const double time, const double dt,
136  const typename Tempus::StepperNewmarkImplicitDFormModifierXBase<double>::MODIFIER_TYPE modType)
137  {
138  switch(modType) {
139  case StepperNewmarkImplicitDFormModifierXBase<double>::X_BEGIN_STEP:
140  {
141  testX_BEGIN_STEP = true;
142  testX = get_ele(*(x), 0);
143  break;
144  }
145  case StepperNewmarkImplicitDFormModifierXBase<double>::X_BEFORE_SOLVE:
146  {
147  testX_BEFORE_SOLVE = true;
148  testDt = dt;
149  break;
150  }
151  case StepperNewmarkImplicitDFormModifierXBase<double>::X_AFTER_SOLVE:
152  {
153  testX_AFTER_SOLVE = true;
154  testTime = time;
155  testX = get_ele(*(x), 0);
156  break;
157  }
158  case StepperNewmarkImplicitDFormModifierXBase<double>::X_END_STEP:
159  {
160  testX_END_STEP = true;
161  testX = get_ele(*(x), 0);
162  break;
163  }
164  default:
165  TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error,
166  "Error - unknown action location.\n");
167  }
168  }
169 
170  bool testX_BEGIN_STEP;
171  bool testX_BEFORE_SOLVE;
172  bool testX_AFTER_SOLVE;
173  bool testX_END_STEP;
174  double testX;
175  double testXDot;
176  double testDt;
177  double testTime;
178 };
179 
180 // ************************************************************
181 // ************************************************************
182 TEUCHOS_UNIT_TEST(NewmarkImplicitDForm, Default_Construction)
183 {
184  auto model = rcp(new Tempus_Test::HarmonicOscillatorModel<double>());
185 
186  // Default construction.
187  auto stepper = rcp(new Tempus::StepperNewmarkImplicitDForm<double>());
188  stepper->setModel(model);
189  stepper->initialize();
190  TEUCHOS_TEST_FOR_EXCEPT(!stepper->isInitialized());
191 
192  auto modifier = rcp(new Tempus::StepperNewmarkImplicitDFormModifierDefault<double>());
193  auto modifierX = rcp(new Tempus::StepperNewmarkImplicitDFormModifierXDefault<double>());
194 
195 
196  // Default values for construction.
197  auto solver = rcp(new Thyra::NOXNonlinearSolver());
198  solver->setParameterList(Tempus::defaultSolverParameters());
199 
200  bool useFSAL = stepper->getUseFSAL();
201  std::string ICConsistency = stepper->getICConsistency();
202  bool ICConsistencyCheck = stepper->getICConsistencyCheck();
203  bool zeroInitialGuess = stepper->getZeroInitialGuess();
204  std::string schemeName = "Average Acceleration";
205  double beta = 0.25;
206  double gamma = 0.5;
207 
208 
209  // Test the set functions.
210  stepper->setAppAction(modifier); stepper->initialize(); TEUCHOS_TEST_FOR_EXCEPT(!stepper->isInitialized());
211  stepper->setAppAction(modifierX); stepper->initialize(); TEUCHOS_TEST_FOR_EXCEPT(!stepper->isInitialized());
212  stepper->setSolver(solver); stepper->initialize(); TEUCHOS_TEST_FOR_EXCEPT(!stepper->isInitialized());
213  stepper->setUseFSAL(useFSAL); stepper->initialize(); TEUCHOS_TEST_FOR_EXCEPT(!stepper->isInitialized());
214  stepper->setICConsistency(ICConsistency); stepper->initialize(); TEUCHOS_TEST_FOR_EXCEPT(!stepper->isInitialized());
215  stepper->setICConsistencyCheck(ICConsistencyCheck); stepper->initialize(); TEUCHOS_TEST_FOR_EXCEPT(!stepper->isInitialized());
216  stepper->setZeroInitialGuess(zeroInitialGuess); stepper->initialize(); TEUCHOS_TEST_FOR_EXCEPT(!stepper->isInitialized());
217 
218  stepper->setSchemeName(schemeName); stepper->initialize(); TEUCHOS_TEST_FOR_EXCEPT(!stepper->isInitialized());
219  stepper->setBeta(beta); stepper->initialize(); TEUCHOS_TEST_FOR_EXCEPT(!stepper->isInitialized());
220  stepper->setGamma(gamma); stepper->initialize(); TEUCHOS_TEST_FOR_EXCEPT(!stepper->isInitialized());
221 
222 
223  // Full argument list construction.
224  stepper = rcp(new Tempus::StepperNewmarkImplicitDForm<double>(
225  model, solver, useFSAL,
226  ICConsistency, ICConsistencyCheck, zeroInitialGuess,
227  schemeName, beta, gamma, modifier));
228  TEUCHOS_TEST_FOR_EXCEPT(!stepper->isInitialized());
229 
230  // Test stepper properties.
231  TEUCHOS_ASSERT(stepper->getOrder() == 2);
232 }
233 
234 
235 // ************************************************************
236 // ************************************************************
237 TEUCHOS_UNIT_TEST(NewmarkImplicitDForm, StepperFactory_Construction)
238 {
239  auto model = rcp(new Tempus_Test::HarmonicOscillatorModel<double>());
240  testFactoryConstruction("Newmark Implicit d-Form", model);
241 }
242 
243 TEUCHOS_UNIT_TEST(NewmarkImplicitDForm, AppAction_Modifier)
244 {
245  using Teuchos::RCP;
246  using Teuchos::getParametersFromXmlFile;
247  using Teuchos::sublist;
248  using Teuchos::ParameterList;
249 
250  double dt = 1.0;
251 
252  // Read params from .xml file
253  RCP<ParameterList> pList = getParametersFromXmlFile(
254  "Tempus_NewmarkImplicitAForm_HarmonicOscillator_Damped_SecondOrder.xml");
255  RCP<ParameterList> pl = sublist(pList, "Tempus", true);
256 
257  // Setup the HarmonicOscillatorModel
258  RCP<ParameterList> hom_pl = sublist(pList, "HarmonicOscillatorModel", true);
259  RCP<const Thyra::ModelEvaluator<double> > model =
260  Teuchos::rcp(new Tempus_Test::HarmonicOscillatorModel<double>(hom_pl));
261 
262  // Setup Stepper for field solve ----------------------------
263  RCP<Tempus::StepperNewmarkImplicitDForm<double> > stepper =
264  Tempus::createStepperNewmarkImplicitDForm(model, Teuchos::null);
265 
266  auto modifier = rcp(new StepperNewmarkImplicitDFormModifierTest());
267  stepper->setAppAction(modifier);
268  stepper->initialize();
269 
270  // Setup TimeStepControl ------------------------------------
271  RCP<Tempus::TimeStepControl<double> > timeStepControl =
272  Teuchos::rcp(new Tempus::TimeStepControl<double>());
273  ParameterList tscPL = pl->sublist("Default Integrator")
274  .sublist("Time Step Control");
275  timeStepControl->setInitIndex(tscPL.get<int> ("Initial Time Index"));
276  timeStepControl->setInitTime (tscPL.get<double>("Initial Time"));
277  timeStepControl->setFinalTime(dt);
278  timeStepControl->setInitTimeStep(dt);
279  timeStepControl->initialize();
280 
281  // Setup initial condition SolutionState --------------------
282  using Teuchos::rcp_const_cast;
283  auto inArgsIC = model->getNominalValues();
284  RCP<Thyra::VectorBase<double> > icX =
285  rcp_const_cast<Thyra::VectorBase<double> > (inArgsIC.get_x());
286  RCP<Thyra::VectorBase<double> > icXDot =
287  rcp_const_cast<Thyra::VectorBase<double> > (inArgsIC.get_x_dot());
288  RCP<Thyra::VectorBase<double> > icXDotDot =
289  rcp_const_cast<Thyra::VectorBase<double> > (inArgsIC.get_x_dot_dot());
290  RCP<Tempus::SolutionState<double> > icState =
291  Tempus::createSolutionStateX(icX, icXDot, icXDotDot);
292  icState->setTime (timeStepControl->getInitTime());
293  icState->setIndex (timeStepControl->getInitIndex());
294  icState->setTimeStep(0.0);
295  icState->setOrder (stepper->getOrder());
296  icState->setSolutionStatus(Tempus::Status::PASSED); // ICs are passing.
297 
298  // Setup SolutionHistory ------------------------------------
299  RCP<Tempus::SolutionHistory<double> > solutionHistory =
300  Teuchos::rcp(new Tempus::SolutionHistory<double>());
301  solutionHistory->setName("Forward States");
302  solutionHistory->setStorageType(Tempus::STORAGE_TYPE_STATIC);
303  solutionHistory->setStorageLimit(2);
304  solutionHistory->addState(icState);
305 
306  // Setup Integrator -----------------------------------------
307  RCP<Tempus::IntegratorBasic<double> > integrator =
308  Tempus::createIntegratorBasic<double>();
309  integrator->setStepper(stepper);
310  integrator->setTimeStepControl(timeStepControl);
311  integrator->setSolutionHistory(solutionHistory);
312  integrator->initialize();
313 
314  // Integrate to timeMax
315  bool integratorStatus = integrator->advanceTime();
316  TEST_ASSERT(integratorStatus)
317 
318 
319  // Testing that each ACTION_LOCATION has been called.
320  TEST_COMPARE(modifier->testBEGIN_STEP, ==, true);
321  TEST_COMPARE(modifier->testBEFORE_SOLVE, ==, true);
322  TEST_COMPARE(modifier->testAFTER_SOLVE, ==, true);
323  TEST_COMPARE(modifier->testEND_STEP, ==, true);
324 
325  // Testing that values can be set through the Modifier.
326  auto x = integrator->getX();
327  auto Dt = integrator->getTime();
328  TEST_FLOATING_EQUALITY(modifier->testDt, Dt, 1.0e-14);
329  TEST_FLOATING_EQUALITY(modifier->testCurrentValue, get_ele(*(x), 0), 1.0e-14);
330  TEST_COMPARE(modifier->testName, ==, stepper->getStepperName());
331 }
332 
333 TEUCHOS_UNIT_TEST(NewmarkImplicitDForm, AppAction_ModifierX)
334 {
335  using Teuchos::RCP;
336  using Teuchos::getParametersFromXmlFile;
337  using Teuchos::sublist;
338  using Teuchos::ParameterList;
339 
340  double dt = 1.0;
341 
342  // Read params from .xml file
343  RCP<ParameterList> pList = getParametersFromXmlFile(
344  "Tempus_NewmarkImplicitAForm_HarmonicOscillator_Damped_SecondOrder.xml");
345  RCP<ParameterList> pl = sublist(pList, "Tempus", true);
346 
347  // Setup the HarmonicOscillatorModel
348  RCP<ParameterList> hom_pl = sublist(pList, "HarmonicOscillatorModel", true);
349  RCP<const Thyra::ModelEvaluator<double> > model =
350  Teuchos::rcp(new Tempus_Test::HarmonicOscillatorModel<double>(hom_pl));
351 
352  // Setup Stepper for field solve ----------------------------
353  RCP<Tempus::StepperNewmarkImplicitDForm<double> > stepper =
354  Tempus::createStepperNewmarkImplicitDForm(model, Teuchos::null);
355 
356  auto modifierX = rcp(new StepperNewmarkImplicitDFormModifierXTest());
357  stepper->setAppAction(modifierX);
358  stepper->initialize();
359 
360  // Setup TimeStepControl ------------------------------------
361  RCP<Tempus::TimeStepControl<double> > timeStepControl =
362  Teuchos::rcp(new Tempus::TimeStepControl<double>());
363  ParameterList tscPL = pl->sublist("Default Integrator")
364  .sublist("Time Step Control");
365  timeStepControl->setInitIndex(tscPL.get<int> ("Initial Time Index"));
366  timeStepControl->setInitTime (tscPL.get<double>("Initial Time"));
367  timeStepControl->setFinalTime(dt);
368  timeStepControl->setInitTimeStep(dt);
369  timeStepControl->initialize();
370 
371  // Setup initial condition SolutionState --------------------
372  using Teuchos::rcp_const_cast;
373  auto inArgsIC = model->getNominalValues();
374  RCP<Thyra::VectorBase<double> > icX =
375  rcp_const_cast<Thyra::VectorBase<double> > (inArgsIC.get_x());
376  RCP<Thyra::VectorBase<double> > icXDot =
377  rcp_const_cast<Thyra::VectorBase<double> > (inArgsIC.get_x_dot());
378  RCP<Thyra::VectorBase<double> > icXDotDot =
379  rcp_const_cast<Thyra::VectorBase<double> > (inArgsIC.get_x_dot_dot());
380  RCP<Tempus::SolutionState<double> > icState =
381  Tempus::createSolutionStateX(icX, icXDot, icXDotDot);
382  icState->setTime (timeStepControl->getInitTime());
383  icState->setIndex (timeStepControl->getInitIndex());
384  icState->setTimeStep(0.0);
385  icState->setOrder (stepper->getOrder());
386  icState->setSolutionStatus(Tempus::Status::PASSED); // ICs are passing.
387 
388  // Setup SolutionHistory ------------------------------------
389  RCP<Tempus::SolutionHistory<double> > solutionHistory =
390  Teuchos::rcp(new Tempus::SolutionHistory<double>());
391  solutionHistory->setName("Forward States");
392  solutionHistory->setStorageType(Tempus::STORAGE_TYPE_STATIC);
393  solutionHistory->setStorageLimit(2);
394  solutionHistory->addState(icState);
395 
396  // Setup Integrator -----------------------------------------
397  RCP<Tempus::IntegratorBasic<double> > integrator =
398  Tempus::createIntegratorBasic<double>();
399  integrator->setStepper(stepper);
400  integrator->setTimeStepControl(timeStepControl);
401  integrator->setSolutionHistory(solutionHistory);
402  integrator->initialize();
403 
404  // Integrate to timeMax
405  bool integratorStatus = integrator->advanceTime();
406  TEST_ASSERT(integratorStatus)
407 
408  // Testing that each ACTION_LOCATION has been called.
409  TEST_COMPARE(modifierX->testX_BEGIN_STEP, ==, true);
410  TEST_COMPARE(modifierX->testX_BEFORE_SOLVE, ==, true);
411  TEST_COMPARE(modifierX->testX_AFTER_SOLVE, ==, true);
412  TEST_COMPARE(modifierX->testX_END_STEP, ==, true);
413 
414 
415  // Testing that values can be set through the Modifier.
416  auto Dt = integrator->getTime();
417  TEST_FLOATING_EQUALITY(modifierX->testDt, Dt, 1.0e-14);
418 
419  const auto x = integrator->getX();
420  TEST_FLOATING_EQUALITY(modifierX->testX, get_ele(*(x), 0), 1.0e-14);
421 }
422 
423 } // namespace Tempus_Test
Tempus::createSolutionStateX
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.
Definition: Tempus_SolutionState_impl.hpp:479
Tempus_Unit_Test::testFactoryConstruction
void testFactoryConstruction(std::string stepperType, const Teuchos::RCP< const Thyra::ModelEvaluator< double > > &model)
Unit test utility for Stepper construction through StepperFactory.
Definition: Tempus_UnitTest_Utils.hpp:45
Tempus::createStepperNewmarkImplicitDForm
Teuchos::RCP< StepperNewmarkImplicitDForm< Scalar > > createStepperNewmarkImplicitDForm(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &model, Teuchos::RCP< Teuchos::ParameterList > pl)
Nonmember constructor - ModelEvaluator and ParameterList.
Definition: Tempus_StepperNewmarkImplicitDForm_impl.hpp:518
Teuchos::ParameterList::get
T & get(const std::string &name, T def_value)
TEUCHOS_TEST_FOR_EXCEPTION
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
TEST_COMPARE
#define TEST_COMPARE(v1, comp, v2)
Tempus::StepperNewmarkImplicitDFormModifierXBase< double >::modify
virtual void modify(Teuchos::RCP< Thyra::VectorBase< double > >, const double, const double, const MODIFIER_TYPE modType)=0
Modify solution based on the MODIFIER_TYPE.
TEST_FLOATING_EQUALITY
#define TEST_FLOATING_EQUALITY(v1, v2, tol)
TEST_ASSERT
#define TEST_ASSERT(v1)
Tempus::defaultSolverParameters
Teuchos::RCP< Teuchos::ParameterList > defaultSolverParameters()
Returns the default solver ParameterList for implicit Steppers.
Definition: Tempus_Stepper_impl.hpp:414
Tempus_StepperNewmarkImplicitDFormModifierDefault.hpp
Tempus_StepperNewmarkImplicitDFormModifierXBase.hpp
Tempus_Test::HarmonicOscillatorModel
Consider the ODE: where is a constant, is a constant damping parameter, is a constant forcing par...
Definition: HarmonicOscillatorModel_decl.hpp:52
Tempus::StepperNewmarkImplicitDFormAppAction::ACTION_LOCATION
ACTION_LOCATION
Indicates the location of application action (see algorithm).
Definition: Tempus_StepperNewmarkImplicitDFormAppAction.hpp:39
Teuchos::rcp
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
Tempus_StepperNewmarkImplicitDFormModifierBase.hpp
Teuchos_UnitTestHarness.hpp
Tempus_UnitTest_Utils.hpp
Tempus_StepperNewmarkImplicitDFormModifierXDefault.hpp
Tempus_Unit_Test::TEUCHOS_UNIT_TEST
TEUCHOS_UNIT_TEST(BackwardEuler, Default_Construction)
Definition: Tempus_UnitTest_BackwardEuler.cpp:49
Tempus::TimeStepControl
TimeStepControl manages the time step size. There several mechanisms that effect the time step size a...
Definition: Tempus_Integrator.hpp:26
Tempus::SolutionHistory
SolutionHistory is basically a container of SolutionStates. SolutionHistory maintains a collection of...
Definition: Tempus_Integrator.hpp:25
Teuchos_XMLParameterListHelpers.hpp
Tempus::STORAGE_TYPE_STATIC
Keep a fix number of states.
Definition: Tempus_SolutionHistory_decl.hpp:23
Teuchos::ParameterList::sublist
ParameterList & sublist(const std::string &name, bool mustAlreadyExist=false, const std::string &docString="")
Teuchos_TimeMonitor.hpp
TEUCHOS_ASSERT
#define TEUCHOS_ASSERT(assertion_test)
Tempus::StepperNewmarkImplicitDFormModifierXBase::MODIFIER_TYPE
MODIFIER_TYPE
Indicates the location of application action (see algorithm).
Definition: Tempus_StepperNewmarkImplicitDFormModifierXBase.hpp:104
TEUCHOS_TEST_FOR_EXCEPT
#define TEUCHOS_TEST_FOR_EXCEPT(throw_exception_test)
Generated by   doxygen 1.8.5