Tempus  Version of the Day
Time Integration
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
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"
18 #include "Tempus_TimeStepControlStrategyIntegralController.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 {
42  auto tscs = rcp(new Tempus::TimeStepControlStrategyIntegralController<double>());
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 {
82  auto tscs = rcp(new Tempus::TimeStepControlStrategyIntegralController<double>(
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 =
140  rcp(new Tempus::TimeStepControlStrategyIntegralController<double>(
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());
155  Tempus::Status status = Tempus::Status::WORKING;
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  // Mock Integrator
172 
173  // -- First Time Step
174  solutionHistory->initWorkingState();
175  auto currentState = solutionHistory->getCurrentState();
176  auto workingState = solutionHistory->getWorkingState();
177 
178  tsc->setNextTimeStep(solutionHistory, status);
179 
180  // First time step should cause no change to dt because
181  // internal relative errors = 1.
182  TEST_FLOATING_EQUALITY(workingState->getTimeStep(), 1.0, 1.0e-14);
183 
184  // Mock takeStep
185  double errN = 0.1;
186  workingState->setErrorRel(errN);
187  workingState->setSolutionStatus(Tempus::Status::PASSED);
188 
189 
190 
191  // -- Second Time Step
192  solutionHistory->initWorkingState();
193  currentState = solutionHistory->getCurrentState();
194  workingState = solutionHistory->getWorkingState();
195  double dt = workingState->getTimeStep();
196 
197  tsc->setNextTimeStep(solutionHistory, status);
198 
199  double p = order - 1.0;
200  double dtNew = dt*safetyFactor*std::pow(errN, -KI/p);
201  TEST_FLOATING_EQUALITY(workingState->getTimeStep(), dtNew, 1.0e-14);
202 
203  // Mock takeStep
204  errN = 0.2;
205  double errNm1 = 0.1;
206  workingState->setErrorRel(errN);
207  workingState->setSolutionStatus(Tempus::Status::PASSED);
208 
209 
210 
211  // -- Third Time Step
212  solutionHistory->initWorkingState();
213  currentState = solutionHistory->getCurrentState();
214  workingState = solutionHistory->getWorkingState();
215  dt = workingState->getTimeStep();
216 
217  tsc->setNextTimeStep(solutionHistory, status);
218 
219  dtNew = dt*safetyFactor*std::pow(errN, -KI/p)
220  *std::pow(errNm1, KP/p);
221  TEST_FLOATING_EQUALITY(workingState->getTimeStep(), dtNew, 1.0e-14);
222 
223  // Mock takeStep
224  errN = 0.3;
225  errNm1 = 0.2;
226  double errNm2 = 0.1;
227  workingState->setErrorRel(errN);
228  workingState->setSolutionStatus(Tempus::Status::PASSED);
229 
230 
231 
232  // -- Fourth Time Step
233  solutionHistory->initWorkingState();
234  currentState = solutionHistory->getCurrentState();
235  workingState = solutionHistory->getWorkingState();
236  dt = workingState->getTimeStep();
237 
238  tsc->setNextTimeStep(solutionHistory, status);
239 
240  dtNew = dt*safetyFactor*std::pow(errN, -KI/p)
241  *std::pow(errNm1, KP/p)
242  *std::pow(errNm2, -KD/p);
243  TEST_FLOATING_EQUALITY(workingState->getTimeStep(), dtNew, 1.0e-14);
244 
245 }
246 
247 
248 // ************************************************************
249 // ************************************************************
250 TEUCHOS_UNIT_TEST(TimeStepControlStrategyIntegralController, getValidParameters)
251 {
252  auto tscs = rcp(new Tempus::TimeStepControlStrategyIntegralController<double>());
253 
254  auto pl = tscs->getValidParameters();
255 
256  TEST_COMPARE ( pl->get<std::string>("Strategy Type"), ==,"Integral Controller");
257  TEST_COMPARE ( pl->get<std::string>("Controller Type"), ==,"PID");
258  TEST_FLOATING_EQUALITY( pl->get<double>("KI"), 0.58, 1.0e-14);
259  TEST_FLOATING_EQUALITY( pl->get<double>("KP"), 0.21, 1.0e-14);
260  TEST_FLOATING_EQUALITY( pl->get<double>("KD"), 0.10, 1.0e-14);
261  TEST_FLOATING_EQUALITY( pl->get<double>("Safety Factor"), 0.9, 1.0e-14);
262  TEST_FLOATING_EQUALITY( pl->get<double>("Safety Factor After Step Rejection"), 0.9, 1.0e-14);
263  TEST_FLOATING_EQUALITY( pl->get<double>("Maximum Safety Factor"), 5.0, 1.0e-14);
264  TEST_FLOATING_EQUALITY( pl->get<double>("Minimum Safety Factor"), 0.5, 1.0e-14);
265 
266  { // Ensure that parameters are "used", excluding sublists.
267  std::ostringstream unusedParameters;
268  pl->unused(unusedParameters);
269  TEST_COMPARE ( unusedParameters.str(), ==, "");
270  }
271 }
272 
273 
274 } // namespace Tempus_Test
Tempus_Test::DahlquistTestModel
The classic Dahlquist Test Problem.
Definition: DahlquistTestModel_decl.hpp:35
TEST_COMPARE
#define TEST_COMPARE(v1, comp, v2)
TEST_FLOATING_EQUALITY
#define TEST_FLOATING_EQUALITY(v1, v2, tol)
Teuchos::rcp
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
Tempus::Status
Status
Status for the Integrator, the Stepper and the SolutionState.
Definition: Tempus_Types.hpp:19
Teuchos_UnitTestHarness.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
Teuchos_TimeMonitor.hpp
Tempus_TimeStepControlStrategyIntegralController.hpp
TEUCHOS_TEST_FOR_EXCEPT
#define TEUCHOS_TEST_FOR_EXCEPT(throw_exception_test)
Generated by   doxygen 1.8.5