Tempus  Version of the Day
Time Integration
Tempus_UnitTest_TimeStepControlStrategyIntegralController.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_StepperFactory.hpp"
17 #include "Tempus_TimeStepControl.hpp"
19 
20 #include "../TestModels/SinCosModel.hpp"
21 #include "../TestModels/DahlquistTestModel.hpp"
22 #include "../TestUtils/Tempus_ConvergenceTestUtils.hpp"
23 
24 #include <fstream>
25 #include <vector>
26 
27 namespace Tempus_Unit_Test {
28 
29 using Teuchos::RCP;
30 using Teuchos::rcp;
31 using Teuchos::rcp_const_cast;
32 using Teuchos::rcp_dynamic_cast;
33 using Teuchos::ParameterList;
34 using Teuchos::sublist;
35 using Teuchos::getParametersFromXmlFile;
36 
37 
38 // ************************************************************
39 // ************************************************************
40 TEUCHOS_UNIT_TEST(TimeStepControlStrategyIntegralController, Default_Construction)
41 {
43  TEUCHOS_TEST_FOR_EXCEPT(!tscs->isInitialized());
44 
45  // Test the get functions (i.e., defaults).
46  TEUCHOS_TEST_FOR_EXCEPT(tscs->getStepType() != "Variable");
47  TEUCHOS_TEST_FOR_EXCEPT(tscs->getController() != "PID");
48  TEUCHOS_TEST_FOR_EXCEPT(tscs->getKI() != 0.58);
49  TEUCHOS_TEST_FOR_EXCEPT(tscs->getKP() != 0.21);
50  TEUCHOS_TEST_FOR_EXCEPT(tscs->getKD() != 0.10);
51  TEUCHOS_TEST_FOR_EXCEPT(tscs->getSafetyFactor() != 0.90);
52  TEUCHOS_TEST_FOR_EXCEPT(tscs->getSafetyFactorAfterReject() != 0.90);
53  TEUCHOS_TEST_FOR_EXCEPT(tscs->getFacMax() != 5.0);
54  TEUCHOS_TEST_FOR_EXCEPT(tscs->getFacMin() != 0.5);
55 
56  // Test the set functions.
57  tscs->setController("I"); tscs->initialize(); TEUCHOS_TEST_FOR_EXCEPT(!tscs->isInitialized());
58  tscs->setKI(0.6); tscs->initialize(); TEUCHOS_TEST_FOR_EXCEPT(!tscs->isInitialized());
59  tscs->setKP(0.0); tscs->initialize(); TEUCHOS_TEST_FOR_EXCEPT(!tscs->isInitialized());
60  tscs->setKD(0.0); tscs->initialize(); TEUCHOS_TEST_FOR_EXCEPT(!tscs->isInitialized());
61  tscs->setSafetyFactor(0.8); tscs->initialize(); TEUCHOS_TEST_FOR_EXCEPT(!tscs->isInitialized());
62  tscs->setSafetyFactorAfterReject(0.8); tscs->initialize(); TEUCHOS_TEST_FOR_EXCEPT(!tscs->isInitialized());
63  tscs->setFacMax(4.0); tscs->initialize(); TEUCHOS_TEST_FOR_EXCEPT(!tscs->isInitialized());
64  tscs->setFacMin(0.4); tscs->initialize(); TEUCHOS_TEST_FOR_EXCEPT(!tscs->isInitialized());
65 
66  TEUCHOS_TEST_FOR_EXCEPT(tscs->getStepType() != "Variable");
67  TEUCHOS_TEST_FOR_EXCEPT(tscs->getController() != "I");
68  TEUCHOS_TEST_FOR_EXCEPT(tscs->getKI() != 0.6);
69  TEUCHOS_TEST_FOR_EXCEPT(tscs->getKP() != 0.0);
70  TEUCHOS_TEST_FOR_EXCEPT(tscs->getKD() != 0.0);
71  TEUCHOS_TEST_FOR_EXCEPT(tscs->getSafetyFactor() != 0.8);
72  TEUCHOS_TEST_FOR_EXCEPT(tscs->getSafetyFactorAfterReject() != 0.8);
73  TEUCHOS_TEST_FOR_EXCEPT(tscs->getFacMax() != 4.0);
74  TEUCHOS_TEST_FOR_EXCEPT(tscs->getFacMin() != 0.4);
75 }
76 
77 
78 // ************************************************************
79 // ************************************************************
80 TEUCHOS_UNIT_TEST(TimeStepControlStrategyIntegralController, Full_Construction)
81 {
83  "I", 0.6, 0.0, 0.0, 0.8, 0.8, 4.0, 0.4));
84  TEUCHOS_TEST_FOR_EXCEPT(!tscs->isInitialized());
85 
86  TEUCHOS_TEST_FOR_EXCEPT(tscs->getStepType() != "Variable");
87  TEUCHOS_TEST_FOR_EXCEPT(tscs->getController() != "I");
88  TEUCHOS_TEST_FOR_EXCEPT(tscs->getKI() != 0.6);
89  TEUCHOS_TEST_FOR_EXCEPT(tscs->getKP() != 0.0);
90  TEUCHOS_TEST_FOR_EXCEPT(tscs->getKD() != 0.0);
91  TEUCHOS_TEST_FOR_EXCEPT(tscs->getSafetyFactor() != 0.8);
92  TEUCHOS_TEST_FOR_EXCEPT(tscs->getSafetyFactorAfterReject() != 0.8);
93  TEUCHOS_TEST_FOR_EXCEPT(tscs->getFacMax() != 4.0);
94  TEUCHOS_TEST_FOR_EXCEPT(tscs->getFacMin() != 0.4);
95 }
96 
97 
98 // ************************************************************
99 // ************************************************************
100 TEUCHOS_UNIT_TEST(TimeStepControlStrategyIntegralController, Create_Construction)
101 {
102  auto pl = Tempus::getTimeStepControlStrategyIntegralControllerPL<double>();
103 
104  pl->set<std::string>("Controller Type", "I");
105  pl->set<double>("KI", 0.6);
106  pl->set<double>("KP", 0.0);
107  pl->set<double>("KD", 0.0);
108  pl->set<double>("Safety Factor", 0.8);
109  pl->set<double>("Safety Factor After Step Rejection", 0.8);
110  pl->set<double>("Maximum Safety Factor", 4.0);
111  pl->set<double>("Minimum Safety Factor", 0.4);
112 
113  auto tscs = Tempus::createTimeStepControlStrategyIntegralController<double>(pl);
114 
115  TEUCHOS_TEST_FOR_EXCEPT(tscs->getStepType() != "Variable");
116  TEUCHOS_TEST_FOR_EXCEPT(tscs->getController() != "I");
117  TEUCHOS_TEST_FOR_EXCEPT(tscs->getKI() != 0.6);
118  TEUCHOS_TEST_FOR_EXCEPT(tscs->getKP() != 0.0);
119  TEUCHOS_TEST_FOR_EXCEPT(tscs->getKD() != 0.0);
120  TEUCHOS_TEST_FOR_EXCEPT(tscs->getSafetyFactor() != 0.8);
121  TEUCHOS_TEST_FOR_EXCEPT(tscs->getSafetyFactorAfterReject() != 0.8);
122  TEUCHOS_TEST_FOR_EXCEPT(tscs->getFacMax() != 4.0);
123  TEUCHOS_TEST_FOR_EXCEPT(tscs->getFacMin() != 0.4);
124 }
125 
126 
127 // ************************************************************
128 // ************************************************************
129 TEUCHOS_UNIT_TEST(TimeStepControlStrategyIntegralController, setNextTimeStep)
130 {
131  double KI = 0.5;
132  double KP = 0.25;
133  double KD = 0.15;
134  double safetyFactor = 0.9;
135  double safetyFactorAfterReject = 0.9;
136  double facMax = 5.0;
137  double facMin = 0.5;
138 
139  auto tscs =
141  "PID", KI, KP, KD, safetyFactor, safetyFactorAfterReject,
142  facMax, facMin));
143 
144  // Setup the TimeStepControl --------------------------------
145  auto tsc = rcp(new Tempus::TimeStepControl<double>());
146  tsc->setTimeStepControlStrategy(tscs);
147  tsc->setInitTime(0.0);
148  tsc->setFinalTime(10.0);
149  tsc->setMinTimeStep (0.01);
150  tsc->setInitTimeStep(1.0);
151  tsc->setMaxTimeStep (10.0);
152  tsc->setFinalIndex(100);
153  tsc->initialize();
154  TEUCHOS_TEST_FOR_EXCEPT(!tsc->isInitialized());
156 
157  // Setup the SolutionHistory --------------------------------
158  auto model = rcp(new Tempus_Test::DahlquistTestModel<double>());
159  auto inArgsIC = model->getNominalValues();
160  auto icSolution = rcp_const_cast<Thyra::VectorBase<double> >(inArgsIC.get_x());
161  auto icState = Tempus::createSolutionStateX<double>(icSolution);
162  auto solutionHistory = rcp(new Tempus::SolutionHistory<double>());
163 
164  double order = 2.0;
165  solutionHistory->addState(icState);
166  solutionHistory->getCurrentState()->setTimeStep(1.0);
167  solutionHistory->getCurrentState()->setTime(0.0);
168  solutionHistory->getCurrentState()->setIndex(0);
169  solutionHistory->getCurrentState()->setOrder(order);
170 
171 
172  // Mock Integrator
173 
174  // -- First Time Step
175  solutionHistory->initWorkingState();
176  auto currentState = solutionHistory->getCurrentState();
177  auto workingState = solutionHistory->getWorkingState();
178 
179  TEST_FLOATING_EQUALITY(workingState->getErrorRel() , 0.0, 1.0e-14);
180  TEST_FLOATING_EQUALITY(workingState->getErrorRelNm1(), 0.0, 1.0e-14);
181  TEST_FLOATING_EQUALITY(workingState->getErrorRelNm2(), 0.0, 1.0e-14);
182 
183  tsc->setNextTimeStep(solutionHistory, status);
184 
185  // First time step should cause no change to dt because
186  // internal relative errors = 1.
187  TEST_FLOATING_EQUALITY(workingState->getTimeStep(), 1.0, 1.0e-14);
188 
189  // Mock takeStep
190  double errN = 0.1;
191  workingState->setErrorRel(errN);
192  workingState->setSolutionStatus(Tempus::Status::PASSED);
193 
194 
195 
196  // -- Second Time Step
197  solutionHistory->initWorkingState();
198  currentState = solutionHistory->getCurrentState();
199  workingState = solutionHistory->getWorkingState();
200  double dt = workingState->getTimeStep();
201 
202  TEST_FLOATING_EQUALITY(workingState->getErrorRel() , 0.1, 1.0e-14);
203  TEST_FLOATING_EQUALITY(workingState->getErrorRelNm1(), 0.0, 1.0e-14);
204  TEST_FLOATING_EQUALITY(workingState->getErrorRelNm2(), 0.0, 1.0e-14);
205 
206  tsc->setNextTimeStep(solutionHistory, status);
207 
208  double p = order - 1.0;
209  double dtNew = dt*safetyFactor*std::pow(errN, -KI/p);
210  TEST_FLOATING_EQUALITY(workingState->getTimeStep(), dtNew, 1.0e-14);
211 
212  // Mock takeStep
213  errN = 0.2;
214  double errNm1 = 0.1;
215  workingState->setErrorRel(errN);
216  workingState->setSolutionStatus(Tempus::Status::PASSED);
217 
218 
219 
220  // -- Third Time Step
221  solutionHistory->initWorkingState();
222  currentState = solutionHistory->getCurrentState();
223  workingState = solutionHistory->getWorkingState();
224  dt = workingState->getTimeStep();
225 
226  TEST_FLOATING_EQUALITY(workingState->getErrorRel() , 0.2, 1.0e-14);
227  TEST_FLOATING_EQUALITY(workingState->getErrorRelNm1(), 0.1, 1.0e-14);
228  TEST_FLOATING_EQUALITY(workingState->getErrorRelNm2(), 0.0, 1.0e-14);
229 
230  tsc->setNextTimeStep(solutionHistory, status);
231 
232  dtNew = dt*safetyFactor*std::pow(errN, -KI/p)
233  *std::pow(errNm1, KP/p);
234  TEST_FLOATING_EQUALITY(workingState->getTimeStep(), dtNew, 1.0e-14);
235 
236  // Mock takeStep
237  errN = 0.3;
238  errNm1 = 0.2;
239  double errNm2 = 0.1;
240  workingState->setErrorRel(errN);
241  workingState->setSolutionStatus(Tempus::Status::PASSED);
242 
243 
244 
245  // -- Fourth Time Step
246  solutionHistory->initWorkingState();
247  currentState = solutionHistory->getCurrentState();
248  workingState = solutionHistory->getWorkingState();
249  dt = workingState->getTimeStep();
250 
251  TEST_FLOATING_EQUALITY(workingState->getErrorRel() , 0.3, 1.0e-14);
252  TEST_FLOATING_EQUALITY(workingState->getErrorRelNm1(), 0.2, 1.0e-14);
253  TEST_FLOATING_EQUALITY(workingState->getErrorRelNm2(), 0.1, 1.0e-14);
254 
255  tsc->setNextTimeStep(solutionHistory, status);
256 
257  dtNew = dt*safetyFactor*std::pow(errN, -KI/p)
258  *std::pow(errNm1, KP/p)
259  *std::pow(errNm2, -KD/p);
260  TEST_FLOATING_EQUALITY(workingState->getTimeStep(), dtNew, 1.0e-14);
261 
262 }
263 
264 
265 // ************************************************************
266 // ************************************************************
267 TEUCHOS_UNIT_TEST(TimeStepControlStrategyIntegralController, getValidParameters)
268 {
270 
271  auto pl = tscs->getValidParameters();
272 
273  TEST_COMPARE ( pl->get<std::string>("Strategy Type"), ==,"Integral Controller");
274  TEST_COMPARE ( pl->get<std::string>("Controller Type"), ==,"PID");
275  TEST_FLOATING_EQUALITY( pl->get<double>("KI"), 0.58, 1.0e-14);
276  TEST_FLOATING_EQUALITY( pl->get<double>("KP"), 0.21, 1.0e-14);
277  TEST_FLOATING_EQUALITY( pl->get<double>("KD"), 0.10, 1.0e-14);
278  TEST_FLOATING_EQUALITY( pl->get<double>("Safety Factor"), 0.9, 1.0e-14);
279  TEST_FLOATING_EQUALITY( pl->get<double>("Safety Factor After Step Rejection"), 0.9, 1.0e-14);
280  TEST_FLOATING_EQUALITY( pl->get<double>("Maximum Safety Factor"), 5.0, 1.0e-14);
281  TEST_FLOATING_EQUALITY( pl->get<double>("Minimum Safety Factor"), 0.5, 1.0e-14);
282 
283  { // Ensure that parameters are "used", excluding sublists.
284  std::ostringstream unusedParameters;
285  pl->unused(unusedParameters);
286  TEST_COMPARE ( unusedParameters.str(), ==, "");
287  }
288 }
289 
290 
291 } // namespace Tempus_Test
The classic Dahlquist Test Problem.
Status
Status for the Integrator, the Stepper and the SolutionState.
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...