russell_ode 1.13.0

use russell_lab::math::SQRT_6;

/// Default number of steps to use when the automatic stepping is not available
pub const N_EQUAL_STEPS: usize = 10;

// References:
//
// 1. Hairer E, Nørsett, SP, Wanner G (2008) Solving Ordinary Differential Equations I.
//    Non-stiff Problems. Second Revised Edition. Corrected 3rd printing 2008. Springer Series
//    in Computational Mathematics, 528p
// 2. Hairer E, Wanner G (2002) Solving Ordinary Differential Equations II.
//    Stiff and Differential-Algebraic Problems. Second Revised Edition.
//    Corrected 2nd printing 2002. Springer Series in Computational Mathematics, 614p
// 3. Kreyszig, E (2011) Advanced engineering mathematics; in collaboration with Kreyszig H,
//    Edward JN 10th ed 2011, Hoboken, New Jersey, Wiley

// Runge-Kutta -- order 2 ---------------------------------------------------------------------

// Table 1.1 on page 135 of Ref#1 (also known as mid-point)

#[rustfmt::skip]
pub(crate) const RUNGE_KUTTA_2_A: [[f64; 2]; 2] = [
    [0.0, 0.0],
    [1.0 / 2.0, 0.0],
];

#[rustfmt::skip]
pub(crate) const RUNGE_KUTTA_2_B: [f64; 2] = [0.0, 1.0];

#[rustfmt::skip]
pub(crate) const RUNGE_KUTTA_2_C: [f64; 2] = [0.0, 1.0 / 2.0];

// Runge-Kutta -- order 3 ---------------------------------------------------------------------

// Table 1.1 on page 135 of Ref#1 (Note: this method has 4 stages)

#[rustfmt::skip]
pub(crate) const RUNGE_KUTTA_3_A: [[f64; 4]; 4] = [
    [0.0, 0.0, 0.0, 0.0],
    [1.0 / 2.0, 0.0, 0.0, 0.0],
    [0.0, 1.0, 0.0, 0.0],
    [0.0, 0.0, 1.0, 0.0],
];

#[rustfmt::skip]
pub(crate) const RUNGE_KUTTA_3_B: [f64; 4] = [1.0 / 6.0, 2.0 / 3.0, 0.0, 1.0 / 6.0];

#[rustfmt::skip]
pub(crate) const RUNGE_KUTTA_3_C: [f64; 4] = [0.0, 1.0 / 2.0, 1.0, 1.0];

// Heun -- order 3 ----------------------------------------------------------------------------

// Table 1.1 on page 135 of Ref#1

#[rustfmt::skip]
pub(crate) const HEUN_3_A: [[f64; 3]; 3] = [
    [0.0, 0.0, 0.0],
    [1.0 / 3.0, 0.0, 0.0],
    [0.0, 2.0 / 3.0, 0.0],
];

#[rustfmt::skip]
pub(crate) const HEUN_3_B: [f64; 3] = [1.0 / 4.0, 0.0, 3.0 / 4.0];

#[rustfmt::skip]
pub(crate) const HEUN_3_C: [f64; 3] = [0.0, 1.0 / 3.0, 2.0 / 3.0];

// Runge-Kutta -- order 4 ---------------------------------------------------------------------

// Table 1.2 on page 138 of Ref#1 ("The" Runge-Kutta method)

#[rustfmt::skip]
pub(crate) const RUNGE_KUTTA_4_A: [[f64; 4]; 4] = [
    [0.0, 0.0, 0.0, 0.0],
    [1.0 / 2.0, 0.0, 0.0, 0.0],
    [0.0, 1.0 / 2.0, 0.0, 0.0],
    [0.0, 0.0, 1.0, 0.0],
];

#[rustfmt::skip]
pub(crate) const RUNGE_KUTTA_4_B: [f64; 4] = [1.0 / 6.0, 2.0 / 6.0, 2.0 / 6.0, 1.0 / 6.0];

#[rustfmt::skip]
pub(crate) const RUNGE_KUTTA_4_C: [f64; 4] = [0.0, 1.0 / 2.0, 1.0 / 2.0, 1.0];

// Runge-Kutta -- alternative 3/8 rule -- order 4 ---------------------------------------------

// Table 1.2 on page 128 of Ref#1

#[rustfmt::skip]
pub(crate) const RUNGE_KUTTA_ALT_4_A: [[f64; 4]; 4] = [
    [0.0, 0.0, 0.0, 0.0],
    [1.0 / 3.0, 0.0, 0.0, 0.0],
    [-1.0 / 3.0, 1.0, 0.0, 0.0],
    [1.0, -1.0, 1.0, 0.0],
];

#[rustfmt::skip]
pub(crate) const RUNGE_KUTTA_ALT_4_B: [f64; 4] = [1.0 / 8.0, 3.0 / 8.0, 3.0 / 8.0, 1.0 / 8.0];

#[rustfmt::skip]
pub(crate) const RUNGE_KUTTA_ALT_4_C: [f64; 4] = [0.0, 1.0 / 3.0, 2.0 / 3.0, 1.0];

// Modified Euler -- order 2 -- embedded 2(1) -------------------------------------------------

// Table 21.1 on page 903 of Ref#3 (also known as Improved Euler method)

#[rustfmt::skip]
pub(crate) const MODIFIED_EULER_A: [[f64; 2]; 2] = [
    [0.0, 0.0],
    [1.0, 0.0],
];

#[rustfmt::skip]
pub(crate) const MODIFIED_EULER_B: [f64; 2] = [1.0 / 2.0, 1.0 / 2.0];

#[allow(unused)]
#[rustfmt::skip]
pub(crate) const MODIFIED_EULER_BE: [f64; 2] = [1.0, 0.0];

#[rustfmt::skip]
pub(crate) const MODIFIED_EULER_C: [f64; 2] = [0.0, 1.0];

#[rustfmt::skip]
pub(crate) const MODIFIED_EULER_E: [f64; 2] = [-1.0 / 2.0, 1.0 / 2.0];

// Merson -- order 4 -- embedded 4("5") or 4(3) -----------------------------------------------

// Table 4.1 on page 167 of Ref#1
// Note: This is the Merson 4("5") method where "5" means that the order 5 is for linear equations
// with constant coefficients; otherwise the method is of order 3.

#[rustfmt::skip]
pub(crate) const MERSON_4_A: [[f64; 5]; 5] = [
    [0.0, 0.0, 0.0, 0.0, 0.0],
    [1.0 / 3.0, 0.0, 0.0, 0.0, 0.0],
    [1.0 / 6.0, 1.0 / 6.0, 0.0, 0.0, 0.0],
    [1.0 / 8.0, 0.0, 3.0 / 8.0, 0.0, 0.0],
    [1.0 / 2.0, 0.0, -3.0 / 2.0, 2.0, 0.0],
];

#[rustfmt::skip]
pub(crate) const MERSON_4_B: [f64; 5] = [1.0 / 6.0, 0.0, 0.0, 2.0 / 3.0, 1.0 / 6.0];

#[allow(unused)]
#[rustfmt::skip]
pub(crate) const MERSON_4_BE: [f64; 5] = [1.0 / 10.0, 0.0, 3.0 / 10.0, 2.0 / 5.0, 1.0 / 5.0];

#[rustfmt::skip]
pub(crate) const MERSON_4_C: [f64; 5] = [0.0, 1.0 / 3.0, 1.0 / 3.0, 1.0 / 2.0, 1.0];

#[rustfmt::skip]
pub(crate) const MERSON_4_E: [f64; 5] = [1.0 / 15.0, 0.0, -3.0 / 10.0, 4.0 / 15.0, -1.0 / 30.0];

// Zonneveld -- order 4 -- embedded 4(3) ------------------------------------------------------

// Table 4.2 on page 167 of Ref#1

#[rustfmt::skip]
pub(crate) const ZONNEVELD_4_A: [[f64; 5]; 5] = [
    [0.0, 0.0, 0.0, 0.0, 0.0],
    [1.0 / 2.0, 0.0, 0.0, 0.0, 0.0],
    [0.0, 1.0 / 2.0, 0.0, 0.0, 0.0],
    [0.0, 0.0, 1.0, 0.0, 0.0],
    [5.0 / 32.0, 7.0 / 32.0, 13.0 / 32.0, -1.0 / 32.0, 0.0],
];

#[rustfmt::skip]
pub(crate) const ZONNEVELD_4_B: [f64; 5] = [1.0 / 6.0, 1.0 / 3.0, 1.0 / 3.0, 1.0 / 6.0, 0.0];

#[allow(unused)]
#[rustfmt::skip]
pub(crate) const ZONNEVELD_4_BE: [f64; 5] = [-1.0 / 2.0, 7.0 / 3.0, 7.0 / 3.0, 13.0 / 6.0, -16.0 / 3.0];

#[rustfmt::skip]
pub(crate) const ZONNEVELD_4_C: [f64; 5] = [0.0, 1.0 / 2.0, 1.0 / 2.0, 1.0, 3.0 / 4.0];

#[rustfmt::skip]
pub(crate) const ZONNEVELD_4_E: [f64; 5] = [2.0 / 3.0, -2.0, -2.0, -2.0, 16.0 / 3.0];

// Fehlberg -- order 4 -- embedded 4(5) -------------------------------------------------------

// Table 5.1 on page 177 of Ref#1

#[rustfmt::skip]
pub(crate) const FEHLBERG_4_A: [[f64; 6]; 6] = [
    [0.0, 0.0, 0.0, 0.0, 0.0, 0.0],
    [1.0 / 4.0, 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],
];

#[rustfmt::skip]
pub(crate) const FEHLBERG_4_B: [f64; 6] = [25.0 / 216.0, 0.0, 1408.0 / 2565.0, 2197.0 / 4104.0, -1.0 / 5.0, 0.0];

#[allow(unused)]
#[rustfmt::skip]
pub(crate) const FEHLBERG_4_BE: [f64; 6] = [16.0 / 135.0, 0.0, 6656.0 / 12825.0, 28561.0 / 56430.0, -9.0 / 50.0, 2.0 / 55.0];

#[rustfmt::skip]
pub(crate) const FEHLBERG_4_C: [f64; 6] = [0.0, 1.0 / 4.0, 3.0 / 8.0, 12.0 / 13.0, 1.0, 1.0 / 2.0];

#[rustfmt::skip]
pub(crate) const FEHLBERG_4_E: [f64; 6] = [-1.0 / 360.0, 0.0, 128.0 / 4275.0, 2197.0 / 75240.0, -1.0 / 50.0, -2.0 / 55.0];

// Dormand-Prince -- order 5 -- embedded 5(4) FSAL --------------------------------------------

// Table 5.2 on page 178 of Ref#1

#[rustfmt::skip]
pub(crate) const DORMAND_PRINCE_5_A: [[f64; 7]; 7] = [
    [0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0],
    [1.0 / 5.0, 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],
];

#[rustfmt::skip]
pub(crate) const DORMAND_PRINCE_5_B: [f64; 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];

#[allow(unused)]
#[rustfmt::skip]
pub(crate) const DORMAND_PRINCE_5_BE: [f64; 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];

#[rustfmt::skip]
pub(crate) const DORMAND_PRINCE_5_C: [f64; 7] = [0.0, 1.0 / 5.0, 3.0 / 10.0, 4.0 / 5.0, 8.0 / 9.0, 1.0, 1.0];

#[rustfmt::skip]
pub(crate) const DORMAND_PRINCE_5_E: [f64; 7] = [71.0 / 57600.0, 0.0, -71.0 / 16695.0, 71.0 / 1920.0, -17253.0 / 339200.0, 22.0 / 525.0, -1.0 / 40.0];

#[rustfmt::skip]
pub(crate) const DORMAND_PRINCE_5_D: [[f64; 7]; 1] = [[
    -12715105075.0 / 11282082432.0,  // D1
    0.00000000000000000000000000,    // D2
    87487479700.0 / 32700410799.0,   // D3
    -10690763975.0 / 1880347072.0,   // D4
    701980252875.0 / 199316789632.0, // D5
    -1453857185.0 / 822651844.0,     // D6
    69997945.0 / 29380423.0,         // D7
]];

// Verner -- order 6 -- embedded -- 6(5) ------------------------------------------------------

// Table 5.4 on page 181 of Ref#1
// Hairer-Wanner: "[...] Fehlberg's methods suffer from the fact that they give identically zero
// error estimates for quadrature problems y' = f(x). The first high order embedded formulas
// which avoid this drawback were constructed by Verner (1978). [...]"

#[rustfmt::skip]
pub(crate) const VERNER_6_A: [[f64; 8]; 8] = [
    [0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0],
    [1.0 / 6.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0],
    [4.0 / 75.0, 16.0 / 75.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0],
    [5.0 / 6.0, -8.0 / 3.0, 5.0 / 2.0, 0.0, 0.0, 0.0, 0.0, 0.0],
    [-165.0 / 64.0, 55.0 / 6.0, -425.0 / 64.0, 85.0 / 96.0, 0.0, 0.0, 0.0, 0.0],
    [12.0 / 5.0, -8.0, 4015.0 / 612.0, -11.0 / 36.0, 88.0 / 255.0, 0.0, 0.0, 0.0],
    [-8263.0 / 15000.0, 124.0 / 75.0, -643.0 / 680.0, -81.0 / 250.0, 2484.0 / 10625.0, 0.0, 0.0, 0.0],
    [3501.0 / 1720.0, -300.0 / 43.0, 297275.0 / 52632.0, -319.0 / 2322.0, 24068.0 / 84065.0, 0.0, 3850.0 / 26703.0, 0.0],
];

#[rustfmt::skip]
pub(crate) const VERNER_6_B: [f64; 8] = [3.0 / 40.0, 0.0, 875.0 / 2244.0, 23.0 / 72.0, 264.0 / 1955.0, 0.0, 125.0 / 11592.0, 43.0 / 616.0];

#[allow(unused)]
#[rustfmt::skip]
pub(crate) const VERNER_6_BE: [f64; 8] = [13.0 / 160.0, 0.0, 2375.0 / 5984.0, 5.0 / 16.0, 12.0 / 85.0, 3.0 / 44.0, 0.0, 0.0];

#[rustfmt::skip]
pub(crate) const VERNER_6_C: [f64; 8] = [0.0, 1.0 / 6.0, 4.0 / 15.0, 2.0 / 3.0, 5.0 / 6.0, 1.0, 1.0 / 15.0, 1.0];

#[rustfmt::skip]
pub(crate) const VERNER_6_E: [f64; 8] = [-1.0 / 160.0, 0.0, -125.0 / 17952.0, 1.0 / 144.0, -12.0 / 1955.0, -3.0 / 44.0, 125.0 / 11592.0, 43.0 / 616.0];

// Fehlberg -- order 7 -- embedded -- 7(8) ----------------------------------------------------

// Table 5.3 on page 180 of Ref#1

#[rustfmt::skip]
pub(crate) const FEHLBERG_7_A: [[f64; 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],
];

#[rustfmt::skip]
pub(crate) const FEHLBERG_7_B: [f64; 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];

#[allow(unused)]
#[rustfmt::skip]
pub(crate) const FEHLBERG_7_BE: [f64; 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];

#[rustfmt::skip]
pub(crate) const FEHLBERG_7_C: [f64; 13] = [0.0, 2.0 / 27.0, 1.0 / 9.0, 1.0 / 6.0, 5.0 / 12.0, 1.0 / 2.0, 5.0 / 6.0, 1.0 / 6.0, 2.0 / 3.0, 1.0 / 3.0, 1.0, 0.0, 1.0];

#[rustfmt::skip]
pub(crate) const FEHLBERG_7_E: [f64; 13] = [41.0 / 840.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 41.0 / 840.0, -41.0 / 840.0, -41.0 / 840.0];

// Dormand-Prince -- order 8 -- embedded -- 8(5,3) --------------------------------------------

// The coefficients here are taken from dop853.f, available on Hairer's website

#[rustfmt::skip]
pub(crate) const DORMAND_PRINCE_8_A: [[f64; 12]; 12] = [
    [0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0],
    [5.26001519587677318785587544488e-2, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0],
    [1.97250569845378994544595329183e-2, 5.91751709536136983633785987549e-2, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0],
    [2.95875854768068491816892993775e-2, 0.0, 8.87627564304205475450678981324e-2, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0],
    [2.41365134159266685502369798665e-1, 0.0, -8.84549479328286085344864962717e-1, 9.24834003261792003115737966543e-1, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0],
    [3.7037037037037037037037037037e-2, 0.0, 0.0, 1.70828608729473871279604482173e-1, 1.25467687566822425016691814123e-1, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0],
    [3.71093750000000000000000000000e-2, 0.0, 0.0, 1.70252211019544039314978060272e-1, 6.02165389804559606850219397283e-2, -1.75781250000000000000000000000e-2, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0],
    [3.70920001185047927108779319836e-2, 0.0, 0.0, 1.70383925712239993810214054705e-1, 1.07262030446373284651809199168e-1, -1.53194377486244017527936158236e-2, 8.27378916381402288758473766002e-3, 0.0, 0.0, 0.0, 0.0, 0.0],
    [6.24110958716075717114429577812e-1, 0.0, 0.0, -3.36089262944694129406857109825, -8.68219346841726006818189891453e-1, 2.75920996994467083049415600797e1, 2.01540675504778934086186788979e1, -4.34898841810699588477366255144e1, 0.0, 0.0, 0.0, 0.0],
    [4.77662536438264365890433908527e-1, 0.0, 0.0, -2.48811461997166764192642586468, -5.90290826836842996371446475743e-1, 2.12300514481811942347288949897e1, 1.52792336328824235832596922938e1, -3.32882109689848629194453265587e1, -2.03312017085086261358222928593e-2, 0.0, 0.0, 0.0],
    [-9.3714243008598732571704021658e-1, 0.0, 0.0, 5.18637242884406370830023853209, 1.09143734899672957818500254654, -8.14978701074692612513997267357, -1.85200656599969598641566180701e1, 2.27394870993505042818970056734e1, 2.49360555267965238987089396762, -3.0467644718982195003823669022, 0.0, 0.0],
    [2.27331014751653820792359768449, 0.0, 0.0, -1.05344954667372501984066689879e1, -2.00087205822486249909675718444, -1.79589318631187989172765950534e1, 2.79488845294199600508499808837e1, -2.85899827713502369474065508674, -8.87285693353062954433549289258, 1.23605671757943030647266201528e1, 6.43392746015763530355970484046e-1, 0.0],
];

#[rustfmt::skip]
pub(crate) const DORMAND_PRINCE_8_B: [f64; 12] = [
    5.42937341165687622380535766363e-2,
    0.0,
    0.0,
    0.0,
    0.0,
    4.45031289275240888144113950566,
    1.89151789931450038304281599044,
    -5.8012039600105847814672114227,
    3.1116436695781989440891606237e-1,
    -1.52160949662516078556178806805e-1,
    2.01365400804030348374776537501e-1,
    4.47106157277725905176885569043e-2,
];

#[rustfmt::skip]
pub(crate) const DORMAND_PRINCE_8_C: [f64; 12] = [
    0.0,
    2.0 * (6.0 - SQRT_6) / 135.0,
    (6.0 - SQRT_6) / 45.0,
    (6.0 - SQRT_6) / 30.0,
    (6.0 + SQRT_6) / 30.0,
    1.0 / 3.0,
    1.0 / 4.0,
    4.0 / 13.0,
    127.0 / 195.0,
    3.0 / 5.0,
    6.0 / 7.0,
    1.0,
];

#[rustfmt::skip]
pub(crate) const DORMAND_PRINCE_8_E: [f64; 12] = [
    0.1312004499419488073250102996e-01,
    0.0,
    0.0,
    0.0,
    0.0,
    -0.1225156446376204440720569753e+01,
    -0.4957589496572501915214079952e+00,
    0.1664377182454986536961530415e+01,
    -0.3503288487499736816886487290e+00,
    0.3341791187130174790297318841e+00,
    0.8192320648511571246570742613e-01,
    -0.2235530786388629525884427845e-01,
];

#[rustfmt::skip]
pub(crate) const DORMAND_PRINCE_8_D: [[f64; 16]; 4] = [
    [
        -0.84289382761090128651353491142e+01 , // d41
        0.0                                  , // d42
        0.0                                  , // d43
        0.0                                  , // d44
        0.0                                  , // d45
        0.56671495351937776962531783590e+00  , // d46
        -0.30689499459498916912797304727e+01 , // d47
        0.23846676565120698287728149680e+01  , // d48
        0.21170345824450282767155149946e+01  , // d49
        -0.87139158377797299206789907490e+00 , // d410
        0.22404374302607882758541771650e+01  , // d411
        0.63157877876946881815570249290e+00  , // d412
        -0.88990336451333310820698117400e-01 , // d413
        0.18148505520854727256656404962e+02  , // d414
        -0.91946323924783554000451984436e+01 , // d415
        -0.44360363875948939664310572000e+01 , // d416
    ],
    [
        0.10427508642579134603413151009e+02  , // d51
        0.0                                  , // d52
        0.0                                  , // d53
        0.0                                  , // d54
        0.0                                  , // d55
        0.24228349177525818288430175319e+03  , // d56
        0.16520045171727028198505394887e+03  , // d57
        -0.37454675472269020279518312152e+03 , // d58
        -0.22113666853125306036270938578e+02 , // d59
        0.77334326684722638389603898808e+01  , // d510
        -0.30674084731089398182061213626e+02 , // d511
        -0.93321305264302278729567221706e+01 , // d512
        0.15697238121770843886131091075e+02  , // d513
        -0.31139403219565177677282850411e+02 , // d514
        -0.93529243588444783865713862664e+01 , // d515
        0.35816841486394083752465898540e+02  , // d516
    ],
    [
        0.19985053242002433820987653617e+02  , // d61
        0.0                                  , // d62
        0.0                                  , // d63
        0.0                                  , // d64
        0.0                                  , // d65
        -0.38703730874935176555105901742e+03 , // d66
        -0.18917813819516756882830838328e+03 , // d67
        0.52780815920542364900561016686e+03  , // d68
        -0.11573902539959630126141871134e+02 , // d69
        0.68812326946963000169666922661e+01  , // d610
        -0.10006050966910838403183860980e+01 , // d611
        0.77771377980534432092869265740e+00  , // d612
        -0.27782057523535084065932004339e+01 , // d613
        -0.60196695231264120758267380846e+02 , // d614
        0.84320405506677161018159903784e+02  , // d615
        0.11992291136182789328035130030e+02  , // d616
    ],
    [
        -0.25693933462703749003312586129e+02 , // d71
        0.0                                  , // d72
        0.0                                  , // d73
        0.0                                  , // d74
        0.0                                  , // d75
        -0.15418974869023643374053993627e+03 , // d76
        -0.23152937917604549567536039109e+03 , // d77
        0.35763911791061412378285349910e+03  , // d78
        0.93405324183624310003907691704e+02  , // d79
        -0.37458323136451633156875139351e+02 , // d710
        0.10409964950896230045147246184e+03  , // d711
        0.29840293426660503123344363579e+02  , // d712
        -0.43533456590011143754432175058e+02 , // d713
        0.96324553959188282948394950600e+02  , // d714
        -0.39177261675615439165231486172e+02 , // d715
        -0.14972683625798562581422125276e+03 , // d716
    ],
];

#[rustfmt::skip]
pub(crate) const DORMAND_PRINCE_8_AD: [[f64; 15]; 3] = [
    [
        5.61675022830479523392909219681e-2  , // 14,1
        0.0                                 , // 14.2
        0.0                                 , // 14.3
        0.0                                 , // 14.4
        0.0                                 , // 14.5
        0.0                                 , // 14.6
        2.53500210216624811088794765333e-1  , // 14,7
        -2.46239037470802489917441475441e-1 , // 14,8
        -1.24191423263816360469010140626e-1 , // 14,9
        1.5329179827876569731206322685e-1   , // 14,10
        8.20105229563468988491666602057e-3  , // 14,11
        7.56789766054569976138603589584e-3  , // 14,12
        -8.2980000000000000000000000000e-3  , // 14,13
        0.0                                 , // 14.14
        0.0                                 , // 14.15
    ],
    [
        3.18346481635021405060768473261e-2  , // 15,1
        0.0                                 , // 15.2
        0.0                                 , // 15.3
        0.0                                 , // 15.4
        0.0                                 , // 15.5
        2.83009096723667755288322961402e-2  , // 15,6
        5.35419883074385676223797384372e-2  , // 15,7
        -5.49237485713909884646569340306e-2 , // 15,8
        0.0                                 , // 15.9
        0.0                                 , // 15.10
        -1.08347328697249322858509316994e-4 , // 15,11
        3.82571090835658412954920192323e-4  , // 15,12
        -3.40465008687404560802977114492e-4 , // 15,13
        1.41312443674632500278074618366e-1  , // 15,14
        0.0                                 , // 16.15
    ],
    [
        -4.28896301583791923408573538692e-1 , // 16,1
        0.0                                 , // 16.2
        0.0                                 , // 16.3
        0.0                                 , // 16.4
        0.0                                 , // 16.5
        -4.69762141536116384314449447206    , // 16,6
        7.68342119606259904184240953878     , // 16,7
        4.06898981839711007970213554331     , // 16,8
        3.56727187455281109270669543021e-1  , // 16,9
        0.0                                 , // 16.10
        0.0                                 , // 16.11
        0.0                                 , // 16.12
        -1.39902416515901462129418009734e-3 , // 16,13
        2.9475147891527723389556272149      , // 16,14
        -9.15095847217987001081870187138    , // 16,15
    ],
];

#[rustfmt::skip]
pub(crate) const DORMAND_PRINCE_8_CD: [f64; 3] = [
    0.1,       // c14
    0.2,       // c15
    7.0 / 9.0, // c16
];

pub(crate) const DORMAND_PRINCE_8_BHH1: f64 = 0.244094488188976377952755905512e+00;
pub(crate) const DORMAND_PRINCE_8_BHH2: f64 = 0.733846688281611857341361741547e+00;
pub(crate) const DORMAND_PRINCE_8_BHH3: f64 = 0.220588235294117647058823529412e-01;

////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////

#[cfg(test)]
mod tests {
    use super::*;
    use russell_lab::approx_eq;

    #[test]
    #[rustfmt::skip]
    fn ee_constants_are_correct() {
        // E = B - Be
        for i in 0..MODIFIED_EULER_B  .len() { approx_eq(MODIFIED_EULER_E   [i], MODIFIED_EULER_B    [i] - MODIFIED_EULER_BE  [i], 1e-15); }
        for i in 0..MERSON_4_B        .len() { approx_eq(MERSON_4_E         [i], MERSON_4_B          [i] - MERSON_4_BE        [i], 1e-15); }
        for i in 0..ZONNEVELD_4_B     .len() { approx_eq(ZONNEVELD_4_E      [i], ZONNEVELD_4_B       [i] - ZONNEVELD_4_BE     [i], 1e-15); }
        for i in 0..FEHLBERG_4_B      .len() { approx_eq(FEHLBERG_4_E       [i], FEHLBERG_4_B        [i] - FEHLBERG_4_BE      [i], 1e-15); }
        for i in 0..DORMAND_PRINCE_5_B.len() { approx_eq(DORMAND_PRINCE_5_E [i], DORMAND_PRINCE_5_B  [i] - DORMAND_PRINCE_5_BE[i], 1e-15); }
        for i in 0..VERNER_6_B        .len() { approx_eq(VERNER_6_E         [i], VERNER_6_B          [i] - VERNER_6_BE        [i], 1e-15); }
        for i in 0..FEHLBERG_7_B      .len() { approx_eq(FEHLBERG_7_E       [i], FEHLBERG_7_B        [i] - FEHLBERG_7_BE      [i], 1e-15); }
    }
}