File butcher_tableau.hpp¶
File List > astrea > astro > astro > propagation > numerical > butcher_tableau.hpp
Go to the documentation of this file
#pragma once
#include <cstddef>
namespace astrea {
namespace astro {
namespace RK45 { // Traditional Runge-Kutta 5(4) 6 stage method
static const std::size_t nStages = 6;
static const double c[6] = { 0.0, 0.2, 0.3, 0.6, 1.0, 0.875 };
static const double b[6] = { 37.0 / 378.0, 0.0, 250.0 / 621.0, 125.0 / 594.0, 0.0, 512.0 / 1771.0 };
static const double bhat[6] = { 2825.0 / 27648.0, 0.0, 18575.0 / 48384.0, 13525.0 / 55296.0, 277.0 / 14336.0, 0.25 };
static const double a[6][6] = { { 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 },
{ 0.2, 0.0, 0.0, 0.0, 0.0, 0.0 },
{ 3.0 / 40.0, 9.0 / 40.0, 0.0, 0.0, 0.0, 0.0 },
{ 0.3, -0.9, 1.2, 0.0, 0.0, 0.0 },
{ -11.0 / 54.0, 2.5, -70.0 / 27.0, 35.0 / 27.0, 0.0, 0.0 },
{ 1631.0 / 55296.0, 175.0 / 512.0, 575.0 / 13824.0, 44275.0 / 110592.0, 253.0 / 4096.0, 0.0 } };
} // namespace RK45
namespace RKF45 { // Runge-Kutta-Fehlberg 5(4) 6 stage method
static const std::size_t nStages = 6;
static const double c[6] = { 0.0, 0.25, 3.0 / 8.0, 12.0 / 13.0, 1.0, 0.5 };
static const double b[6] = { 25.0 / 216.0, 0.0, 1408.0 / 2565.0, 2197.0 / 4104.0, -0.2, 0.0 };
static const double bhat[6] = { 16.0 / 135.0, 0.0, 6656.0 / 12825.0, 28561.0 / 56430.0, -9.0 / 50.0, 2.0 / 55.0 };
static const double a[6][6] = { { 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 },
{ 0.25, 0.0, 0.0, 0.0, 0.0, 0.0 },
{ 3.0 / 32.0, 9.0 / 32.0, 0.0, 0.0, 0.0, 0.0 },
{ 1932.0 / 2197.0, -7200.0 / 2197.0, 7296.0 / 2197.0, 0.0, 0.0, 0.0 },
{ 439.0 / 216.0, -8.0, 3680.0 / 513.0, -845.0 / 4104.0, 0.0, 0.0 },
{ -8.0 / 27.0, 2.0, -3544.0 / 2565.0, 1859.0 / 4104.0, -11.0 / 40.0, 0.0 } };
} // namespace RKF45
namespace DOP45 { // Dormand-Prince 5(4) 7 stage method
static const std::size_t nStages = 7;
static const double c[7] = { 0.0, 0.2, 0.3, 0.8, 8.0 / 9.0, 1.0, 1.0 };
static const double b[7] = { 5179.0 / 57600.0, 0.0, 7571.0 / 16695.0, 393.0 / 640.0, -92097.0 / 339200.0,
187.0 / 2100.0, 1.0 / 40.0 };
static const double bhat[7] = { 35.0 / 384.0, 0.0, 500.0 / 1113.0, 125.0 / 192.0, -2187.0 / 6784.0, 11.0 / 84.0, 0.0 };
static const double a[7][7] = { { 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 },
{ 0.2, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 },
{ 3.0 / 40.0, 9.0 / 40.0, 0.0, 0.0, 0.0, 0.0, 0.0 },
{ 44.0 / 45.0, -56.0 / 15.0, 32.0 / 9.0, 0.0, 0.0, 0.0, 0.0 },
{ 19372.0 / 6561.0, -25360.0 / 2187.0, 64448.0 / 6561.0, -212.0 / 729.0, 0.0, 0.0, 0.0 },
{ 9017.0 / 3168.0, -355.0 / 33.0, 46732.0 / 5247.0, 49.0 / 176.0, -5103.0 / 18656.0, 0.0, 0.0 },
{ 35.0 / 384.0, 0.0, 500.0 / 1113.0, 125.0 / 192.0, -2187.0 / 6784.0, 11.0 / 84.0, 0.0 } };
} // namespace DOP45
namespace RKF78 { // Runge-Kutta-Fehlberg 7(8) 13 stage method (this one is so fucking fast)
static const std::size_t nStages = 13;
static const double c[13] = { 0.0, 2.0 / 27.0, 1.0 / 9.0, 1.0 / 6.0, 5.0 / 12.0, 0.5, 5.0 / 6.0,
1.0 / 6.0, 2.0 / 3.0, 1.0 / 3.0, 1.0, 0.0, 1.0 };
static const double b[13] = { 41.0 / 840.0, 0.0, 0.0, 0.0, 0.0, 34.0 / 105.0, 9.0 / 35.0,
9.0 / 35.0, 9.0 / 280.0, 9.0 / 280.0, 41.0 / 840.0, 0.0, 0.0 };
static const double bhat[13] = { 0.0, 0.0, 0.0, 0.0, 0.0, 34.0 / 105.0, 9.0 / 35.0,
9.0 / 35.0, 9.0 / 280.0, 9.0 / 280.0, 0.0, 41.0 / 840.0, 41.0 / 840.0 };
static const double a[13][13] = {
{ 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 },
{ 2.0 / 27.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 },
{ 1.0 / 36.0, 1.0 / 12.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 },
{ 1.0 / 24.0, 0.0, 1.0 / 8.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 },
{ 5.0 / 12.0, 0.0, -25.0 / 16.0, 25.0 / 16.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 },
{ 1.0 / 20.0, 0.0, 0.0, 1.0 / 4.0, 1.0 / 5.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 },
{ -25.0 / 108.0, 0.0, 0.0, 125.0 / 108.0, -65.0 / 27.0, 125.0 / 54.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 },
{ 31.0 / 300.0, 0.0, 0.0, 0.0, 61.0 / 225.0, -2.0 / 9.0, 13.0 / 900.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 },
{ 2.0, 0.0, 0.0, -53.0 / 6.0, 704.0 / 45.0, -107.0 / 9.0, 67.0 / 90.0, 3.0, 0.0, 0.0, 0.0, 0.0, 0.0 },
{ -91.0 / 108.0, 0.0, 0.0, 23.0 / 108.0, -976.0 / 135.0, 311.0 / 54.0, -19.0 / 60.0, 17.0 / 6.0, -1.0 / 12.0, 0.0, 0.0, 0.0, 0.0 },
{ 2383.0 / 4100.0, 0.0, 0.0, -341.0 / 164.0, 4496.0 / 1025.0, -301.0 / 82.0, 2133.0 / 4100.0, 45.0 / 82.0, 45.0 / 164.0, 18.0 / 41.0, 0.0, 0.0, 0.0 },
{ 3.0 / 205.0, 0.0, 0.0, 0.0, 0.0, -6.0 / 41.0, -3.0 / 205.0, -3.0 / 41.0, 3.0 / 41.0, 6.0 / 41.0, 0.0, 0.0, 0.0 },
{ -1777.0 / 4100.0, 0.0, 0.0, -341.0 / 164.0, 4496.0 / 1025.0, -289.0 / 82.0, 2193.0 / 4100.0, 51.0 / 82.0, 33.0 / 164.0, 12.0 / 41.0, 0.0, 1.0, 0.0 }
};
} // namespace RKF78
namespace DOP78 { // Dormand-Prince 7(8) 13 stage method
static const std::size_t nStages = 13;
static const double c[13] = { 0.0,
1.0 / 18.0,
1.0 / 12.0,
1.0 / 8.0,
5.0 / 16.0,
3.0 / 8.0,
59.0 / 400.0,
93.0 / 200.0,
5490023248.0 / 9719169821.0,
13.0 / 20.0,
1201146811.0 / 1299019798.0,
1.0,
1.0 };
static const double b[13] = { 13451932.0 / 455176623.0,
0.0,
0.0,
0.0,
0.0,
-808719846.0 / 976000145.0,
1757004468.0 / 5645159321.0,
656045339.0 / 265891186.0,
-3867574721.0 / 1518517206.0,
465885868.0 / 322736535.0,
53011238.0 / 667516719.0,
2.0 / 45.0,
0.0 };
static const double bhat[13] = { 14005451.0 / 335480064.0,
0.0,
0.0,
0.0,
0.0,
-59238493.0 / 1068277825.0,
181606767.0 / 758867731.0,
561292985.0 / 797845732.0,
-1041891430.0 / 1371343529.0,
760417239.0 / 1151165299.0,
118820643.0 / 751138087.0,
-528747749.0 / 2220607170.0,
1.0 / 4.0 };
static const double a[13][13] = {
{ 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 },
{ 1.0 / 18.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 },
{ 1.0 / 48.0, 1.0 / 16.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 },
{ 1.0 / 32.0, 0.0, 3.0 / 32.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 },
{ 5.0 / 16.0, 0.0, -75.0 / 64.0, 75.0 / 64.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 },
{ 3.0 / 80.0, 0.0, 0.0, 3.0 / 16.0, 3.0 / 20.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 },
{ 29443841.0 / 614563906.0, 0.0, 0.0, 77736538.0 / 692538347.0, -28693883.0 / 1125000000.0, 23124283.0 / 1800000000.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 },
{ 16016141.0 / 946692911.0, 0.0, 0.0, 61564180.0 / 158732637.0, 22789713.0 / 633445777.0, 545815736.0 / 2771057229.0, -180193667.0 / 1043307555.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 },
{ 39632708.0 / 573591083.0,
0.0,
0.0,
-433636366.0 / 683701615.0,
-421739975.0 / 2616292301.0,
100302831.0 / 723423059.0,
790204164.0 / 839813087.0,
800635310.0 / 3783071287.0,
0.0,
0.0,
0.0,
0.0,
0.0 },
{ 246121993.0 / 1340847787.0,
0.0,
0.0,
-37695042795.0 / 15268766246.0,
-309121744.0 / 1061227803.0,
-12992083.0 / 490766935.0,
6005943493.0 / 2108947869.0,
393006217.0 / 1396673457.0,
123872331.0 / 1001029789.0,
0.0,
0.0,
0.0,
0.0 },
{ -1028468189.0 / 846180014.0,
0.0,
0.0,
8478235783.0 / 508512852.0,
1311729495.0 / 1432422823.0,
-10304129995.0 / 1701304382.0,
-48777925059.0 / 3047939560.0,
15336726248.0 / 1032824649.0,
-45442868181.0 / 3398467696.0,
3065993473.0 / 597172653.0,
0.0,
0.0,
0.0 },
{ 185892177.0 / 718116043.0,
0.0,
0.0,
-3185094517.0 / 667107341.0,
-477755414.0 / 1098053517.0,
-703635378.0 / 230739211.0,
5731566787.0 / 1027545527.0,
5232866602.0 / 850066563.0,
-4093664535.0 / 808688257.0,
3962137247.0 / 1805957418.0,
65686358.0 / 487910083.0,
0.0,
0.0 },
{ 403863854.0 / 491063109.0,
0.0,
0.0,
-5068492393.0 / 434740067.0,
-411421997.0 / 543043805.0,
652783627.0 / 914296604.0,
11173962825.0 / 925320556.0,
-13158990841.0 / 6184727034.0,
3936647629.0 / 1978049680.0,
-160528059.0 / 685178525.0,
248638103.0 / 1413531060.0,
0.0,
0.0 }
};
} // namespace DOP78
} // namespace astro
} // namespace astrea