9 #ifndef Tempus_StepperRKButcherTableau_hpp
10 #define Tempus_StepperRKButcherTableau_hpp
14 #pragma clang system_header
17 #include "Tempus_config.hpp"
18 #include "Tempus_StepperExplicitRK.hpp"
19 #include "Tempus_StepperDIRK.hpp"
42 template<
class Scalar>
66 std::string ICConsistency,
67 bool ICConsistencyCheck,
74 this->
setup(appModel, useFSAL, ICConsistency,
75 ICConsistencyCheck, useEmbedded, stepperRKAppAction);
80 std::ostringstream Description;
85 return Description.str();
113 template<
class Scalar>
122 if (pl != Teuchos::null) {
124 pl->
get<std::string>(
"Stepper Type", stepper->getStepperType());
127 stepperType != stepper->getStepperType() &&
128 stepperType !=
"RK1", std::logic_error,
129 " ParameterList 'Stepper Type' (='" + stepperType +
"')\n"
130 " does not match type for this Stepper (='" + stepper->getStepperType() +
131 "')\n or one of its aliases ('RK1').\n");
134 pl->
set<std::string>(
"Stepper Type", stepper->getStepperType());
137 stepper->setStepperRKValues(pl);
139 if (model != Teuchos::null) {
140 stepper->setModel(model);
141 stepper->initialize();
167 template<
class Scalar>
191 std::string ICConsistency,
192 bool ICConsistencyCheck,
199 this->
setup(appModel, useFSAL, ICConsistency,
200 ICConsistencyCheck, useEmbedded, stepperRKAppAction);
205 std::ostringstream Description;
207 <<
"\"The\" Runge-Kutta Method (explicit):\n"
208 <<
"Solving Ordinary Differential Equations I:\n"
209 <<
"Nonstiff Problems, 2nd Revised Edition\n"
210 <<
"E. Hairer, S.P. Norsett, G. Wanner\n"
211 <<
"Table 1.2, pg 138\n"
212 <<
"c = [ 0 1/2 1/2 1 ]'\n"
217 <<
"b = [ 1/6 1/3 1/3 1/6 ]'";
218 return Description.str();
224 const Scalar one = ST::one();
225 const Scalar zero = ST::zero();
226 const Scalar onehalf = one/(2*one);
227 const Scalar onesixth = one/(6*one);
228 const Scalar onethird = one/(3*one);
236 A(0,0) = zero; A(0,1) = zero; A(0,2) = zero; A(0,3) = zero;
237 A(1,0) = onehalf; A(1,1) = zero; A(1,2) = zero; A(1,3) = zero;
238 A(2,0) = zero; A(2,1) = onehalf; A(2,2) = zero; A(2,3) = zero;
239 A(3,0) = zero; A(3,1) = zero; A(3,2) = one; A(3,3) = zero;
242 b(0) = onesixth; b(1) = onethird; b(2) = onethird; b(3) = onesixth;
245 c(0) = zero; c(1) = onehalf; c(2) = onehalf; c(3) = one;
257 template<
class Scalar>
264 stepper->setStepperRKValues(pl);
266 if (model != Teuchos::null) {
267 stepper->setModel(model);
268 stepper->initialize();
299 template<
class Scalar>
323 std::string ICConsistency,
324 bool ICConsistencyCheck,
331 this->
setup(appModel, useFSAL, ICConsistency,
332 ICConsistencyCheck, useEmbedded, stepperRKAppAction);
337 std::ostringstream Description;
339 <<
"P. Bogacki and L.F. Shampine.\n"
340 <<
"A 3(2) pair of Runge–Kutta formulas.\n"
341 <<
"Applied Mathematics Letters, 2(4):321 – 325, 1989.\n"
342 <<
"c = [ 0 1/2 3/4 1 ]'\n"
346 <<
" [ 2/9 1/3 4/9 0 ]\n"
347 <<
"b = [ 2/9 1/3 4/9 0 ]'\n"
348 <<
"bstar = [ 7/24 1/4 1/3 1/8 ]'";
349 return Description.str();
366 const Scalar one = ST::one();
367 const Scalar zero = ST::zero();
368 const Scalar onehalf = one/(2*one);
369 const Scalar onethird = one/(3*one);
370 const Scalar threefourths = (3*one)/(4*one);
371 const Scalar twoninths = (2*one)/(9*one);
372 const Scalar fourninths = (4*one)/(9*one);
375 A(0,0) = zero; A(0,1) = zero; A(0,2) = zero; A(0,3) = zero;
376 A(1,0) = onehalf; A(1,1) = zero; A(1,2) = zero; A(1,3) = zero;
377 A(2,0) = zero; A(2,1) =threefourths; A(2,2) = zero; A(2,3) = zero;
378 A(3,0) =twoninths; A(3,1) = onethird; A(3,2) =fourninths; A(3,3) = zero;
381 b(0) = A(3,0); b(1) = A(3,1); b(2) = A(3,2); b(3) = A(3,3);
384 c(0) = zero; c(1) = onehalf; c(2) = threefourths; c(3) = one;
387 bstar(0) = as<Scalar>(7*one/(24*one));
388 bstar(1) = as<Scalar>(1*one/(4*one));
389 bstar(2) = as<Scalar>(1*one/(3*one));
390 bstar(3) = as<Scalar>(1*one/(8*one));
401 template<
class Scalar>
408 stepper->setStepperRKValues(pl);
410 if (model != Teuchos::null) {
411 stepper->setModel(model);
412 stepper->initialize();
445 template<
class Scalar>
469 std::string ICConsistency,
470 bool ICConsistencyCheck,
477 this->
setup(appModel, useFSAL, ICConsistency,
478 ICConsistencyCheck, useEmbedded, stepperRKAppAction);
483 std::ostringstream Description;
485 <<
"Solving Ordinary Differential Equations I:\n"
486 <<
"Nonstiff Problems, 2nd Revised Edition\n"
487 <<
"E. Hairer, S.P. Norsett, G. Wanner\n"
488 <<
"Table 4.1, pg 167\n"
489 <<
"c = [ 0 1/3 1/3 1/2 1 ]'\n"
492 <<
" [ 1/6 1/6 0 ]\n"
493 <<
" [ 1/8 0 3/8 0 ]\n"
494 <<
" [ 1/2 0 -3/2 2 0 ]\n"
495 <<
"b = [ 1/6 0 0 2/3 1/6 ]'\n"
496 <<
"bstar = [ 1/10 0 3/10 2/5 1/5 ]'";
497 return Description.str();
513 const Scalar one = ST::one();
514 const Scalar zero = ST::zero();
517 A(1,0) = as<Scalar>(one/(3*one));;
519 A(2,0) = as<Scalar>(one/(6*one));;
520 A(2,1) = as<Scalar>(one/(6*one));;
522 A(3,0) = as<Scalar>(one/(8*one));;
523 A(3,2) = as<Scalar>(3*one/(8*one));;
525 A(4,0) = as<Scalar>(one/(2*one));;
526 A(4,2) = as<Scalar>(-3*one/(2*one));;
530 b(0) = as<Scalar>(one/(6*one));
531 b(3) = as<Scalar>(2*one/(3*one));
532 b(4) = as<Scalar>(one/(6*one));
536 c(1) = as<Scalar>(1*one/(3*one));
537 c(2) = as<Scalar>(1*one/(3*one));
538 c(3) = as<Scalar>(1*one/(2*one));
542 bstar(0) = as<Scalar>(1*one/(10*one));
543 bstar(2) = as<Scalar>(3*one/(10*one));
544 bstar(3) = as<Scalar>(2*one/(5*one));
545 bstar(4) = as<Scalar>(1*one/(5*one));
556 template<
class Scalar>
563 stepper->setStepperRKValues(pl);
565 if (model != Teuchos::null) {
566 stepper->setModel(model);
567 stepper->initialize();
597 template<
class Scalar>
617 std::string ICConsistency,
618 bool ICConsistencyCheck,
625 this->
setup(appModel, useFSAL, ICConsistency,
626 ICConsistencyCheck, useEmbedded, stepperRKAppAction);
631 std::ostringstream Description;
633 <<
"Solving Ordinary Differential Equations I:\n"
634 <<
"Nonstiff Problems, 2nd Revised Edition\n"
635 <<
"E. Hairer, S.P. Norsett, G. Wanner\n"
636 <<
"Table 1.2, pg 138\n"
637 <<
"c = [ 0 1/3 2/3 1 ]'\n"
642 <<
"b = [ 1/8 3/8 3/8 1/8 ]'";
643 return Description.str();
658 const Scalar one = ST::one();
659 const Scalar zero = ST::zero();
660 const Scalar onethird = as<Scalar>(one/(3*one));
661 const Scalar twothirds = as<Scalar>(2*one/(3*one));
662 const Scalar oneeighth = as<Scalar>(one/(8*one));
663 const Scalar threeeighths = as<Scalar>(3*one/(8*one));
666 A(0,0) = zero; A(0,1) = zero; A(0,2) = zero; A(0,3) = zero;
667 A(1,0) = onethird; A(1,1) = zero; A(1,2) = zero; A(1,3) = zero;
668 A(2,0) = -onethird; A(2,1) = one; A(2,2) = zero; A(2,3) = zero;
669 A(3,0) = one; A(3,1) = -one; A(3,2) = one; A(3,3) = zero;
672 b(0) =oneeighth; b(1) =threeeighths; b(2) =threeeighths; b(3) =oneeighth;
675 c(0) = zero; c(1) = onethird; c(2) = twothirds; c(3) = one;
687 template<
class Scalar>
694 stepper->setStepperRKValues(pl);
696 if (model != Teuchos::null) {
697 stepper->setModel(model);
698 stepper->initialize();
728 template<
class Scalar>
752 std::string ICConsistency,
753 bool ICConsistencyCheck,
760 this->
setup(appModel, useFSAL, ICConsistency,
761 ICConsistencyCheck, useEmbedded, stepperRKAppAction);
766 std::ostringstream Description;
768 <<
"Solving Ordinary Differential Equations I:\n"
769 <<
"Nonstiff Problems, 2nd Revised Edition\n"
770 <<
"E. Hairer, S.P. Norsett, G. Wanner\n"
771 <<
"Table 1.1, pg 135\n"
772 <<
"c = [ 0 1/2 1 1 ]'\n"
777 <<
"b = [ 1/6 2/3 0 1/6 ]'";
778 return Description.str();
790 const Scalar one = ST::one();
791 const Scalar onehalf = one/(2*one);
792 const Scalar onesixth = one/(6*one);
793 const Scalar twothirds = 2*one/(3*one);
794 const Scalar zero = ST::zero();
797 A(0,0) = zero; A(0,1) = zero; A(0,2) = zero; A(0,3) = zero;
798 A(1,0) = onehalf; A(1,1) = zero; A(1,2) = zero; A(1,3) = zero;
799 A(2,0) = zero; A(2,1) = one; A(2,2) = zero; A(2,3) = zero;
800 A(3,0) = zero; A(3,1) = zero; A(3,2) = one; A(3,3) = zero;
803 b(0) = onesixth; b(1) = twothirds; b(2) = zero; b(3) = onesixth;
806 c(0) = zero; c(1) = onehalf; c(2) = one; c(3) = one;
818 template<
class Scalar>
825 stepper->setStepperRKValues(pl);
827 if (model != Teuchos::null) {
828 stepper->setModel(model);
829 stepper->initialize();
858 template<
class Scalar>
870 this->
setStepperName(
"RK Explicit 5 Stage 3rd order by Kinnmark and Gray");
871 this->
setStepperType(
"RK Explicit 5 Stage 3rd order by Kinnmark and Gray");
882 std::string ICConsistency,
883 bool ICConsistencyCheck,
887 this->
setStepperName(
"RK Explicit 5 Stage 3rd order by Kinnmark and Gray");
888 this->
setStepperType(
"RK Explicit 5 Stage 3rd order by Kinnmark and Gray");
890 this->
setup(appModel, useFSAL, ICConsistency,
891 ICConsistencyCheck, useEmbedded,stepperRKAppAction);
896 std::ostringstream Description;
898 <<
"Kinnmark & Gray 5 stage, 3rd order scheme \n"
899 <<
"Modified by P. Ullrich. From the prim_advance_mod.F90 \n"
900 <<
"routine in the HOMME atmosphere model code.\n"
901 <<
"c = [ 0 1/5 1/5 1/3 2/3 ]'\n"
905 <<
" [ 0 0 1/3 0 ]\n"
906 <<
" [ 0 0 0 2/3 0 ]\n"
907 <<
"b = [ 1/4 0 0 0 3/4 ]'";
908 return Description.str();
921 const Scalar one = ST::one();
922 const Scalar onefifth = one/(5*one);
923 const Scalar onefourth = one/(4*one);
924 const Scalar onethird = one/(3*one);
925 const Scalar twothirds = 2*one/(3*one);
926 const Scalar threefourths = 3*one/(4*one);
927 const Scalar zero = ST::zero();
930 A(0,0) = zero; A(0,1) = zero; A(0,2) = zero; A(0,3) = zero; A(0,4) = zero;
931 A(1,0) = onefifth; A(1,1) = zero; A(1,2) = zero; A(1,3) = zero; A(1,4) = zero;
932 A(2,0) = zero; A(2,1) = onefifth; A(2,2) = zero; A(2,3) = zero; A(2,4) = zero;
933 A(3,0) = zero; A(3,1) = zero; A(3,2) = onethird; A(3,3) = zero; A(3,4) = zero;
934 A(4,0) = zero; A(4,1) = zero; A(4,2) = zero; A(4,3) = twothirds; A(4,4) = zero;
937 b(0) =onefourth; b(1) =zero; b(2) =zero; b(3) =zero; b(4) =threefourths;
940 c(0) =zero; c(1) =onefifth; c(2) =onefifth; c(3) =onethird; c(4) =twothirds;
952 template<
class Scalar>
959 stepper->setStepperRKValues(pl);
961 if (model != Teuchos::null) {
962 stepper->setModel(model);
963 stepper->initialize();
988 template<
class Scalar>
1012 std::string ICConsistency,
1013 bool ICConsistencyCheck,
1020 this->
setup(appModel, useFSAL, ICConsistency,
1021 ICConsistencyCheck, useEmbedded, stepperRKAppAction);
1026 std::ostringstream Description;
1028 <<
"c = [ 0 1/2 1 ]'\n"
1032 <<
"b = [ 1/6 4/6 1/6 ]'";
1033 return Description.str();
1041 const Scalar one = ST::one();
1042 const Scalar two = Teuchos::as<Scalar>(2*one);
1043 const Scalar zero = ST::zero();
1044 const Scalar onehalf = one/(2*one);
1045 const Scalar onesixth = one/(6*one);
1046 const Scalar foursixth = 4*one/(6*one);
1054 A(0,0) = zero; A(0,1) = zero; A(0,2) = zero;
1055 A(1,0) = onehalf; A(1,1) = zero; A(1,2) = zero;
1056 A(2,0) = -one; A(2,1) = two; A(2,2) = zero;
1059 b(0) = onesixth; b(1) = foursixth; b(2) = onesixth;
1062 c(0) = zero; c(1) = onehalf; c(2) = one;
1074 template<
class Scalar>
1081 stepper->setStepperRKValues(pl);
1083 if (model != Teuchos::null) {
1084 stepper->setModel(model);
1085 stepper->initialize();
1121 template<
class Scalar>
1145 std::string ICConsistency,
1146 bool ICConsistencyCheck,
1153 this->
setup(appModel, useFSAL, ICConsistency,
1154 ICConsistencyCheck, useEmbedded, stepperRKAppAction);
1159 std::ostringstream Description;
1161 <<
"This Stepper is known as 'RK Explicit 3 Stage 3rd order TVD' or 'SSPERK33'.\n"
1162 <<
"Sigal Gottlieb and Chi-Wang Shu\n"
1163 <<
"`Total Variation Diminishing Runge-Kutta Schemes'\n"
1164 <<
"Mathematics of Computation\n"
1165 <<
"Volume 67, Number 221, January 1998, pp. 73-85\n"
1166 <<
"c = [ 0 1 1/2 ]'\n"
1169 <<
" [ 1/4 1/4 0 ]\n"
1170 <<
"b = [ 1/6 1/6 4/6 ]'\n"
1171 <<
"This is also written in the following set of updates.\n"
1172 <<
"u1 = u^n + dt L(u^n)\n"
1173 <<
"u2 = 3 u^n/4 + u1/4 + dt L(u1)/4\n"
1174 <<
"u^(n+1) = u^n/3 + 2 u2/2 + 2 dt L(u2)/3";
1175 return Description.str();
1184 const Scalar one = ST::one();
1185 const Scalar zero = ST::zero();
1186 const Scalar onehalf = one/(2*one);
1187 const Scalar onefourth = one/(4*one);
1188 const Scalar onesixth = one/(6*one);
1189 const Scalar foursixth = 4*one/(6*one);
1198 A(0,0) = zero; A(0,1) = zero; A(0,2) = zero;
1199 A(1,0) = one; A(1,1) = zero; A(1,2) = zero;
1200 A(2,0) = onefourth; A(2,1) = onefourth; A(2,2) = zero;
1203 b(0) = onesixth; b(1) = onesixth; b(2) = foursixth;
1206 c(0) = zero; c(1) = one; c(2) = onehalf;
1209 bstar(0) = as<Scalar>(0.291485418878409);
1210 bstar(1) = as<Scalar>(0.291485418878409);
1211 bstar(2) = as<Scalar>(0.417029162243181);
1225 template<
class Scalar>
1234 if (pl != Teuchos::null) {
1236 pl->
get<std::string>(
"Stepper Type", stepper->getStepperType());
1239 stepperType != stepper->getStepperType() &&
1240 stepperType !=
"SSPERK33" && stepperType !=
"SSPRK3", std::logic_error,
1241 " ParameterList 'Stepper Type' (='" + stepperType +
"')\n"
1242 " does not match type for this Stepper (='" + stepper->getStepperType() +
1243 "')\n or one of its aliases ('SSPERK33' or 'SSPRK3').\n");
1246 pl->
set<std::string>(
"Stepper Type", stepper->getStepperType());
1249 stepper->setStepperRKValues(pl);
1251 if (model != Teuchos::null) {
1252 stepper->setModel(model);
1253 stepper->initialize();
1282 template<
class Scalar>
1306 std::string ICConsistency,
1307 bool ICConsistencyCheck,
1314 this->
setup(appModel, useFSAL, ICConsistency,
1315 ICConsistencyCheck, useEmbedded, stepperRKAppAction);
1320 std::ostringstream Description;
1322 <<
"Solving Ordinary Differential Equations I:\n"
1323 <<
"Nonstiff Problems, 2nd Revised Edition\n"
1324 <<
"E. Hairer, S.P. Norsett, G. Wanner\n"
1325 <<
"Table 1.1, pg 135\n"
1326 <<
"c = [ 0 1/3 2/3 ]'\n"
1330 <<
"b = [ 1/4 0 3/4 ]'";
1331 return Description.str();
1339 const Scalar one = ST::one();
1340 const Scalar zero = ST::zero();
1341 const Scalar onethird = one/(3*one);
1342 const Scalar twothirds = 2*one/(3*one);
1343 const Scalar onefourth = one/(4*one);
1344 const Scalar threefourths = 3*one/(4*one);
1352 A(0,0) = zero; A(0,1) = zero; A(0,2) = zero;
1353 A(1,0) = onethird; A(1,1) = zero; A(1,2) = zero;
1354 A(2,0) = zero; A(2,1) = twothirds; A(2,2) = zero;
1357 b(0) = onefourth; b(1) = zero; b(2) = threefourths;
1360 c(0) = zero; c(1) = onethird; c(2) = twothirds;
1372 template<
class Scalar>
1379 stepper->setStepperRKValues(pl);
1381 if (model != Teuchos::null) {
1382 stepper->setModel(model);
1383 stepper->initialize();
1411 template<
class Scalar>
1435 std::string ICConsistency,
1436 bool ICConsistencyCheck,
1443 this->
setup(appModel, useFSAL, ICConsistency,
1444 ICConsistencyCheck, useEmbedded, stepperRKAppAction);
1449 std::ostringstream Description;
1451 <<
"Solving Ordinary Differential Equations I:\n"
1452 <<
"Nonstiff Problems, 2nd Revised Edition\n"
1453 <<
"E. Hairer, S.P. Norsett, G. Wanner\n"
1454 <<
"Table 1.1, pg 135\n"
1455 <<
"c = [ 0 1/2 ]'\n"
1459 return Description.str();
1467 const Scalar one = ST::one();
1468 const Scalar zero = ST::zero();
1469 const Scalar onehalf = one/(2*one);
1477 A(0,0) = zero; A(0,1) = zero;
1478 A(1,0) = onehalf; A(1,1) = zero;
1481 b(0) = zero; b(1) = one;
1484 c(0) = zero; c(1) = onehalf;
1496 template<
class Scalar>
1503 stepper->setStepperRKValues(pl);
1505 if (model != Teuchos::null) {
1506 stepper->setModel(model);
1507 stepper->initialize();
1531 template<
class Scalar>
1555 std::string ICConsistency,
1556 bool ICConsistencyCheck,
1563 this->
setup(appModel, useFSAL, ICConsistency,
1564 ICConsistencyCheck, useEmbedded, stepperRKAppAction);
1569 std::ostringstream Description;
1571 <<
"This Stepper is known as 'RK Explicit Ralston' or 'RK2'.\n"
1572 <<
"c = [ 0 2/3 ]'\n"
1575 <<
"b = [ 1/4 3/4 ]'";
1576 return Description.str();
1584 const Scalar one = ST::one();
1585 const Scalar zero = ST::zero();
1587 const int NumStages = 2;
1588 const int order = 2;
1594 A(0,0) = zero; A(0,1) = zero; A(1,1) = zero;
1595 A(1,0) = (2*one)/(3*one);
1598 b(0) = (1*one)/(4*one);
1599 b(1) = (3*one)/(4*one);
1603 c(1) = (2*one)/(3*one);
1616 template<
class Scalar>
1625 if (pl != Teuchos::null) {
1627 pl->
get<std::string>(
"Stepper Type", stepper->getStepperType());
1630 stepperType != stepper->getStepperType() &&
1631 stepperType !=
"RK2", std::logic_error,
1632 " ParameterList 'Stepper Type' (='" + stepperType +
"')\n"
1633 " does not match type for this Stepper (='" + stepper->getStepperType() +
1634 "')\n or one of its aliases ('RK2').\n");
1637 pl->
set<std::string>(
"Stepper Type", stepper->getStepperType());
1640 stepper->setStepperRKValues(pl);
1642 if (model != Teuchos::null) {
1643 stepper->setModel(model);
1644 stepper->initialize();
1669 template<
class Scalar>
1693 std::string ICConsistency,
1694 bool ICConsistencyCheck,
1701 this->
setup(appModel, useFSAL, ICConsistency,
1702 ICConsistencyCheck, useEmbedded, stepperRKAppAction);
1707 std::ostringstream Description;
1709 <<
"This Stepper is known as 'RK Explicit Trapezoidal' or 'Heuns Method' or 'SSPERK22'.\n"
1713 <<
"b = [ 1/2 1/2 ]\n"
1714 <<
"bstart = [ 3/4 1/4 ]'";
1715 return Description.str();
1724 const Scalar one = ST::one();
1725 const Scalar zero = ST::zero();
1726 const Scalar onehalf = one/(2*one);
1735 A(0,0) = zero; A(0,1) = zero;
1736 A(1,0) = one; A(1,1) = zero;
1739 b(0) = onehalf; b(1) = onehalf;
1742 c(0) = zero; c(1) = one;
1745 bstar(0) = as<Scalar>(3*one/(4*one));
1746 bstar(1) = as<Scalar>(1*one/(4*one));
1760 template<
class Scalar>
1769 if (pl != Teuchos::null) {
1771 pl->
get<std::string>(
"Stepper Type", stepper->getStepperType());
1774 stepperType != stepper->getStepperType() &&
1775 stepperType !=
"Heuns Method" && stepperType !=
"SSPERK22" &&
1776 stepperType !=
"SSPRK2", std::logic_error,
1777 " ParameterList 'Stepper Type' (='" + stepperType +
"')\n"
1778 " does not match type for this Stepper (='" + stepper->getStepperType() +
1779 "')\n or one of its aliases ('Heuns Method' or 'SSPERK22' or 'SSPRK2').\n");
1782 pl->
set<std::string>(
"Stepper Type", stepper->getStepperType());
1785 stepper->setStepperRKValues(pl);
1787 if (model != Teuchos::null) {
1788 stepper->setModel(model);
1789 stepper->initialize();
1813 template<
class Scalar>
1832 std::string ICConsistency,
1833 bool ICConsistencyCheck,
1840 this->
setup(appModel, useFSAL, ICConsistency,
1841 ICConsistencyCheck, useEmbedded, stepperRKAppAction);
1846 std::ostringstream Description;
1848 <<
"Strong Stability Preserving Explicit RK (stage=5, order=4)\n"
1850 return Description.str();
1860 const int NumStages = 5;
1861 const int order = 4;
1862 const Scalar sspcoef = 1.5082;
1867 const Scalar zero = ST::zero();
1870 A(0,0) = A(0,1) = A(0,2) = A(0,3) = A(0,4) = zero;
1872 A(1,0) = as<Scalar>(0.391752226571889);
1873 A(1,1) = A(1,2) = A(1,3) = A(0,4) = zero;
1875 A(2,0) = as<Scalar>(0.217669096261169);
1876 A(2,1) = as<Scalar>(0.368410593050372);
1877 A(2,2) = A(2,3) = A(2,4) = zero;
1879 A(3,0) = as<Scalar>(0.082692086657811);
1880 A(3,1) = as<Scalar>(0.139958502191896);
1881 A(3,2) = as<Scalar>(0.251891774271693);
1882 A(3,3) = A(2,4) = zero;
1884 A(4,0) = as<Scalar>(0.067966283637115);
1885 A(4,1) = as<Scalar>(0.115034698504632);
1886 A(4,2) = as<Scalar>(0.207034898597385);
1887 A(4,3) = as<Scalar>(0.544974750228520);
1891 b(0) = as<Scalar>(0.146811876084786);
1892 b(1) = as<Scalar>(0.248482909444976);
1893 b(2) = as<Scalar>(0.104258830331980);
1894 b(3) = as<Scalar>(0.274438900901350);
1895 b(4) = as<Scalar>(0.226007483236908);
1900 c(2) = A(2,0) + A(2,1);
1901 c(3) = A(3,0) + A(3,1) + A(3,2);
1902 c(4) = A(4,0) + A(4,1) + A(4,2) + A(4,3);
1905 bstar(0) = as<Scalar>(0.130649104813131);
1906 bstar(1) = as<Scalar>(0.317716031201302);
1907 bstar(2) = as<Scalar>(0.000000869337261);
1908 bstar(3) = as<Scalar>(0.304581512634772);
1909 bstar(4) = as<Scalar>(0.247052482013534);
1914 this->
tableau_->setTVDCoeff(sspcoef);
1921 template<
class Scalar>
1928 stepper->setStepperRKValues(pl);
1930 if (model != Teuchos::null) {
1931 stepper->setModel(model);
1932 stepper->initialize();
1969 template<
class Scalar>
1988 std::string ICConsistency,
1989 bool ICConsistencyCheck,
2002 this->
setTableau(A,b,c,order,orderMin,orderMax,bstar);
2005 this->
tableau_->isImplicit() ==
true, std::logic_error,
2006 "Error - General ERK received an implicit Butcher Tableau!\n");
2008 this->
setup(appModel, useFSAL, ICConsistency,
2009 ICConsistencyCheck, useEmbedded, stepperRKAppAction);
2014 std::stringstream Description;
2016 <<
"The format of the Butcher Tableau parameter list is\n"
2017 <<
" <Parameter name=\"A\" type=\"string\" value=\"# # # ;\n"
2020 <<
" <Parameter name=\"b\" type=\"string\" value=\"# # #\"/>\n"
2021 <<
" <Parameter name=\"c\" type=\"string\" value=\"# # #\"/>\n\n"
2022 <<
"Note the number of stages is implicit in the number of entries.\n"
2023 <<
"The number of stages must be consistent.\n"
2025 <<
"Default tableau is RK4 (order=4):\n"
2026 <<
"c = [ 0 1/2 1/2 1 ]'\n"
2031 <<
"b = [ 1/6 1/3 1/3 1/6 ]'";
2032 return Description.str();
2037 if (this->
tableau_ == Teuchos::null) {
2040 auto t = stepper->getTableau();
2043 t->A(),t->b(),t->c(),
2044 t->order(),t->orderMin(),t->orderMax(),
2078 std::ostringstream Astream;
2079 Astream.precision(15);
2080 for(
int i = 0; i < A.
numRows(); i++) {
2081 for(
int j = 0; j < A.
numCols()-1; j++) {
2082 Astream << A(i,j) <<
" ";
2084 Astream << A(i,A.
numCols()-1);
2085 if ( i != A.
numRows()-1 ) Astream <<
"; ";
2087 tableauPL->
set<std::string>(
"A", Astream.str());
2089 std::ostringstream bstream;
2090 bstream.precision(15);
2091 for(
int i = 0; i < b.
length()-1; i++) {
2092 bstream << b(i) <<
" ";
2094 bstream << b(b.
length()-1);
2095 tableauPL->
set<std::string>(
"b", bstream.str());
2097 std::ostringstream cstream;
2098 cstream.precision(15);
2099 for(
int i = 0; i < c.
length()-1; i++) {
2100 cstream << c(i) <<
" ";
2102 cstream << c(c.
length()-1);
2103 tableauPL->
set<std::string>(
"c", cstream.str());
2105 tableauPL->
set<
int>(
"order", this->
tableau_->order());
2107 if ( bstar.
length() == 0 ) {
2108 tableauPL->
set(
"bstar",
"");
2110 std::ostringstream bstarstream;
2111 bstarstream.precision(15);
2112 for(
int i = 0; i < bstar.
length()-1; i++) {
2113 bstarstream << bstar(i) <<
" ";
2115 bstarstream << bstar(bstar.
length()-1);
2116 tableauPL->
set<std::string>(
"bstar", bstarstream.str());
2119 pl->
set(
"Tableau", *tableauPL);
2128 template<
class Scalar>
2135 stepper->setStepperRKValues(pl);
2137 if (pl != Teuchos::null) {
2139 auto t = stepper->createTableau(pl);
2140 stepper->setTableau( t->A(),t->b(),t->c(),
2141 t->order(),t->orderMin(),t->orderMax(),
2146 stepper->getTableau()->isImplicit() ==
true, std::logic_error,
2147 "Error - General ERK received an implicit Butcher Tableau!\n");
2149 if (model != Teuchos::null) {
2150 stepper->setModel(model);
2151 stepper->initialize();
2174 template<
class Scalar>
2199 std::string ICConsistency,
2200 bool ICConsistencyCheck,
2202 bool zeroInitialGuess,
2208 this->
setup(appModel, solver, useFSAL, ICConsistency, ICConsistencyCheck,
2209 useEmbedded, zeroInitialGuess, stepperRKAppAction);
2214 std::ostringstream Description;
2219 return Description.str();
2227 const Scalar sspcoef = std::numeric_limits<Scalar>::max();
2247 this->
tableau_->setTVDCoeff(sspcoef);
2254 template<
class Scalar>
2261 stepper->setStepperDIRKValues(pl);
2263 if (model != Teuchos::null) {
2264 stepper->setModel(model);
2265 stepper->initialize();
2297 template<
class Scalar>
2310 const Scalar one = ST::one();
2311 gammaDefault_ = Teuchos::as<Scalar>((2*one-ST::squareroot(2*one))/(2*one));
2327 std::string ICConsistency,
2328 bool ICConsistencyCheck,
2330 bool zeroInitialGuess,
2332 Scalar gamma = Scalar(0.2928932188134524))
2335 const Scalar one = ST::one();
2336 gammaDefault_ = Teuchos::as<Scalar>((2*one-ST::squareroot(2*one))/(2*one));
2342 this->
setup(appModel, solver, useFSAL, ICConsistency, ICConsistencyCheck,
2343 useEmbedded, zeroInitialGuess, stepperRKAppAction);
2357 std::ostringstream Description;
2359 <<
"Computer Methods for ODEs and DAEs\n"
2360 <<
"U. M. Ascher and L. R. Petzold\n"
2362 <<
"gamma = (2+-sqrt(2))/2\n"
2363 <<
"c = [ gamma 1 ]'\n"
2364 <<
"A = [ gamma 0 ]\n"
2365 <<
" [ 1-gamma gamma ]\n"
2366 <<
"b = [ 1-gamma gamma ]'";
2367 return Description.str();
2374 pl->template set<double>(
"gamma", this->
getGamma(),
2375 "The default value is gamma = (2-sqrt(2))/2. "
2376 "This will produce an L-stable 2nd order method with the stage "
2377 "times within the timestep. Other values of gamma will still "
2378 "produce an L-stable scheme, but will only be 1st order accurate.");
2393 const Scalar one = ST::one();
2394 const Scalar zero = ST::zero();
2397 A(0,0) =
gamma_; A(0,1) = zero;
2398 A(1,0) = Teuchos::as<Scalar>( one -
gamma_ ); A(1,1) =
gamma_;
2401 b(0) = Teuchos::as<Scalar>( one -
gamma_ ); b(1) =
gamma_;
2404 c(0) =
gamma_; c(1) = one;
2421 template<
class Scalar>
2428 stepper->setStepperDIRKValues(pl);
2430 if (pl != Teuchos::null)
2431 stepper->setGamma(pl->
get<
double>(
"gamma", 0.2928932188134524));
2433 if (model != Teuchos::null) {
2434 stepper->setModel(model);
2435 stepper->initialize();
2470 template<
class Scalar>
2495 std::string ICConsistency,
2496 bool ICConsistencyCheck,
2498 bool zeroInitialGuess,
2505 this->
setup(appModel, solver, useFSAL, ICConsistency, ICConsistencyCheck,
2506 useEmbedded, zeroInitialGuess, stepperRKAppAction);
2511 std::ostringstream Description;
2513 <<
"Implicit-explicit Runge-Kutta schemes and applications to\n"
2514 <<
"hyperbolic systems with relaxation\n"
2515 <<
"L Pareschi, G Russo\n"
2516 <<
"Journal of Scientific computing, 2005 - Springer\n"
2518 <<
"gamma = 1/(2+sqrt(2))\n"
2519 <<
"c = [ gamma (1-gamma) 1/2 ]'\n"
2520 <<
"A = [ gamma 0 0 ]\n"
2521 <<
" [ 1-2gamma gamma 0 ]\n"
2522 <<
" [ 1/2-gamma 0 gamma ]\n"
2523 <<
"b = [ 1/6 1/6 2/3 ]'";
2524 return Description.str();
2533 const int NumStages = 3;
2534 const int order = 2;
2535 const Scalar sspcoef = 1.0529;
2539 const Scalar one = ST::one();
2540 const Scalar zero = ST::zero();
2541 const Scalar gamma = as<Scalar>(one - ( one / ST::squareroot(2*one) ) );
2544 A(0,0) = A(1,1) = A(2,2) = gamma;
2545 A(0,1) = A(0,2) = A(1,2) = A(2,1) = zero;
2546 A(1,0) = as<Scalar>(one - 2*gamma);
2547 A(2,0) = as<Scalar>( ( one/ (2.*one)) - gamma );
2550 b(0) = b(1) = ( one / (6*one) );
2551 b(2) = (2*one)/(3*one);
2556 c(2) = one / (2*one);
2561 this->
tableau_->setTVDCoeff(sspcoef);
2569 template<
class Scalar>
2576 stepper->setStepperDIRKValues(pl);
2578 if (model != Teuchos::null) {
2579 stepper->setModel(model);
2580 stepper->initialize();
2616 template<
class Scalar>
2630 const Scalar one = ST::one();
2631 gammaDefault_ = as<Scalar>((3*one+ST::squareroot(3*one))/(6*one));
2649 std::string ICConsistency,
2650 bool ICConsistencyCheck,
2652 bool zeroInitialGuess,
2654 std::string gammaType =
"3rd Order A-stable",
2655 Scalar gamma = Scalar(0.7886751345948128))
2659 const Scalar one = ST::one();
2660 gammaDefault_ = as<Scalar>((3*one+ST::squareroot(3*one))/(6*one));
2668 this->
setup(appModel, solver, useFSAL, ICConsistency, ICConsistencyCheck,
2669 useEmbedded, zeroInitialGuess, stepperRKAppAction);
2675 !(gammaType ==
"3rd Order A-stable" ||
2676 gammaType ==
"2nd Order L-stable" ||
2677 gammaType ==
"gamma"), std::logic_error,
2678 "gammaType needs to be '3rd Order A-stable', '2nd Order L-stable' "
2701 std::ostringstream Description;
2703 <<
"Solving Ordinary Differential Equations I:\n"
2704 <<
"Nonstiff Problems, 2nd Revised Edition\n"
2705 <<
"E. Hairer, S. P. Norsett, and G. Wanner\n"
2706 <<
"Table 7.2, pg 207\n"
2707 <<
"gamma = (3+sqrt(3))/6 -> 3rd order and A-stable\n"
2708 <<
"gamma = (2-sqrt(2))/2 -> 2nd order and L-stable\n"
2709 <<
"c = [ gamma 1-gamma ]'\n"
2710 <<
"A = [ gamma 0 ]\n"
2711 <<
" [ 1-2*gamma gamma ]\n"
2712 <<
"b = [ 1/2 1/2 ]'";
2713 return Description.str();
2721 pl->template set<std::string>(
"Gamma Type", this->
getGammaType(),
2722 "Valid values are '3rd Order A-stable' ((3+sqrt(3))/6.), "
2723 "'2nd Order L-stable' ((2-sqrt(2))/2) and 'gamma' (user defined). "
2724 "The default value is '3rd Order A-stable'.");
2725 pl->template set<double>(
"gamma", this->
getGamma(),
2726 "Equal to (3+sqrt(3))/6 if 'Gamma Type' = '3rd Order A-stable', or "
2727 "(2-sqrt(2))/2 if 'Gamma Type' = '2nd Order L-stable', or "
2728 "user-defined gamma value if 'Gamma Type = 'gamma'. "
2729 "The default value is gamma = (3+sqrt(3))/6, which matches "
2730 "the default 'Gamma Type' = '3rd Order A-stable'.");
2745 const Scalar one = ST::one();
2746 const Scalar zero = ST::zero();
2752 }
else if (
gammaType_ ==
"2nd Order L-stable") {
2754 gamma_ = as<Scalar>( (2*one - ST::squareroot(2*one))/(2*one) );
2760 A(0,0) =
gamma_; A(0,1) = zero;
2764 b(0) = as<Scalar>( one/(2*one) ); b(1) = as<Scalar>( one/(2*one) );
2783 template<
class Scalar>
2790 stepper->setStepperDIRKValues(pl);
2792 if (pl != Teuchos::null) {
2793 stepper->setGammaType(
2794 pl->
get<std::string>(
"Gamma Type",
"3rd Order A-stable"));
2795 stepper->setGamma(pl->
get<
double>(
"gamma", 0.7886751345948128));
2798 if (model != Teuchos::null) {
2799 stepper->setModel(model);
2800 stepper->initialize();
2828 template<
class Scalar>
2853 std::string ICConsistency,
2854 bool ICConsistencyCheck,
2856 bool zeroInitialGuess,
2862 this->
setup(appModel, solver, useFSAL, ICConsistency, ICConsistencyCheck,
2863 useEmbedded, zeroInitialGuess, stepperRKAppAction);
2868 std::ostringstream Description;
2870 <<
"Hammer & Hollingsworth method\n"
2871 <<
"Solving Ordinary Differential Equations I:\n"
2872 <<
"Nonstiff Problems, 2nd Revised Edition\n"
2873 <<
"E. Hairer, S. P. Norsett, and G. Wanner\n"
2874 <<
"Table 7.1, pg 205\n"
2875 <<
"c = [ 0 2/3 ]'\n"
2878 <<
"b = [ 1/4 3/4 ]'";
2879 return Description.str();
2892 const Scalar one = ST::one();
2893 const Scalar zero = ST::zero();
2896 A(0,0) = zero; A(0,1) = zero;
2897 A(1,0) = as<Scalar>( one/(3*one) ); A(1,1) = as<Scalar>( one/(3*one) );
2900 b(0) = as<Scalar>( one/(4*one) ); b(1) = as<Scalar>( 3*one/(4*one) );
2903 c(0) = zero; c(1) = as<Scalar>( 2*one/(3*one) );
2915 template<
class Scalar>
2922 stepper->setStepperDIRKValues(pl);
2924 if (model != Teuchos::null) {
2925 stepper->setModel(model);
2926 stepper->initialize();
2961 template<
class Scalar>
2990 std::string ICConsistency,
2991 bool ICConsistencyCheck,
2993 bool zeroInitialGuess,
2995 Scalar theta = Scalar(0.5))
3004 this->
setup(appModel, solver, useFSAL, ICConsistency, ICConsistencyCheck,
3005 useEmbedded, zeroInitialGuess, stepperRKAppAction);
3012 "'theta' can not be zero, as it makes this stepper explicit. \n"
3013 "Try using the 'RK Forward Euler' stepper.\n");
3023 std::ostringstream Description;
3025 <<
"Non-standard finite-difference methods\n"
3026 <<
"in dynamical systems, P. Kama,\n"
3027 <<
"Dissertation, University of Pretoria, pg. 49.\n"
3028 <<
"Comment: Generalized Implicit Midpoint Method\n"
3029 <<
"c = [ theta ]'\n"
3030 <<
"A = [ theta ]\n"
3032 return Description.str();
3040 pl->template set<double>(
"theta",
getTheta(),
3041 "Valid values are 0 <= theta <= 1, where theta = 0 "
3042 "implies Forward Euler, theta = 1/2 implies implicit midpoint "
3043 "method (default), and theta = 1 implies Backward Euler. "
3044 "For theta != 1/2, this method is first-order accurate, "
3045 "and with theta = 1/2, it is second-order accurate. "
3046 "This method is A-stable, but becomes L-stable with theta=1.");
3081 template<
class Scalar>
3088 stepper->setStepperDIRKValues(pl);
3090 if (pl != Teuchos::null) {
3091 stepper->setTheta(pl->
get<
double>(
"theta", 0.5));
3094 if (model != Teuchos::null) {
3095 stepper->setModel(model);
3096 stepper->initialize();
3130 template<
class Scalar>
3159 std::string ICConsistency,
3160 bool ICConsistencyCheck,
3162 bool zeroInitialGuess,
3164 Scalar theta = Scalar(0.5))
3173 this->
setup(appModel, solver, useFSAL, ICConsistency, ICConsistencyCheck,
3174 useEmbedded, zeroInitialGuess, stepperRKAppAction);
3181 "'theta' can not be zero, as it makes this stepper explicit. \n"
3182 "Try using the 'RK Forward Euler' stepper.\n");
3192 std::ostringstream Description;
3194 <<
"Computer Methods for ODEs and DAEs\n"
3195 <<
"U. M. Ascher and L. R. Petzold\n"
3199 <<
" [ 1-theta theta ]\n"
3200 <<
"b = [ 1-theta theta ]'";
3201 return Description.str();
3209 pl->template set<double>(
"theta",
getTheta(),
3210 "Valid values are 0 < theta <= 1, where theta = 0 "
3211 "implies Forward Euler, theta = 1/2 implies trapezoidal "
3212 "method (default), and theta = 1 implies Backward Euler. "
3213 "For theta != 1/2, this method is first-order accurate, "
3214 "and with theta = 1/2, it is second-order accurate. "
3215 "This method is A-stable, but becomes L-stable with theta=1.");
3227 const Scalar one = ST::one();
3228 const Scalar zero = ST::zero();
3236 A(0,0) = zero; A(0,1) = zero;
3237 A(1,0) = Teuchos::as<Scalar>( one -
theta_ ); A(1,1) =
theta_;
3240 b(0) = Teuchos::as<Scalar>( one -
theta_ );
3264 template<
class Scalar>
3271 stepper->setStepperDIRKValues(pl);
3273 if (pl != Teuchos::null) {
3274 stepper->setTheta(pl->
get<
double>(
"theta", 0.5));
3277 if (model != Teuchos::null) {
3278 stepper->setModel(model);
3279 stepper->initialize();
3307 template<
class Scalar>
3332 std::string ICConsistency,
3333 bool ICConsistencyCheck,
3335 bool zeroInitialGuess,
3341 this->
setup(appModel, solver, useFSAL, ICConsistency, ICConsistencyCheck,
3342 useEmbedded, zeroInitialGuess, stepperRKAppAction);
3347 std::ostringstream Description;
3349 <<
"Also known as Crank-Nicolson Method.\n"
3353 <<
"b = [ 1/2 1/2 ]'";
3354 return Description.str();
3364 const Scalar one = ST::one();
3365 const Scalar zero = ST::zero();
3366 const Scalar onehalf = ST::one()/(2*ST::one());
3374 A(0,0) = zero; A(0,1) = zero;
3375 A(1,0) = onehalf; A(1,1) = onehalf;
3397 template<
class Scalar>
3406 if (pl != Teuchos::null) {
3408 pl->
get<std::string>(
"Stepper Type", stepper->getStepperType());
3411 stepperType != stepper->getStepperType() &&
3412 stepperType !=
"RK Crank-Nicolson", std::logic_error,
3413 " ParameterList 'Stepper Type' (='" + stepperType +
"')\n"
3414 " does not match type for this Stepper (='" + stepper->getStepperType() +
3415 "')\n or one of its aliases ('RK Crank-Nicolson').\n");
3418 pl->
set<std::string>(
"Stepper Type", stepper->getStepperType());
3421 stepper->setStepperDIRKValues(pl);
3423 if (model != Teuchos::null) {
3424 stepper->setModel(model);
3425 stepper->initialize();
3459 template<
class Scalar>
3484 std::string ICConsistency,
3485 bool ICConsistencyCheck,
3487 bool zeroInitialGuess,
3493 this->
setup(appModel, solver, useFSAL, ICConsistency, ICConsistencyCheck,
3494 useEmbedded, zeroInitialGuess, stepperRKAppAction);
3499 std::ostringstream Description;
3502 <<
"Solving Ordinary Differential Equations II:\n"
3503 <<
"Stiff and Differential-Algebraic Problems,\n"
3504 <<
"2nd Revised Edition\n"
3505 <<
"E. Hairer and G. Wanner\n"
3506 <<
"Table 5.2, pg 72\n"
3507 <<
"Solving Ordinary Differential Equations I:\n"
3508 <<
"Nonstiff Problems, 2nd Revised Edition\n"
3509 <<
"E. Hairer, S. P. Norsett, and G. Wanner\n"
3510 <<
"Table 7.1, pg 205\n"
3514 return Description.str();
3526 const Scalar onehalf = ST::one()/(2*ST::one());
3527 const Scalar one = ST::one();
3550 template<
class Scalar>
3557 stepper->setStepperDIRKValues(pl);
3559 if (model != Teuchos::null) {
3560 stepper->setModel(model);
3561 stepper->initialize();
3589 template<
class Scalar>
3609 std::string ICConsistency,
3610 bool ICConsistencyCheck,
3612 bool zeroInitialGuess,
3618 this->
setup(appModel, solver, useFSAL, ICConsistency, ICConsistencyCheck,
3619 useEmbedded, zeroInitialGuess, stepperRKAppAction);
3624 std::ostringstream Description;
3626 <<
"Strong Stability Preserving Diagonally-Implicit RK (stage=2, order=2)\n"
3628 <<
"c = [ 1/4 3/4 ]'\n"
3631 <<
"b = [ 1/2 1/2 ]\n" << std::endl;
3632 return Description.str();
3641 const int NumStages = 2;
3642 const int order = 2;
3647 const Scalar one = ST::one();
3648 const Scalar zero = ST::zero();
3649 const Scalar onehalf = one/(2*one);
3650 const Scalar onefourth = one/(4*one);
3653 A(0,0) = A(1,1) = onefourth;
3658 b(0) = b(1) = onehalf;
3662 c(1) = A(1,0) + A(1,1);
3674 template<
class Scalar>
3681 stepper->setStepperDIRKValues(pl);
3683 if (model != Teuchos::null) {
3684 stepper->setModel(model);
3685 stepper->initialize();
3714 template<
class Scalar>
3734 std::string ICConsistency,
3735 bool ICConsistencyCheck,
3737 bool zeroInitialGuess,
3743 this->
setup(appModel, solver, useFSAL, ICConsistency, ICConsistencyCheck,
3744 useEmbedded, zeroInitialGuess, stepperRKAppAction);
3749 std::ostringstream Description;
3751 <<
"Strong Stability Preserving Diagonally-Implicit RK (stage=3, order=2)\n"
3753 <<
"c = [ 1/6 1/2 5/6 ]'\n"
3756 <<
" [ 1/3 1/3 1/6 ]\n"
3757 <<
"b = [ 1/3 1/3 1/3 ]\n" << std::endl;
3758 return Description.str();
3768 const int NumStages = 3;
3769 const int order = 2;
3774 const Scalar one = ST::one();
3775 const Scalar zero = ST::zero();
3776 const Scalar onethird = one/(3*one);
3777 const Scalar onesixth = one/(6*one);
3780 A(0,0) = A(1,1) = A(2,2) = onesixth;
3781 A(1,0) = A(2,0) = A(2,1) = onethird;
3782 A(0,1) = A(0,2) = A(1,2) = zero;
3785 b(0) = b(1) = b(2) = onethird;
3789 c(1) = A(1,0) + A(1,1);
3790 c(2) = A(2,0) + A(2,1) + A(2,2);
3802 template<
class Scalar>
3809 stepper->setStepperDIRKValues(pl);
3811 if (model != Teuchos::null) {
3812 stepper->setModel(model);
3813 stepper->initialize();
3841 template<
class Scalar>
3861 std::string ICConsistency,
3862 bool ICConsistencyCheck,
3864 bool zeroInitialGuess,
3870 this->
setup(appModel, solver, useFSAL, ICConsistency, ICConsistencyCheck,
3871 useEmbedded, zeroInitialGuess, stepperRKAppAction);
3876 std::ostringstream Description;
3878 <<
"Strong Stability Preserving Diagonally-Implicit RK (stage=2, order=3)\n"
3879 <<
"SSP-Coef = 1 + sqrt( 3 )\n"
3880 <<
"c = [ 1/(3 + sqrt( 3 )) (1/6)(3 + sqrt( 3 )) ] '\n"
3881 <<
"A = [ 1/(3 + sqrt( 3 )) ] \n"
3882 <<
" [ 1/sqrt( 3 ) 1/(3 + sqrt( 3 )) ] \n"
3883 <<
"b = [ 1/2 1/2 ] \n" << std::endl;
3884 return Description.str();
3894 const int NumStages = 2;
3895 const int order = 3;
3896 const Scalar sspcoef = 2.7321;
3901 const Scalar one = ST::one();
3902 const Scalar zero = ST::zero();
3903 const Scalar onehalf = one/(2*one);
3904 const Scalar rootthree = ST::squareroot(3*one);
3907 A(0,0) = A(1,1) = one/(3*one + rootthree);
3908 A(1,0) = one/rootthree;
3912 b(0) = b(1) = onehalf;
3916 c(1) = A(1,0) + A(1,1);
3921 this->
tableau_->setTVDCoeff(sspcoef);
3928 template<
class Scalar>
3935 stepper->setStepperDIRKValues(pl);
3937 if (model != Teuchos::null) {
3938 stepper->setModel(model);
3939 stepper->initialize();
3968 template<
class Scalar>
3988 std::string ICConsistency,
3989 bool ICConsistencyCheck,
3991 bool zeroInitialGuess,
3997 this->
setup(appModel, solver, useFSAL, ICConsistency, ICConsistencyCheck,
3998 useEmbedded, zeroInitialGuess, stepperRKAppAction);
4003 std::ostringstream Description;
4005 <<
"Strong Stability Preserving Diagonally-Implicit RK (stage=3, order=3)\n"
4006 <<
"SSP-Coef = 2 + 2 sqrt(2)\n"
4007 <<
"c = [ 1/( 4 + 2 sqrt(2) 1/2 (1/4)(2 + sqrt(2) ] '\n"
4008 <<
"A = [ 1/( 4 + 2 sqrt(2) ] \n"
4009 <<
" [ 1/(2 sqrt(2) 1/( 4 + 2 sqrt(2) ] \n"
4010 <<
" [ 1/(2 sqrt(2) 1/(2 sqrt(2) 1/( 4 + 2 sqrt(2) ] \n"
4011 <<
"b = [ 1/3 1/3 1/3 ] \n"
4013 return Description.str();
4023 const int NumStages = 3;
4024 const int order = 3;
4025 const Scalar sspcoef= 4.8284;
4030 const Scalar one = ST::one();
4031 const Scalar zero = ST::zero();
4032 const Scalar onethird = one/(3*one);
4033 const Scalar rootwo = ST::squareroot(2*one);
4036 A(0,0) = A(1,1) = A(2,2) = one / (4*one + 2*rootwo);
4037 A(1,0) = A(2,0) = A(2,1) = one / (2*rootwo);
4038 A(0,1) = A(0,2) = A(1,2) = zero;
4041 b(0) = b(1) = b(2) = onethird;
4045 c(1) = A(1,0) + A(1,1);
4046 c(2) = A(2,0) + A(2,1) + A(2,2);
4051 this->
tableau_->setTVDCoeff(sspcoef);
4058 template<
class Scalar>
4065 stepper->setStepperDIRKValues(pl);
4067 if (model != Teuchos::null) {
4068 stepper->setModel(model);
4069 stepper->initialize();
4097 template<
class Scalar>
4122 std::string ICConsistency,
4123 bool ICConsistencyCheck,
4125 bool zeroInitialGuess,
4131 this->
setup(appModel, solver, useFSAL, ICConsistency, ICConsistencyCheck,
4132 useEmbedded, zeroInitialGuess, stepperRKAppAction);
4137 std::ostringstream Description;
4140 <<
"Solving Ordinary Differential Equations II:\n"
4141 <<
"Stiff and Differential-Algebraic Problems,\n"
4142 <<
"2nd Revised Edition\n"
4143 <<
"E. Hairer and G. Wanner\n"
4144 <<
"Table 5.3, pg 73\n"
4148 return Description.str();
4160 const Scalar one = ST::one();
4161 const Scalar zero = ST::zero();
4169 this->
getStepperType(),A,b,c,order,order,order,emptyBStar,
false));
4176 template<
class Scalar>
4183 stepper->setStepperDIRKValues(pl);
4185 if (model != Teuchos::null) {
4186 stepper->setModel(model);
4187 stepper->initialize();
4217 template<
class Scalar>
4229 this->
setStepperName(
"RK Implicit 2 Stage 2nd order Lobatto IIIB");
4230 this->
setStepperType(
"RK Implicit 2 Stage 2nd order Lobatto IIIB");
4242 std::string ICConsistency,
4243 bool ICConsistencyCheck,
4245 bool zeroInitialGuess,
4248 this->
setStepperName(
"RK Implicit 2 Stage 2nd order Lobatto IIIB");
4249 this->
setStepperType(
"RK Implicit 2 Stage 2nd order Lobatto IIIB");
4251 this->
setup(appModel, solver, useFSAL, ICConsistency, ICConsistencyCheck,
4252 useEmbedded, zeroInitialGuess, stepperRKAppAction);
4257 std::ostringstream Description;
4260 <<
"Solving Ordinary Differential Equations II:\n"
4261 <<
"Stiff and Differential-Algebraic Problems,\n"
4262 <<
"2nd Revised Edition\n"
4263 <<
"E. Hairer and G. Wanner\n"
4264 <<
"Table 5.9, pg 76\n"
4266 <<
"A = [ 1/2 0 ]\n"
4268 <<
"b = [ 1/2 1/2 ]'";
4269 return Description.str();
4282 const Scalar zero = ST::zero();
4283 const Scalar one = ST::one();
4286 A(0,0) = as<Scalar>( one/(2*one) );
4288 A(1,0) = as<Scalar>( one/(2*one) );
4292 b(0) = as<Scalar>( one/(2*one) );
4293 b(1) = as<Scalar>( one/(2*one) );
4302 this->
getStepperType(),A,b,c,order,order,order,emptyBStar,
false));
4303 this->tableau_->setTVD(
true);
4304 this->tableau_->setTVDCoeff(2.0);
4312 template<
class Scalar>
4319 stepper->setStepperDIRKValues(pl);
4321 if (model != Teuchos::null) {
4322 stepper->setModel(model);
4323 stepper->initialize();
4356 template<
class Scalar>
4381 std::string ICConsistency,
4382 bool ICConsistencyCheck,
4384 bool zeroInitialGuess,
4390 this->
setup(appModel, solver, useFSAL, ICConsistency, ICConsistencyCheck,
4391 useEmbedded, zeroInitialGuess, stepperRKAppAction);
4396 std::ostringstream Description;
4399 <<
"Solving Ordinary Differential Equations II:\n"
4400 <<
"Stiff and Differential-Algebraic Problems,\n"
4401 <<
"2nd Revised Edition\n"
4402 <<
"E. Hairer and G. Wanner\n"
4404 <<
"c = [ 1/4 3/4 11/20 1/2 1 ]'\n"
4407 <<
" [ 17/50 -1/25 1/4 ]\n"
4408 <<
" [ 371/1360 -137/2720 15/544 1/4 ]\n"
4409 <<
" [ 25/24 -49/48 125/16 -85/12 1/4 ]\n"
4410 <<
"b = [ 25/24 -49/48 125/16 -85/12 1/4 ]'";
4412 return Description.str();
4425 const Scalar zero = ST::zero();
4426 const Scalar one = ST::one();
4427 const Scalar onequarter = as<Scalar>( one/(4*one) );
4430 A(0,0) = onequarter;
4436 A(1,0) = as<Scalar>( one / (2*one) );
4437 A(1,1) = onequarter;
4442 A(2,0) = as<Scalar>( 17*one/(50*one) );
4443 A(2,1) = as<Scalar>( -one/(25*one) );
4444 A(2,2) = onequarter;
4448 A(3,0) = as<Scalar>( 371*one/(1360*one) );
4449 A(3,1) = as<Scalar>( -137*one/(2720*one) );
4450 A(3,2) = as<Scalar>( 15*one/(544*one) );
4451 A(3,3) = onequarter;
4454 A(4,0) = as<Scalar>( 25*one/(24*one) );
4455 A(4,1) = as<Scalar>( -49*one/(48*one) );
4456 A(4,2) = as<Scalar>( 125*one/(16*one) );
4457 A(4,3) = as<Scalar>( -85*one/(12*one) );
4458 A(4,4) = onequarter;
4461 b(0) = as<Scalar>( 25*one/(24*one) );
4462 b(1) = as<Scalar>( -49*one/(48*one) );
4463 b(2) = as<Scalar>( 125*one/(16*one) );
4464 b(3) = as<Scalar>( -85*one/(12*one) );
4478 c(1) = as<Scalar>( 3*one/(4*one) );
4479 c(2) = as<Scalar>( 11*one/(20*one) );
4480 c(3) = as<Scalar>( one/(2*one) );
4493 template<
class Scalar>
4500 stepper->setStepperDIRKValues(pl);
4502 if (model != Teuchos::null) {
4503 stepper->setModel(model);
4504 stepper->initialize();
4536 template<
class Scalar>
4561 std::string ICConsistency,
4562 bool ICConsistencyCheck,
4564 bool zeroInitialGuess,
4570 this->
setup(appModel, solver, useFSAL, ICConsistency, ICConsistencyCheck,
4571 useEmbedded, zeroInitialGuess, stepperRKAppAction);
4576 std::ostringstream Description;
4579 <<
"Solving Ordinary Differential Equations II:\n"
4580 <<
"Stiff and Differential-Algebraic Problems,\n"
4581 <<
"2nd Revised Edition\n"
4582 <<
"E. Hairer and G. Wanner\n"
4584 <<
"gamma = (1/sqrt(3))*cos(pi/18)+1/2\n"
4585 <<
"delta = 1/(6*(2*gamma-1)^2)\n"
4586 <<
"c = [ gamma 1/2 1-gamma ]'\n"
4587 <<
"A = [ gamma ]\n"
4588 <<
" [ 1/2-gamma gamma ]\n"
4589 <<
" [ 2*gamma 1-4*gamma gamma ]\n"
4590 <<
"b = [ delta 1-2*delta delta ]'";
4591 return Description.str();
4604 const Scalar zero = ST::zero();
4605 const Scalar one = ST::one();
4606 const Scalar pi = as<Scalar>(4*one)*std::atan(one);
4607 const Scalar gamma = as<Scalar>( one/ST::squareroot(3*one)*std::cos(pi/(18*one))+one/(2*one) );
4608 const Scalar delta = as<Scalar>( one/(6*one*std::pow(2*gamma-one,2*one)) );
4615 A(1,0) = as<Scalar>( one/(2*one) - gamma );
4619 A(2,0) = as<Scalar>( 2*gamma );
4620 A(2,1) = as<Scalar>( one - 4*gamma );
4625 b(1) = as<Scalar>( one-2*delta );
4630 c(1) = as<Scalar>( one/(2*one) );
4631 c(2) = as<Scalar>( one - gamma );
4643 template<
class Scalar>
4650 stepper->setStepperDIRKValues(pl);
4652 if (model != Teuchos::null) {
4653 stepper->setModel(model);
4654 stepper->initialize();
4687 template<
class Scalar>
4712 std::string ICConsistency,
4713 bool ICConsistencyCheck,
4715 bool zeroInitialGuess,
4721 this->
setup(appModel, solver, useFSAL, ICConsistency, ICConsistencyCheck,
4722 useEmbedded, zeroInitialGuess, stepperRKAppAction);
4727 std::ostringstream Description;
4729 <<
"Solving Ordinary Differential Equations II:\n"
4730 <<
"Stiff and Differential-Algebraic Problems,\n"
4731 <<
"2nd Revised Edition\n"
4732 <<
"E. Hairer and G. Wanner\n"
4734 <<
"c = [ (6-sqrt(6))/10 ]\n"
4735 <<
" [ (6+9*sqrt(6))/35 ]\n"
4737 <<
" [ (4-sqrt(6))/10 ]\n"
4738 <<
" [ (4+sqrt(6))/10 ]\n"
4739 <<
"A = [ A1 A2 A3 A4 A5 ]\n"
4740 <<
" A1 = [ (6-sqrt(6))/10 ]\n"
4741 <<
" [ (-6+5*sqrt(6))/14 ]\n"
4742 <<
" [ (888+607*sqrt(6))/2850 ]\n"
4743 <<
" [ (3153-3082*sqrt(6))/14250 ]\n"
4744 <<
" [ (-32583+14638*sqrt(6))/71250 ]\n"
4746 <<
" [ (6-sqrt(6))/10 ]\n"
4747 <<
" [ (126-161*sqrt(6))/1425 ]\n"
4748 <<
" [ (3213+1148*sqrt(6))/28500 ]\n"
4749 <<
" [ (-17199+364*sqrt(6))/142500 ]\n"
4752 <<
" [ (6-sqrt(6))/10 ]\n"
4753 <<
" [ (-267+88*sqrt(6))/500 ]\n"
4754 <<
" [ (1329-544*sqrt(6))/2500 ]\n"
4758 <<
" [ (6-sqrt(6))/10 ]\n"
4759 <<
" [ (-96+131*sqrt(6))/625 ]\n"
4764 <<
" [ (6-sqrt(6))/10 ]\n"
4768 <<
" [ (16-sqrt(6))/36 ]\n"
4769 <<
" [ (16+sqrt(6))/36 ]'";
4770 return Description.str();
4783 const Scalar zero = ST::zero();
4784 const Scalar one = ST::one();
4785 const Scalar sqrt6 = ST::squareroot(as<Scalar>(6*one));
4786 const Scalar gamma = as<Scalar>( (6*one - sqrt6) / (10*one) );
4795 A(1,0) = as<Scalar>( (-6*one+5*one*sqrt6)/(14*one) );
4801 A(2,0) = as<Scalar>( (888*one+607*one*sqrt6)/(2850*one) );
4802 A(2,1) = as<Scalar>( (126*one-161*one*sqrt6)/(1425*one) );
4807 A(3,0) = as<Scalar>( (3153*one-3082*one*sqrt6)/(14250*one) );
4808 A(3,1) = as<Scalar>( (3213*one+1148*one*sqrt6)/(28500*one) );
4809 A(3,2) = as<Scalar>( (-267*one+88*one*sqrt6)/(500*one) );
4813 A(4,0) = as<Scalar>( (-32583*one+14638*one*sqrt6)/(71250*one) );
4814 A(4,1) = as<Scalar>( (-17199*one+364*one*sqrt6)/(142500*one) );
4815 A(4,2) = as<Scalar>( (1329*one-544*one*sqrt6)/(2500*one) );
4816 A(4,3) = as<Scalar>( (-96*one+131*sqrt6)/(625*one) );
4822 b(2) = as<Scalar>( one/(9*one) );
4823 b(3) = as<Scalar>( (16*one-sqrt6)/(36*one) );
4824 b(4) = as<Scalar>( (16*one+sqrt6)/(36*one) );
4828 c(1) = as<Scalar>( (6*one+9*one*sqrt6)/(35*one) );
4830 c(3) = as<Scalar>( (4*one-sqrt6)/(10*one) );
4831 c(4) = as<Scalar>( (4*one+sqrt6)/(10*one) );
4843 template<
class Scalar>
4850 stepper->setStepperDIRKValues(pl);
4852 if (model != Teuchos::null) {
4853 stepper->setModel(model);
4854 stepper->initialize();
4880 template<
class Scalar>
4905 std::string ICConsistency,
4906 bool ICConsistencyCheck,
4908 bool zeroInitialGuess,
4914 this->
setup(appModel, solver, useFSAL, ICConsistency, ICConsistencyCheck,
4915 useEmbedded, zeroInitialGuess, stepperRKAppAction);
4920 std::ostringstream Description;
4925 <<
"b = [ 1/2 1/2 ]'\n"
4926 <<
"bstar = [ 1 0 ]'";
4927 return Description.str();
4942 const Scalar one = ST::one();
4943 const Scalar zero = ST::zero();
4946 A(0,0) = one; A(0,1) = zero;
4947 A(1,0) = -one; A(1,1) = one;
4950 b(0) = as<Scalar>(one/(2*one));
4951 b(1) = as<Scalar>(one/(2*one));
4970 template<
class Scalar>
4977 stepper->setStepperDIRKValues(pl);
4979 if (model != Teuchos::null) {
4980 stepper->setModel(model);
4981 stepper->initialize();
5021 template<
class Scalar>
5046 std::string ICConsistency,
5047 bool ICConsistencyCheck,
5049 bool zeroInitialGuess,
5061 this->
setTableau(A,b,c,order,orderMin,orderMax,bstar);
5064 this->
tableau_->isImplicit() !=
true, std::logic_error,
5065 "Error - General DIRK did not receive a DIRK Butcher Tableau!\n");
5067 this->
setup(appModel, solver, useFSAL, ICConsistency, ICConsistencyCheck,
5068 useEmbedded, zeroInitialGuess, stepperRKAppAction);
5073 std::stringstream Description;
5075 <<
"The format of the Butcher Tableau parameter list is\n"
5076 <<
" <Parameter name=\"A\" type=\"string\" value=\"# # # ;\n"
5079 <<
" <Parameter name=\"b\" type=\"string\" value=\"# # #\"/>\n"
5080 <<
" <Parameter name=\"c\" type=\"string\" value=\"# # #\"/>\n\n"
5081 <<
"Note the number of stages is implicit in the number of entries.\n"
5082 <<
"The number of stages must be consistent.\n"
5084 <<
"Default tableau is 'SDIRK 2 Stage 2nd order':\n"
5085 <<
" Computer Methods for ODEs and DAEs\n"
5086 <<
" U. M. Ascher and L. R. Petzold\n"
5088 <<
" gamma = (2-sqrt(2))/2\n"
5089 <<
" c = [ gamma 1 ]'\n"
5090 <<
" A = [ gamma 0 ]\n"
5091 <<
" [ 1-gamma gamma ]\n"
5092 <<
" b = [ 1-gamma gamma ]'";
5093 return Description.str();
5098 if (this->
tableau_ == Teuchos::null) {
5101 auto t = stepper->getTableau();
5104 t->A(),t->b(),t->c(),
5105 t->order(),t->orderMin(),t->orderMax(),
5138 std::ostringstream Astream;
5139 Astream.precision(15);
5140 for(
int i = 0; i < A.
numRows(); i++) {
5141 for(
int j = 0; j < A.
numCols()-1; j++) {
5142 Astream << A(i,j) <<
" ";
5144 Astream << A(i,A.
numCols()-1);
5145 if ( i != A.
numRows()-1 ) Astream <<
"; ";
5147 tableauPL->
set<std::string>(
"A", Astream.str());
5149 std::ostringstream bstream;
5150 bstream.precision(15);
5151 for(
int i = 0; i < b.
length()-1; i++) {
5152 bstream << b(i) <<
" ";
5154 bstream << b(b.
length()-1);
5155 tableauPL->
set<std::string>(
"b", bstream.str());
5157 std::ostringstream cstream;
5158 cstream.precision(15);
5159 for(
int i = 0; i < c.
length()-1; i++) {
5160 cstream << c(i) <<
" ";
5162 cstream << c(c.
length()-1);
5163 tableauPL->
set<std::string>(
"c", cstream.str());
5165 tableauPL->
set<
int>(
"order", this->
tableau_->order());
5167 if ( bstar.
length() == 0 ) {
5168 tableauPL->
set(
"bstar",
"");
5170 std::ostringstream bstarstream;
5171 bstarstream.precision(15);
5172 for(
int i = 0; i < bstar.
length()-1; i++) {
5173 bstarstream << bstar(i) <<
" ";
5175 bstarstream << bstar(bstar.
length()-1);
5176 tableauPL->
set<std::string>(
"bstar", bstarstream.str());
5179 pl->
set(
"Tableau", *tableauPL);
5188 template<
class Scalar>
5195 stepper->setStepperDIRKValues(pl);
5197 if (pl != Teuchos::null) {
5199 auto t = stepper->createTableau(pl);
5200 stepper->setTableau( t->A(),t->b(),t->c(),
5201 t->order(),t->orderMin(),t->orderMax(),
5206 stepper->getTableau()->isDIRK() !=
true, std::logic_error,
5207 "Error - General DIRK did not receive a DIRK Butcher Tableau!\n");
5209 if (model != Teuchos::null) {
5210 stepper->setModel(model);
5211 stepper->initialize();
5221 #endif // Tempus_StepperRKButcherTableau_hpp
std::string getDescription() const
StepperERK_3Stage3rdOrderTVD(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &appModel, bool useFSAL, std::string ICConsistency, bool ICConsistencyCheck, bool useEmbedded, const Teuchos::RCP< StepperRKAppAction< Scalar > > &stepperRKAppAction)
StepperERK_Merson45(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &appModel, bool useFSAL, std::string ICConsistency, bool ICConsistencyCheck, bool useEmbedded, const Teuchos::RCP< StepperRKAppAction< Scalar > > &stepperRKAppAction)
Teuchos::RCP< StepperEDIRK_TrapezoidalRule< Scalar > > createStepperEDIRK_TrapezoidalRule(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &model, Teuchos::RCP< Teuchos::ParameterList > pl)
Nonmember constructor - ModelEvaluator and ParameterList.
General Implicit Runge-Kutta Butcher Tableau.
StepperERK_3Stage3rdOrder()
Default constructor.
StepperEDIRK_TrapezoidalRule(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &appModel, const Teuchos::RCP< Thyra::NonlinearSolverBase< Scalar > > &solver, bool useFSAL, std::string ICConsistency, bool ICConsistencyCheck, bool useEmbedded, bool zeroInitialGuess, const Teuchos::RCP< StepperRKAppAction< Scalar > > &stepperRKAppAction)
StepperERK_3Stage3rdOrderTVD()
Default constructor.
StepperERK_3Stage3rdOrderHeun(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &appModel, bool useFSAL, std::string ICConsistency, bool ICConsistencyCheck, bool useEmbedded, const Teuchos::RCP< StepperRKAppAction< Scalar > > &stepperRKAppAction)
RK Implicit 2 Stage 2nd order Lobatto IIIB.
virtual std::string getDefaultICConsistency() const
std::string getStepperType() const
Get the stepper type. The stepper type is used as an identifier for the stepper, and can only be set ...
StepperERK_Ralston()
Default constructor.
Teuchos::RCP< StepperSDIRK_SSPDIRK33< Scalar > > createStepperSDIRK_SSPDIRK33(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &model, Teuchos::RCP< Teuchos::ParameterList > pl)
Nonmember constructor - ModelEvaluator and ParameterList.
std::string getDescription() const
Explicit Runge-Kutta time stepper.
virtual void setupDefault()
Default setup for constructor.
Strong Stability Preserving Diagonally-Implicit RK Butcher Tableau.
Teuchos::RCP< Teuchos::ParameterList > getValidParametersBasicDIRK() const
StepperSDIRK_SSPDIRK33(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &appModel, const Teuchos::RCP< Thyra::NonlinearSolverBase< Scalar > > &solver, bool useFSAL, std::string ICConsistency, bool ICConsistencyCheck, bool useEmbedded, bool zeroInitialGuess, const Teuchos::RCP< StepperRKAppAction< Scalar > > &stepperRKAppAction)
StepperDIRK_BackwardEuler(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &appModel, const Teuchos::RCP< Thyra::NonlinearSolverBase< Scalar > > &solver, bool useFSAL, std::string ICConsistency, bool ICConsistencyCheck, bool useEmbedded, bool zeroInitialGuess, const Teuchos::RCP< StepperRKAppAction< Scalar > > &stepperRKAppAction)
StepperSDIRK_3Stage2ndOrder()
Default constructor.
Teuchos::RCP< StepperEDIRK_2Stage3rdOrder< Scalar > > createStepperEDIRK_2Stage3rdOrder(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &model, Teuchos::RCP< Teuchos::ParameterList > pl)
Nonmember constructor - ModelEvaluator and ParameterList.
StepperERK_SSPERK54(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &appModel, bool useFSAL, std::string ICConsistency, bool ICConsistencyCheck, bool useEmbedded, const Teuchos::RCP< StepperRKAppAction< Scalar > > &stepperRKAppAction)
Teuchos::RCP< StepperSDIRK_5Stage4thOrder< Scalar > > createStepperSDIRK_5Stage4thOrder(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &model, Teuchos::RCP< Teuchos::ParameterList > pl)
Nonmember constructor - ModelEvaluator and ParameterList.
StepperSDIRK_ImplicitMidpoint()
Default constructor.
StepperEDIRK_2StageTheta(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &appModel, const Teuchos::RCP< Thyra::NonlinearSolverBase< Scalar > > &solver, bool useFSAL, std::string ICConsistency, bool ICConsistencyCheck, bool useEmbedded, bool zeroInitialGuess, const Teuchos::RCP< StepperRKAppAction< Scalar > > &stepperRKAppAction, Scalar theta=Scalar(0.5))
Backward Euler Runge-Kutta Butcher Tableau.
StepperSDIRK_3Stage4thOrder(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &appModel, const Teuchos::RCP< Thyra::NonlinearSolverBase< Scalar > > &solver, bool useFSAL, std::string ICConsistency, bool ICConsistencyCheck, bool useEmbedded, bool zeroInitialGuess, const Teuchos::RCP< StepperRKAppAction< Scalar > > &stepperRKAppAction)
StepperERK_Midpoint()
Default constructor.
Teuchos::RCP< RKButcherTableau< Scalar > > tableau_
T & get(const std::string &name, T def_value)
StepperSDIRK_2Stage2ndOrder()
Default constructor.
Teuchos::RCP< StepperERK_SSPERK54< Scalar > > createStepperERK_SSPERK54(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &model, Teuchos::RCP< Teuchos::ParameterList > pl)
Nonmember constructor - ModelEvaluator and ParameterList.
StepperSDIRK_ImplicitMidpoint(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &appModel, const Teuchos::RCP< Thyra::NonlinearSolverBase< Scalar > > &solver, bool useFSAL, std::string ICConsistency, bool ICConsistencyCheck, bool useEmbedded, bool zeroInitialGuess, const Teuchos::RCP< StepperRKAppAction< Scalar > > &stepperRKAppAction)
General Explicit Runge-Kutta Butcher Tableau.
std::string getDescription() const
StepperDIRK_General(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &appModel, const Teuchos::RCP< Thyra::NonlinearSolverBase< Scalar > > &solver, bool useFSAL, std::string ICConsistency, bool ICConsistencyCheck, bool useEmbedded, bool zeroInitialGuess, const Teuchos::RCP< StepperRKAppAction< Scalar > > &stepperRKAppAction, const Teuchos::SerialDenseMatrix< int, Scalar > &A, const Teuchos::SerialDenseVector< int, Scalar > &b, const Teuchos::SerialDenseVector< int, Scalar > &c, const int order, const int orderMin, const int orderMax, const Teuchos::SerialDenseVector< int, Scalar > &bstar)
ParameterList & set(std::string const &name, T const &value, std::string const &docString="", RCP< const ParameterEntryValidator > const &validator=null)
Teuchos::RCP< StepperERK_3Stage3rdOrderHeun< Scalar > > createStepperERK_3Stage3rdOrderHeun(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &model, Teuchos::RCP< Teuchos::ParameterList > pl)
Nonmember constructor - ModelEvaluator and ParameterList.
RK Explicit 5 Stage 3rd order by Kinnmark and Gray.
Teuchos::RCP< StepperERK_Merson45< Scalar > > createStepperERK_Merson45(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &model, Teuchos::RCP< Teuchos::ParameterList > pl)
Nonmember constructor - ModelEvaluator and ParameterList.
Teuchos::RCP< StepperEDIRK_2StageTheta< Scalar > > createStepperEDIRK_2StageTheta(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &model, Teuchos::RCP< Teuchos::ParameterList > pl)
Nonmember constructor - ModelEvaluator and ParameterList.
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
StepperSDIRK_SSPDIRK23(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &appModel, const Teuchos::RCP< Thyra::NonlinearSolverBase< Scalar > > &solver, bool useFSAL, std::string ICConsistency, bool ICConsistencyCheck, bool useEmbedded, bool zeroInitialGuess, const Teuchos::RCP< StepperRKAppAction< Scalar > > &stepperRKAppAction)
std::string getDescription() const
StepperSDIRK_3Stage2ndOrder(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &appModel, const Teuchos::RCP< Thyra::NonlinearSolverBase< Scalar > > &solver, bool useFSAL, std::string ICConsistency, bool ICConsistencyCheck, bool useEmbedded, bool zeroInitialGuess, const Teuchos::RCP< StepperRKAppAction< Scalar > > &stepperRKAppAction)
void setStepperName(std::string s)
Set the stepper name.
StepperDIRK_1StageTheta()
Default constructor.
std::string getDescription() const
StepperEDIRK_2Stage3rdOrder()
Default constructor.
Strong Stability Preserving Explicit RK Butcher Tableau.
RK Explicit 4 Stage 3rd order by Runge.
std::string getDescription() const
StepperERK_5Stage3rdOrderKandG(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &appModel, bool useFSAL, std::string ICConsistency, bool ICConsistencyCheck, bool useEmbedded, const Teuchos::RCP< StepperRKAppAction< Scalar > > &stepperRKAppAction)
std::string gammaTypeDefault_
std::string getDescription() const
StepperEDIRK_2Stage3rdOrder(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &appModel, const Teuchos::RCP< Thyra::NonlinearSolverBase< Scalar > > &solver, bool useFSAL, std::string ICConsistency, bool ICConsistencyCheck, bool useEmbedded, bool zeroInitialGuess, const Teuchos::RCP< StepperRKAppAction< Scalar > > &stepperRKAppAction)
std::string getDescription() const
Teuchos::RCP< StepperERK_3_8Rule< Scalar > > createStepperERK_3_8Rule(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &model, Teuchos::RCP< Teuchos::ParameterList > pl)
Nonmember constructor - ModelEvaluator and ParameterList.
Teuchos::RCP< StepperERK_BogackiShampine32< Scalar > > createStepperERK_BogackiShampine32(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &model, Teuchos::RCP< Teuchos::ParameterList > pl)
Nonmember constructor - ModelEvaluator and ParameterList.
Teuchos::RCP< StepperSDIRK_5Stage5thOrder< Scalar > > createStepperSDIRK_5Stage5thOrder(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &model, Teuchos::RCP< Teuchos::ParameterList > pl)
Nonmember constructor - ModelEvaluator and ParameterList.
Teuchos::RCP< StepperERK_3Stage3rdOrder< Scalar > > createStepperERK_3Stage3rdOrder(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &model, Teuchos::RCP< Teuchos::ParameterList > pl)
Nonmember constructor - ModelEvaluator and ParameterList.
StepperDIRK_2Stage2ndOrderLobattoIIIB(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &appModel, const Teuchos::RCP< Thyra::NonlinearSolverBase< Scalar > > &solver, bool useFSAL, std::string ICConsistency, bool ICConsistencyCheck, bool useEmbedded, bool zeroInitialGuess, const Teuchos::RCP< StepperRKAppAction< Scalar > > &stepperRKAppAction)
StepperERK_5Stage3rdOrderKandG()
Default constructor.
RK Trapezoidal Rule (A.K.A. RK Crank-Nicolson)
Teuchos::RCP< StepperDIRK_2Stage2ndOrderLobattoIIIB< Scalar > > createStepperDIRK_2Stage2ndOrderLobattoIIIB(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &model, Teuchos::RCP< Teuchos::ParameterList > pl)
Nonmember constructor - ModelEvaluator and ParameterList.
Teuchos::RCP< StepperERK_ForwardEuler< Scalar > > createStepperERK_ForwardEuler(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &model, Teuchos::RCP< Teuchos::ParameterList > pl)
Nonmember constructor - ModelEvaluator and ParameterList.
Teuchos::RCP< const Teuchos::ParameterList > getValidParameters() const
Teuchos::RCP< StepperSDIRK_21Pair< Scalar > > createStepperSDIRK_21Pair(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &model, Teuchos::RCP< Teuchos::ParameterList > pl)
Nonmember constructor - ModelEvaluator and ParameterList.
StepperERK_BogackiShampine32()
Default constructor.
bool isInitialized_
True if stepper's member data is initialized.
std::string getDescription() const
std::string getDescription() const
Teuchos::RCP< StepperSDIRK_SSPDIRK23< Scalar > > createStepperSDIRK_SSPDIRK23(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &model, Teuchos::RCP< Teuchos::ParameterList > pl)
Nonmember constructor - ModelEvaluator and ParameterList.
StepperERK_4Stage3rdOrderRunge(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &appModel, bool useFSAL, std::string ICConsistency, bool ICConsistencyCheck, bool useEmbedded, const Teuchos::RCP< StepperRKAppAction< Scalar > > &stepperRKAppAction)
Teuchos::RCP< StepperERK_3Stage3rdOrderTVD< Scalar > > createStepperERK_3Stage3rdOrderTVD(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &model, Teuchos::RCP< Teuchos::ParameterList > pl)
Nonmember constructor - ModelEvaluator and ParameterList.
void setGamma(Scalar gamma)
Teuchos::RCP< StepperSDIRK_3Stage4thOrder< Scalar > > createStepperSDIRK_3Stage4thOrder(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &model, Teuchos::RCP< Teuchos::ParameterList > pl)
Nonmember constructor - ModelEvaluator and ParameterList.
Explicit RK 3/8th Rule Butcher Tableau.
std::string getDescription() const
StepperSDIRK_2Stage3rdOrder(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &appModel, const Teuchos::RCP< Thyra::NonlinearSolverBase< Scalar > > &solver, bool useFSAL, std::string ICConsistency, bool ICConsistencyCheck, bool useEmbedded, bool zeroInitialGuess, const Teuchos::RCP< StepperRKAppAction< Scalar > > &stepperRKAppAction, std::string gammaType="3rd Order A-stable", Scalar gamma=Scalar(0.7886751345948128))
bool isParameter(const std::string &name) const
StepperERK_Midpoint(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &appModel, bool useFSAL, std::string ICConsistency, bool ICConsistencyCheck, bool useEmbedded, const Teuchos::RCP< StepperRKAppAction< Scalar > > &stepperRKAppAction)
void setTheta(Scalar theta)
Diagonally Implicit Runge-Kutta (DIRK) time stepper.
void setTheta(Scalar theta)
std::string getDescription() const
StepperERK_4Stage4thOrder(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &appModel, bool useFSAL, std::string ICConsistency, bool ICConsistencyCheck, bool useEmbedded, const Teuchos::RCP< StepperRKAppAction< Scalar > > &stepperRKAppAction)
StepperERK_Trapezoidal()
Default constructor.
Strong Stability Preserving Diagonally-Implicit RK Butcher Tableau.
StepperERK_Merson45()
Default constructor.
std::string getDescription() const
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
StepperSDIRK_5Stage4thOrder()
Default constructor.
Application Action for StepperRKBase.
Teuchos::RCP< StepperDIRK_1Stage1stOrderRadauIA< Scalar > > createStepperDIRK_1Stage1stOrderRadauIA(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &model, Teuchos::RCP< Teuchos::ParameterList > pl)
Nonmember constructor - ModelEvaluator and ParameterList.
bool useFSAL_
Use First-Same-As-Last (FSAL) principle.
std::string getDescription() const
virtual void setupDefault()
Default setup for constructor.
Explicit RK Bogacki-Shampine Butcher Tableau.
StepperSDIRK_5Stage5thOrder(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &appModel, const Teuchos::RCP< Thyra::NonlinearSolverBase< Scalar > > &solver, bool useFSAL, std::string ICConsistency, bool ICConsistencyCheck, bool useEmbedded, bool zeroInitialGuess, const Teuchos::RCP< StepperRKAppAction< Scalar > > &stepperRKAppAction)
EDIRK 2 Stage Theta Method.
std::string getDescription() const
void setICConsistencyCheck(bool c)
std::string getDescription() const
Teuchos::RCP< StepperSDIRK_3Stage2ndOrder< Scalar > > createStepperSDIRK_3Stage2ndOrder(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &model, Teuchos::RCP< Teuchos::ParameterList > pl)
Nonmember constructor - ModelEvaluator and ParameterList.
Teuchos::RCP< StepperERK_5Stage3rdOrderKandG< Scalar > > createStepperERK_5Stage3rdOrderKandG(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &model, Teuchos::RCP< Teuchos::ParameterList > pl)
Nonmember constructor - ModelEvaluator and ParameterList.
OrdinalType length() const
Teuchos::RCP< const Teuchos::ParameterList > getValidParameters() const
StepperSDIRK_21Pair(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &appModel, const Teuchos::RCP< Thyra::NonlinearSolverBase< Scalar > > &solver, bool useFSAL, std::string ICConsistency, bool ICConsistencyCheck, bool useEmbedded, bool zeroInitialGuess, const Teuchos::RCP< StepperRKAppAction< Scalar > > &stepperRKAppAction)
std::string getDescription() const
Teuchos::RCP< StepperDIRK_BackwardEuler< Scalar > > createStepperDIRK_BackwardEuler(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &model, Teuchos::RCP< Teuchos::ParameterList > pl)
Nonmember constructor - ModelEvaluator and ParameterList.
StepperERK_Trapezoidal(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &appModel, bool useFSAL, std::string ICConsistency, bool ICConsistencyCheck, bool useEmbedded, const Teuchos::RCP< StepperRKAppAction< Scalar > > &stepperRKAppAction)
StepperDIRK_1Stage1stOrderRadauIA(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &appModel, const Teuchos::RCP< Thyra::NonlinearSolverBase< Scalar > > &solver, bool useFSAL, std::string ICConsistency, bool ICConsistencyCheck, bool useEmbedded, bool zeroInitialGuess, const Teuchos::RCP< StepperRKAppAction< Scalar > > &stepperRKAppAction)
std::string getDescription() const
std::string getGammaType() const
Strong Stability Preserving Diagonally-Implicit RK Butcher Tableau.
std::string getDescription() const
Teuchos::RCP< StepperSDIRK_2Stage3rdOrder< Scalar > > createStepperSDIRK_2Stage3rdOrder(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &model, Teuchos::RCP< Teuchos::ParameterList > pl)
Nonmember constructor - ModelEvaluator and ParameterList.
Teuchos::RCP< Teuchos::ParameterList > getValidParametersBasicERK() const
std::string getDescription() const
StepperSDIRK_SSPDIRK22(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &appModel, const Teuchos::RCP< Thyra::NonlinearSolverBase< Scalar > > &solver, bool useFSAL, std::string ICConsistency, bool ICConsistencyCheck, bool useEmbedded, bool zeroInitialGuess, const Teuchos::RCP< StepperRKAppAction< Scalar > > &stepperRKAppAction)
StepperSDIRK_21Pair()
Default constructor.
std::string getDescription() const
StepperSDIRK_2Stage3rdOrder()
Default constructor.
std::string getDescription() const
std::string getDescription() const
StepperDIRK_1Stage1stOrderRadauIA()
Default constructor.
void setTableau(const Teuchos::SerialDenseMatrix< int, Scalar > &A, const Teuchos::SerialDenseVector< int, Scalar > &b, const Teuchos::SerialDenseVector< int, Scalar > &c, const int order, const int orderMin, const int orderMax, const Teuchos::SerialDenseVector< int, Scalar > &bstar=Teuchos::SerialDenseVector< int, Scalar >())
Teuchos::RCP< StepperERK_Midpoint< Scalar > > createStepperERK_Midpoint(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &model, Teuchos::RCP< Teuchos::ParameterList > pl)
Nonmember constructor - ModelEvaluator and ParameterList.
RK Explicit 3 Stage 3rd order TVD.
Teuchos::RCP< const Teuchos::ParameterList > getValidParameters() const
std::string getDescription() const
OrdinalType numCols() const
StepperEDIRK_2StageTheta()
Default constructor.
Teuchos::RCP< StepperSDIRK_SSPDIRK22< Scalar > > createStepperSDIRK_SSPDIRK22(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &model, Teuchos::RCP< Teuchos::ParameterList > pl)
Nonmember constructor - ModelEvaluator and ParameterList.
Teuchos::RCP< const Teuchos::ParameterList > getValidParameters() const
StepperERK_ForwardEuler(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &appModel, bool useFSAL, std::string ICConsistency, bool ICConsistencyCheck, bool useEmbedded, const Teuchos::RCP< StepperRKAppAction< Scalar > > &stepperRKAppAction)
std::string getDescription() const
std::string getDescription() const
StepperDIRK_General()
Default constructor.
Teuchos::RCP< StepperERK_Trapezoidal< Scalar > > createStepperERK_Trapezoidal(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &model, Teuchos::RCP< Teuchos::ParameterList > pl)
Nonmember constructor - ModelEvaluator and ParameterList.
Forward Euler Runge-Kutta Butcher Tableau.
StepperDIRK_1StageTheta(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &appModel, const Teuchos::RCP< Thyra::NonlinearSolverBase< Scalar > > &solver, bool useFSAL, std::string ICConsistency, bool ICConsistencyCheck, bool useEmbedded, bool zeroInitialGuess, const Teuchos::RCP< StepperRKAppAction< Scalar > > &stepperRKAppAction, Scalar theta=Scalar(0.5))
std::string getDescription() const
StepperERK_BogackiShampine32(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &appModel, bool useFSAL, std::string ICConsistency, bool ICConsistencyCheck, bool useEmbedded, const Teuchos::RCP< StepperRKAppAction< Scalar > > &stepperRKAppAction)
Teuchos::RCP< StepperSDIRK_ImplicitMidpoint< Scalar > > createStepperSDIRK_ImplicitMidpoint(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &model, Teuchos::RCP< Teuchos::ParameterList > pl)
Nonmember constructor - ModelEvaluator and ParameterList.
StepperERK_3Stage3rdOrder(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &appModel, bool useFSAL, std::string ICConsistency, bool ICConsistencyCheck, bool useEmbedded, const Teuchos::RCP< StepperRKAppAction< Scalar > > &stepperRKAppAction)
Teuchos::RCP< StepperERK_General< Scalar > > createStepperERK_General(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &model, Teuchos::RCP< Teuchos::ParameterList > pl)
Nonmember constructor - ModelEvaluator and ParameterList.
StepperDIRK_BackwardEuler()
Default constructor.
virtual std::string getDescription() const
TypeTo as(const TypeFrom &t)
void setGammaType(std::string gammaType)
StepperERK_3Stage3rdOrderHeun()
Default constructor.
Teuchos::RCP< StepperSDIRK_2Stage2ndOrder< Scalar > > createStepperSDIRK_2Stage2ndOrder(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &model, Teuchos::RCP< Teuchos::ParameterList > pl)
Nonmember constructor - ModelEvaluator and ParameterList.
Explicit RK Merson Butcher Tableau.
std::string getDescription() const
Runge-Kutta 4th order Butcher Tableau.
std::string getDescription() const
Teuchos::RCP< const Teuchos::ParameterList > getValidParameters() const
RK Explicit 3 Stage 3rd order.
StepperERK_4Stage4thOrder()
Default constructor.
void setGamma(Scalar gamma)
std::string getDescription() const
StepperSDIRK_5Stage5thOrder()
Default constructor.
void setTableau(const Teuchos::SerialDenseMatrix< int, Scalar > &A, const Teuchos::SerialDenseVector< int, Scalar > &b, const Teuchos::SerialDenseVector< int, Scalar > &c, const int order, const int orderMin, const int orderMax, const Teuchos::SerialDenseVector< int, Scalar > &bstar=Teuchos::SerialDenseVector< int, Scalar >())
RK Implicit 1 Stage 1st order Radau IA.
StepperERK_Ralston(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &appModel, bool useFSAL, std::string ICConsistency, bool ICConsistencyCheck, bool useEmbedded, const Teuchos::RCP< StepperRKAppAction< Scalar > > &stepperRKAppAction)
StepperSDIRK_3Stage4thOrder()
Default constructor.
virtual void setup(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &wrapperModel, const Teuchos::RCP< Thyra::NonlinearSolverBase< Scalar > > &solver, bool useFSAL, std::string ICConsistency, bool ICConsistencyCheck, bool useEmbedded, bool zeroInitialGuess, const Teuchos::RCP< StepperRKAppAction< Scalar > > &stepperRKAppAction)
Setup for constructor.
StepperSDIRK_5Stage4thOrder(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &appModel, const Teuchos::RCP< Thyra::NonlinearSolverBase< Scalar > > &solver, bool useFSAL, std::string ICConsistency, bool ICConsistencyCheck, bool useEmbedded, bool zeroInitialGuess, const Teuchos::RCP< StepperRKAppAction< Scalar > > &stepperRKAppAction)
Teuchos::RCP< StepperERK_Ralston< Scalar > > createStepperERK_Ralston(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &model, Teuchos::RCP< Teuchos::ParameterList > pl)
Nonmember constructor - ModelEvaluator and ParameterList.
StepperDIRK_2Stage2ndOrderLobattoIIIB()
Default constructor.
std::string getDescription() const
virtual void setUseFSAL(bool a)
StepperSDIRK_2Stage2ndOrder(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &appModel, const Teuchos::RCP< Thyra::NonlinearSolverBase< Scalar > > &solver, bool useFSAL, std::string ICConsistency, bool ICConsistencyCheck, bool useEmbedded, bool zeroInitialGuess, const Teuchos::RCP< StepperRKAppAction< Scalar > > &stepperRKAppAction, Scalar gamma=Scalar(0.2928932188134524))
StepperERK_ForwardEuler()
Default constructor.
Teuchos::RCP< const Teuchos::ParameterList > getValidParameters() const
StepperSDIRK_SSPDIRK32(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &appModel, const Teuchos::RCP< Thyra::NonlinearSolverBase< Scalar > > &solver, bool useFSAL, std::string ICConsistency, bool ICConsistencyCheck, bool useEmbedded, bool zeroInitialGuess, const Teuchos::RCP< StepperRKAppAction< Scalar > > &stepperRKAppAction)
StepperERK_3_8Rule(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &appModel, bool useFSAL, std::string ICConsistency, bool ICConsistencyCheck, bool useEmbedded, const Teuchos::RCP< StepperRKAppAction< Scalar > > &stepperRKAppAction)