54 RCP<ParameterList> pList =
55 getParametersFromXmlFile(
"Tempus_Subcycling_SinCos.xml");
62 RCP<ParameterList> scm_pl = sublist(pList,
"SinCosModel",
true);
63 auto model = rcp(
new SinCosModel<double> (scm_pl));
65 RCP<ParameterList> tempusPL = sublist(pList,
"Tempus",
true);
69 RCP<Tempus::IntegratorBasic<double> > integrator =
70 Tempus::createIntegratorBasic<double>(tempusPL, model);
72 RCP<ParameterList> stepperPL = sublist(tempusPL,
"Demo Stepper",
true);
73 RCP<const ParameterList> defaultPL =
74 integrator->getStepper()->getValidParameters();
76 bool pass = haveSameValuesSorted(*stepperPL, *defaultPL,
true);
79 out <<
"stepperPL -------------- \n" << *stepperPL << std::endl;
80 out <<
"defaultPL -------------- \n" << *defaultPL << std::endl;
87 RCP<Tempus::IntegratorBasic<double> > integrator =
88 Tempus::createIntegratorBasic<double>(model, std::string(
"Forward Euler"));
90 RCP<ParameterList> stepperPL = sublist(tempusPL,
"Demo Stepper",
true);
91 RCP<const ParameterList> defaultPL =
92 integrator->getStepper()->getValidParameters();
94 bool pass = haveSameValuesSorted(*stepperPL, *defaultPL,
true);
96 std::cout << std::endl;
97 std::cout <<
"stepperPL -------------- \n" << *stepperPL << std::endl;
98 std::cout <<
"defaultPL -------------- \n" << *defaultPL << std::endl;
112 auto model = rcp(
new SinCosModel<double>());
113 auto modelME = rcp_dynamic_cast<const Thyra::ModelEvaluator<double> > (model);
118 stepper->setSubcyclingStepper(stepperFE);
120 stepper->setSubcyclingMinTimeStep (0.1);
121 stepper->setSubcyclingInitTimeStep (0.1);
122 stepper->setSubcyclingMaxTimeStep (0.1);
123 stepper->setSubcyclingMaxFailures (10);
124 stepper->setSubcyclingMaxConsecFailures(5);
125 stepper->setSubcyclingScreenOutputIndexInterval(1);
126 stepper->setSubcyclingIntegratorObserver(
128 stepper->setSubcyclingPrintDtChanges (
true);
132 stepper->setSubcyclingTimeStepControlStrategy(subStrategy);
136 timeStepControl->setInitIndex(0);
137 timeStepControl->setFinalIndex(10);
138 timeStepControl->setInitTime (0.0);
139 timeStepControl->setFinalTime(1.0);
140 timeStepControl->setInitTimeStep(dt);
144 strategy->initialize();
145 timeStepControl->setTimeStepControlStrategy(strategy);
147 timeStepControl->initialize();
150 auto inArgsIC = stepper->getModel()->getNominalValues();
151 auto icSolution = rcp_const_cast<Thyra::VectorBase<double> > (inArgsIC.get_x());
153 icState->setTime (timeStepControl->getInitTime());
154 icState->setIndex (timeStepControl->getInitIndex());
155 icState->setTimeStep(0.0);
160 solutionHistory->setName(
"Forward States");
162 solutionHistory->setStorageLimit(2);
163 solutionHistory->addState(icState);
166 stepper->setInitialConditions(solutionHistory);
167 stepper->initialize();
170 RCP<Tempus::IntegratorBasic<double> > integrator =
171 Tempus::createIntegratorBasic<double>();
172 integrator->setStepper(stepper);
173 integrator->setTimeStepControl(timeStepControl);
174 integrator->setSolutionHistory(solutionHistory);
175 integrator->setScreenOutputIndexInterval(1);
177 integrator->initialize();
181 bool integratorStatus = integrator->advanceTime();
182 TEST_ASSERT(integratorStatus)
186 double time = integrator->getTime();
187 double timeFinal = 1.0;
188 TEST_FLOATING_EQUALITY(time, timeFinal, 1.0e-14);
191 RCP<Thyra::VectorBase<double> > x = integrator->getX();
192 RCP<const Thyra::VectorBase<double> > x_exact =
193 model->getExactSolution(time).get_x();
196 RCP<Thyra::VectorBase<double> > xdiff = x->clone_v();
197 Thyra::V_StVpStV(xdiff.ptr(), 1.0, *x_exact, -1.0, *(x));
200 std::cout <<
" Stepper = " << stepper->description() << std::endl;
201 std::cout <<
" =========================" << std::endl;
202 std::cout <<
" Exact solution : " << get_ele(*(x_exact), 0) <<
" "
203 << get_ele(*(x_exact), 1) << std::endl;
204 std::cout <<
" Computed solution: " << get_ele(*(x ), 0) <<
" "
205 << get_ele(*(x ), 1) << std::endl;
206 std::cout <<
" Difference : " << get_ele(*(xdiff ), 0) <<
" "
207 << get_ele(*(xdiff ), 1) << std::endl;
208 std::cout <<
" =========================" << std::endl;
209 TEST_FLOATING_EQUALITY(get_ele(*(x), 0), 0.882508, 1.0e-4 );
210 TEST_FLOATING_EQUALITY(get_ele(*(x), 1), 0.570790, 1.0e-4 );
218 RCP<Tempus::IntegratorBasic<double> > integrator;
219 std::vector<RCP<Thyra::VectorBase<double>>> solutions;
220 std::vector<RCP<Thyra::VectorBase<double>>> solutionsDot;
221 std::vector<double> StepSize;
226 const int nTimeStepSizes = 2;
227 std::string output_file_string =
"Tempus_Subcycling_SinCos";
228 std::string output_file_name = output_file_string +
".dat";
229 std::string err_out_file_name = output_file_string +
"-Error.dat";
231 for (
int n=0; n<nTimeStepSizes; n++) {
236 auto model = rcp(
new SinCosModel<double>());
237 auto modelME=rcp_dynamic_cast<const Thyra::ModelEvaluator<double> > (model);
242 stepper->setSubcyclingStepper(stepperFE);
244 stepper->setSubcyclingMinTimeStep (dt/10.0);
245 stepper->setSubcyclingInitTimeStep (dt/10.0);
246 stepper->setSubcyclingMaxTimeStep (dt);
247 stepper->setSubcyclingMaxFailures (10);
248 stepper->setSubcyclingMaxConsecFailures(5);
249 stepper->setSubcyclingScreenOutputIndexInterval(1);
256 strategy->setMinEta(0.02);
257 strategy->setMaxEta(0.04);
258 strategy->initialize();
259 stepper->setSubcyclingTimeStepControlStrategy(strategy);
263 timeStepControl->setInitIndex(0);
264 timeStepControl->setInitTime (0.0);
265 timeStepControl->setFinalTime(1.0);
266 timeStepControl->setMinTimeStep (dt);
267 timeStepControl->setInitTimeStep(dt);
268 timeStepControl->setMaxTimeStep (dt);
269 timeStepControl->initialize();
272 auto inArgsIC = stepper->getModel()->getNominalValues();
273 auto icSolution = rcp_const_cast<Thyra::VectorBase<double> > (inArgsIC.get_x());
275 icState->setTime (timeStepControl->getInitTime());
276 icState->setIndex (timeStepControl->getInitIndex());
277 icState->setTimeStep(0.0);
282 solutionHistory->setName(
"Forward States");
284 solutionHistory->setStorageLimit(2);
285 solutionHistory->addState(icState);
288 stepper->setInitialConditions(solutionHistory);
289 stepper->initialize();
292 integrator = Tempus::createIntegratorBasic<double>();
293 integrator->setStepper(stepper);
294 integrator->setTimeStepControl(timeStepControl);
295 integrator->setSolutionHistory(solutionHistory);
296 integrator->setScreenOutputIndexInterval(10);
298 integrator->initialize();
302 bool integratorStatus = integrator->advanceTime();
303 TEST_ASSERT(integratorStatus)
306 time = integrator->getTime();
307 double timeFinal = 1.0;
308 TEST_FLOATING_EQUALITY(time, timeFinal, 1.0e-14);
311 RCP<Thyra::VectorBase<double> > x = integrator->getX();
312 RCP<const Thyra::VectorBase<double> > x_exact =
313 model->getExactSolution(time).get_x();
344 StepSize.push_back(dt);
345 auto solution = Thyra::createMember(model->get_x_space());
346 Thyra::copy(*(integrator->getX()),solution.ptr());
347 solutions.push_back(solution);
348 if (n == nTimeStepSizes-1) {
349 StepSize.push_back(0.0);
350 auto solutionExact = Thyra::createMember(model->get_x_space());
351 Thyra::copy(*(model->getExactSolution(time).get_x()),solutionExact.ptr());
352 solutions.push_back(solutionExact);
357 if (nTimeStepSizes > 1) {
359 double xDotSlope = 0.0;
360 std::vector<double> xErrorNorm;
361 std::vector<double> xDotErrorNorm;
362 RCP<Tempus::Stepper<double> > stepper = integrator->getStepper();
366 solutions, xErrorNorm, xSlope,
367 solutionsDot, xDotErrorNorm, xDotSlope);
369 TEST_FLOATING_EQUALITY( xSlope, 1.00137, 0.01 );
371 TEST_FLOATING_EQUALITY( xErrorNorm[0], 0.00387948, 1.0e-4 );
375 Teuchos::TimeMonitor::summarize();
383 RCP<Tempus::IntegratorBasic<double> > integrator;
384 std::vector<RCP<Thyra::VectorBase<double>>> solutions;
385 std::vector<RCP<Thyra::VectorBase<double>>> solutionsDot;
386 std::vector<double> StepSize;
387 std::vector<double> xErrorNorm;
388 std::vector<double> xDotErrorNorm;
389 const int nTimeStepSizes = 4;
392 for (
int n=0; n<nTimeStepSizes; n++) {
396 if (n == nTimeStepSizes-1) dt /= 10.0;
399 auto tmpModel = rcp(
new VanDerPol_IMEX_ExplicitModel<double>());
400 auto pl = Teuchos::rcp_const_cast<Teuchos::ParameterList> (
401 tmpModel->getValidParameters());
402 pl->set(
"Coeff epsilon", 0.1);
403 RCP<const Thyra::ModelEvaluator<double> > explicitModel =
404 rcp(
new VanDerPol_IMEX_ExplicitModel<double>(pl));
405 RCP<const Thyra::ModelEvaluator<double> > implicitModel =
406 rcp(
new VanDerPol_IMEX_ImplicitModel<double>(pl));
413 stepperFE->setUseFSAL(
false);
414 stepperFE->initialize();
415 stepperSC->setSubcyclingStepper(stepperFE);
417 stepperSC->setSubcyclingMinTimeStep (0.00001);
418 stepperSC->setSubcyclingInitTimeStep (dt/10.0);
419 stepperSC->setSubcyclingMaxTimeStep (dt/10.0);
420 stepperSC->setSubcyclingMaxFailures (10);
421 stepperSC->setSubcyclingMaxConsecFailures(5);
422 stepperSC->setSubcyclingScreenOutputIndexInterval(1);
428 strategySC->setMinEta(0.000001);
429 strategySC->setMaxEta(0.01);
430 strategySC->initialize();
431 stepperSC->setSubcyclingTimeStepControlStrategy(strategySC);
438 stepper->addStepper(stepperSC);
439 stepper->addStepper(stepperBE);
443 timeStepControl->setInitIndex(0);
444 timeStepControl->setInitTime (0.0);
446 timeStepControl->setFinalTime(2.0);
447 timeStepControl->setMinTimeStep (0.000001);
448 timeStepControl->setInitTimeStep(dt);
449 timeStepControl->setMaxTimeStep (dt);
460 timeStepControl->initialize();
463 auto inArgsIC = stepper->getModel()->getNominalValues();
464 auto icX = rcp_const_cast<Thyra::VectorBase<double> > (inArgsIC.get_x());
465 auto icXDot = rcp_const_cast<Thyra::VectorBase<double> > (inArgsIC.get_x_dot());
467 icState->setTime (timeStepControl->getInitTime());
468 icState->setIndex (timeStepControl->getInitIndex());
469 icState->setTimeStep(0.0);
470 icState->setOrder (stepper->getOrder());
475 solutionHistory->setName(
"Forward States");
478 solutionHistory->setStorageLimit(3);
479 solutionHistory->addState(icState);
482 stepperSC->setInitialConditions(solutionHistory);
483 stepper->initialize();
486 integrator = Tempus::createIntegratorBasic<double>();
487 integrator->setStepper(stepper);
488 integrator->setTimeStepControl(timeStepControl);
489 integrator->setSolutionHistory(solutionHistory);
490 integrator->setScreenOutputIndexInterval(10);
492 integrator->initialize();
496 bool integratorStatus = integrator->advanceTime();
497 TEST_ASSERT(integratorStatus)
500 time = integrator->getTime();
501 double timeFinal = 2.0;
502 double tol = 100.0 * std::numeric_limits<double>::epsilon();
503 TEST_FLOATING_EQUALITY(time, timeFinal, tol);
506 StepSize.push_back(dt);
507 auto solution = Thyra::createMember(implicitModel->get_x_space());
508 Thyra::copy(*(integrator->getX()),solution.ptr());
509 solutions.push_back(solution);
510 auto solutionDot = Thyra::createMember(implicitModel->get_x_space());
511 Thyra::copy(*(integrator->getXDot()),solutionDot.ptr());
512 solutionsDot.push_back(solutionDot);
516 if ((n == 0) || (n == nTimeStepSizes-1)) {
517 std::string fname =
"Tempus_Subcycling_VanDerPol-Ref.dat";
518 if (n == 0) fname =
"Tempus_Subcycling_VanDerPol.dat";
526 double xDotSlope = 0.0;
527 RCP<Tempus::Stepper<double> > stepper = integrator->getStepper();
531 solutions, xErrorNorm, xSlope,
532 solutionsDot, xDotErrorNorm, xDotSlope);
534 TEST_FLOATING_EQUALITY( xSlope, 1.25708, 0.05 );
535 TEST_FLOATING_EQUALITY( xDotSlope, 1.20230, 0.05 );
536 TEST_FLOATING_EQUALITY( xErrorNorm[0], 0.37156, 1.0e-4 );
537 TEST_FLOATING_EQUALITY( xDotErrorNorm[0], 3.11651, 1.0e-4 );
539 Teuchos::TimeMonitor::summarize();
SolutionHistory is basically a container of SolutionStates. SolutionHistory maintains a collection of...