Algorithms

List of all algorithms (skeleton of each numerical method to solve)

Runge-Kutta methods

A Runge-Kutta method is defined in ponio by its Butcher tableau

\[\begin{split}\begin{array}{c|c} c & A \\ \hline & b^\top \end{array}\end{split}\]

which is stored as a JSON file in database folder.

Runge-Kutta methods are useful to solve a problem like

\[\dot{u} = f(t, u)\]

The method, with the previous Butcher tableau, reads as

\[\begin{split}\begin{aligned} u^{(i)} &= u^n + \Delta t \sum_j a_{ij} k_j, \quad i = 1, \dots, s \\ k_i &= f(t^n + c_i\Delta t, u^{(i)}) \\ u^{n+1} &= u^n + \Delta t \sum_i b_i k_i \end{aligned}\end{split}\]

Explicit methods

When the matrix \(A\) is strictly lower triangular, the Runge-Kutta method is called explicit. The only think you need to provide to solve a problem with this kind of method is the function \(f\). See below for the list of explicit Runge-Kutta methods in ponio.

template<typename value_t>
using ponio::runge_kutta::euler_t = explicit_runge_kutta::explicit_runge_kutta<butcher_euler<value_t>>

Euler method.

see more on ponio

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|c} 0 & 0 \\ \hline & 1 \end{array} \end{split}\]

• stages: 1

• order: 1

• stages order: 1

• stability function:

  \[ z + 1 \]

Template Parameters:

value_t – type of coefficient (doubleby default)

template<typename value_t>
using ponio::runge_kutta::explicit_euler_sub4_t = explicit_runge_kutta::explicit_runge_kutta<butcher_explicit_euler_sub4<value_t>>

Explicit Euler sub4 method.

see more on ponio

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|cccc} 0 & 0 & 0 & 0 & 0 \\ \frac{1}{4} & \frac{1}{4} & 0 & 0 & 0 \\ \frac{1}{2} & \frac{1}{4} & \frac{1}{4} & 0 & 0 \\ \frac{3}{4} & \frac{1}{4} & \frac{1}{4} & \frac{1}{4} & 0 \\ \hline & \frac{1}{4} & \frac{1}{4} & \frac{1}{4} & \frac{1}{4} \end{array} \end{split}\]

• stages: 4

• order: 1

• stages order: 1

• stability function:

  \[ \frac{z^{4}}{256} + \frac{z^{3}}{16} + \frac{3 z^{2}}{8} + z + 1 \]

Template Parameters:

value_t – type of coefficient (doubleby default)

template<typename value_t>
using ponio::runge_kutta::rk56_t = explicit_runge_kutta::explicit_runge_kutta<butcher_rk56<value_t>>

RK56 method.

see more on ponio

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|cccccc} 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ \frac{2}{23} & \frac{2}{23} & 0 & 0 & 0 & 0 & 0 \\ \frac{12}{37} & - \frac{384}{1369} & \frac{828}{1369} & 0 & 0 & 0 & 0 \\ \frac{27}{29} & \frac{5631861}{622340} & - \frac{8039673}{622340} & \frac{24}{5} & 0 & 0 & 0 \\ \frac{199}{200} & \frac{22890428394764947}{1641109248000000} & - \frac{78986676649487}{3964032000000} & \frac{103911467638313}{14784768000000} & - \frac{308153608007}{5544288000000} & 0 & 0 \\ 1 & \frac{1781059432255}{124415232228} & - \frac{12785194207}{625202172} & \frac{904736654489}{125792366742} & - \frac{246740990}{4701689307} & - \frac{1472000000}{264184008767} & 0 \\ \hline & \frac{75317}{773712} & 0 & \frac{1145112371}{2326257360} & \frac{386882707}{156505608} & - \frac{7360000000}{366413327} & \frac{721}{40} \end{array} \end{split}\]

• stages: 6

• order: 5

• stages order: 1

• stability function:

  \[ \frac{13456 z^{6}}{9548775} + \frac{z^{5}}{120} + \frac{z^{4}}{24} + \frac{z^{3}}{6} + \frac{z^{2}}{2} + z + 1 \]

• bibliography: Michael, Akpini K. A. & Richard, Assui K & Gozo, Yoro & Bale, Bailly, , 2019, JMR

Template Parameters:

value_t – type of coefficient (doubleby default)

template<typename value_t>
using ponio::runge_kutta::rk56a_t = explicit_runge_kutta::explicit_runge_kutta<butcher_rk56a<value_t>>

RK56a method.

see more on ponio

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|cccccc} 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ \frac{1}{6} & \frac{1}{6} & 0 & 0 & 0 & 0 & 0 \\ \frac{12}{37} & \frac{12}{1369} & \frac{432}{1369} & 0 & 0 & 0 & 0 \\ \frac{15}{16} & \frac{84995}{28416} & - \frac{66417}{9472} & \frac{119}{24} & 0 & 0 & 0 \\ \frac{74}{75} & \frac{9130103592001}{2191442343750} & - \frac{238171162168}{24349359375} & \frac{2910761155207}{438288468750} & - \frac{14380276736}{365240390625} & 0 & 0 \\ 1 & \frac{4874959019}{1089361548} & - \frac{317812436}{30260043} & \frac{47844943720}{6764346369} & - \frac{3895040}{144757503} & - \frac{111796875}{6026555708} & 0 \\ \hline & \frac{15557}{159840} & 0 & \frac{1182595591}{2401898400} & \frac{3227648}{963765} & - \frac{335390625}{32098832} & \frac{751}{100} \end{array} \end{split}\]

• stages: 6

• order: 5

• stages order: 1

• stability function:

  \[ \frac{53312 z^{6}}{37271025} + \frac{z^{5}}{120} + \frac{z^{4}}{24} + \frac{z^{3}}{6} + \frac{z^{2}}{2} + z + 1 \]

• bibliography: Michael, Akpini K. A. & Richard, Assui K & Gozo, Yoro & Bale, Bailly, , 2019, JMR

Template Parameters:

value_t – type of coefficient (doubleby default)

template<typename value_t>
using ponio::runge_kutta::rk6es_t = explicit_runge_kutta::explicit_runge_kutta<butcher_rk6es<value_t>>

RK6ES method.

see more on ponio

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|ccccccc} 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ \frac{202276644898141}{1000000000000000} & \frac{202276644898141}{1000000000000000} & 0 & 0 & 0 & 0 & 0 & 0 \\ \frac{303414967347211}{1000000000000000} & \frac{758537418368027}{10000000000000000} & \frac{28445153188801}{125000000000000} & 0 & 0 & 0 & 0 & 0 \\ \frac{7}{8} & \frac{13592822172833}{10000000000000} & - \frac{523788570262881}{100000000000000} & \frac{475360348534551}{100000000000000} & 0 & 0 & 0 & 0 \\ \frac{1}{2} & - \frac{160546001129011}{500000000000000} & \frac{82567656396119}{50000000000000} & - \frac{905286676763721}{1000000000000000} & \frac{375127755496799}{5000000000000000} & 0 & 0 & 0 \\ \frac{1}{8} & \frac{73080459837341}{250000000000000} & - \frac{74826938608983}{100000000000000} & \frac{118494168993297}{200000000000000} & - \frac{98886347122859}{2500000000000000} & \frac{56062481246249}{2000000000000000} & 0 & 0 \\ 1 & - \frac{206627618949041}{10000000000000} & \frac{638523209463321}{10000000000000} & - \frac{92689688684611}{1250000000000} & \frac{108014705546673}{125000000000000} & \frac{36263704148237}{2500000000000} & \frac{82962962962963}{5000000000000} & 0 \\ \hline & \frac{142857142857143}{10000000000000000} & 0 & 0 & \frac{270899470899471}{1000000000000000} & \frac{42962962962963}{100000000000000} & \frac{270899470899471}{1000000000000000} & \frac{142857142857143}{10000000000000000} \end{array} \end{split}\]

• stages: 7

• order: 6

• stages order: 0

• stability function:

  \[ \frac{681733599944196256481840469616871820319429391550419341109169405135576788400240780340726058231849169 z^{7}}{6250000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000} + \frac{694444444444442776900226135732387139110736631542963745911023814596592719893682433983637 z^{6}}{500000000000000000000000000000000000000000000000000000000000000000000000000000000000000000} + \frac{83333333333333377197708904510207888651440709904449388217043165275172792917 z^{5}}{10000000000000000000000000000000000000000000000000000000000000000000000000000} + \frac{833333333333329808545102361789100693956364827031005789219211 z^{4}}{20000000000000000000000000000000000000000000000000000000000000} + \frac{1666666666666657993708579763794892360468111301 z^{3}}{10000000000000000000000000000000000000000000000} + \frac{4999999999999986695238095238081 z^{2}}{10000000000000000000000000000000} + \frac{5000000000000003 z}{5000000000000000} + 1 \]

• bibliography: Lawson, J. Douglas, , 1967, SIAM J. Numer. Anal.

Template Parameters:

value_t – type of coefficient (doubleby default)

template<typename value_t>
using ponio::runge_kutta::rk_118_t = explicit_runge_kutta::explicit_runge_kutta<butcher_rk_118<value_t>>

RK (11,8) method.

see more on ponio

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|ccccccccccc} 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ \frac{1}{2} & \frac{1}{2} & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ \frac{1}{2} & \frac{1}{4} & \frac{1}{4} & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ \frac{1}{2} - \frac{\sqrt{21}}{14} & \frac{1}{7} & - \frac{1}{14} + \frac{3 \sqrt{21}}{98} & \frac{3}{7} - \frac{5 \sqrt{21}}{49} & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ \frac{1}{2} - \frac{\sqrt{21}}{14} & \frac{11}{84} - \frac{\sqrt{21}}{84} & 0 & \frac{2}{7} - \frac{4 \sqrt{21}}{63} & \frac{\sqrt{21}}{252} + \frac{1}{12} & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ \frac{1}{2} & \frac{5}{48} - \frac{\sqrt{21}}{48} & 0 & \frac{1}{4} - \frac{\sqrt{21}}{36} & - \frac{77}{120} - \frac{7 \sqrt{21}}{180} & \frac{7 \sqrt{21}}{80} + \frac{63}{80} & 0 & 0 & 0 & 0 & 0 & 0 \\ \frac{\sqrt{21}}{14} + \frac{1}{2} & \frac{\sqrt{21}}{42} + \frac{5}{21} & 0 & - \frac{48}{35} - \frac{92 \sqrt{21}}{315} & \frac{211}{30} + \frac{29 \sqrt{21}}{18} & - \frac{23 \sqrt{21}}{14} - \frac{36}{5} & \frac{13 \sqrt{21}}{35} + \frac{9}{5} & 0 & 0 & 0 & 0 & 0 \\ \frac{\sqrt{21}}{14} + \frac{1}{2} & \frac{1}{14} & 0 & 0 & 0 & \frac{\sqrt{21}}{42} + \frac{1}{9} & \frac{13}{63} + \frac{\sqrt{21}}{21} & \frac{1}{9} & 0 & 0 & 0 & 0 \\ \frac{1}{2} & \frac{1}{32} & 0 & 0 & 0 & \frac{91}{576} + \frac{7 \sqrt{21}}{192} & \frac{11}{72} & - \frac{385}{1152} + \frac{25 \sqrt{21}}{384} & \frac{63}{128} - \frac{13 \sqrt{21}}{128} & 0 & 0 & 0 \\ \frac{1}{2} - \frac{\sqrt{21}}{14} & \frac{1}{14} & 0 & 0 & 0 & \frac{1}{9} & - \frac{733}{2205} + \frac{\sqrt{21}}{15} & \frac{515}{504} - \frac{37 \sqrt{21}}{168} & - \frac{51}{56} + \frac{11 \sqrt{21}}{56} & \frac{132}{245} - \frac{4 \sqrt{21}}{35} & 0 & 0 \\ 1 & 0 & 0 & 0 & 0 & - \frac{7}{3} - \frac{7 \sqrt{21}}{18} & - \frac{28 \sqrt{21}}{45} - \frac{2}{5} & - \frac{91}{24} + \frac{53 \sqrt{21}}{72} & \frac{301}{72} - \frac{53 \sqrt{21}}{72} & \frac{28}{45} + \frac{28 \sqrt{21}}{45} & \frac{7 \sqrt{21}}{18} + \frac{49}{18} & 0 \\ \hline & \frac{1}{20} & 0 & 0 & 0 & 0 & 0 & 0 & \frac{49}{180} & \frac{16}{45} & \frac{49}{180} & \frac{1}{20} \end{array} \end{split}\]

• stages: 11

• order: 8

• stages order: 1

• stability function:

  \[ z^{11} \left(\frac{1}{1382400} - \frac{11 \sqrt{21}}{67737600}\right) + z^{10} \left(- \frac{13}{8467200} + \frac{19 \sqrt{21}}{59270400}\right) + z^{9} \left(- \frac{191}{16934400} + \frac{13 \sqrt{21}}{5644800}\right) + \frac{z^{8}}{40320} + \frac{z^{7}}{5040} + \frac{z^{6}}{720} + \frac{z^{5}}{120} + \frac{z^{4}}{24} + \frac{z^{3}}{6} + \frac{z^{2}}{2} + z + 1 \]

• bibliography: Cooper, G. J. & Verner, J. H., , 1972, SIAM J. Numer. Anal.

Template Parameters:

value_t – type of coefficient (doubleby default)

template<typename value_t>
using ponio::runge_kutta::rk_21_t = explicit_runge_kutta::explicit_runge_kutta<butcher_rk_21<value_t>>

RK (2,1) method.

see more on ponio

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|cc} 0 & 0 & 0 \\ \frac{1}{2} & 1 & 0 \\ \hline & \frac{1}{2} & \frac{1}{2} \end{array} \end{split}\]

• stages: 2

• order: 2

• stages order: 0

• stability function:

  \[ \frac{z^{2}}{2} + z + 1 \]

Template Parameters:

value_t – type of coefficient (doubleby default)

template<typename value_t>
using ponio::runge_kutta::rk_22_midpoint_t = explicit_runge_kutta::explicit_runge_kutta<butcher_rk_22_midpoint<value_t>>

RK (2,2) mid-point method.

see more on ponio

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|cc} 0 & 0 & 0 \\ \frac{1}{2} & \frac{1}{2} & 0 \\ \hline & 0 & 1 \end{array} \end{split}\]

• stages: 2

• order: 2

• stages order: 1

• stability function:

  \[ \frac{z^{2}}{2} + z + 1 \]

Template Parameters:

value_t – type of coefficient (doubleby default)

template<typename value_t>
using ponio::runge_kutta::rk_22_ralston_t = explicit_runge_kutta::explicit_runge_kutta<butcher_rk_22_ralston<value_t>>

RK (2,2) Ralston method.

see more on ponio

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|cc} 0 & 0 & 0 \\ \frac{2}{3} & \frac{2}{3} & 0 \\ \hline & \frac{1}{4} & \frac{3}{4} \end{array} \end{split}\]

• stages: 2

• order: 2

• stages order: 1

• stability function:

  \[ \frac{z^{2}}{2} + z + 1 \]

Template Parameters:

value_t – type of coefficient (doubleby default)

template<typename value_t>
using ponio::runge_kutta::rk_32_best_t = explicit_runge_kutta::explicit_runge_kutta<butcher_rk_32_best<value_t>>

RK (3,2) best method.

see more on ponio

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|ccc} 0 & 0 & 0 & 0 \\ \frac{1}{2} & \frac{1}{2} & 0 & 0 \\ \frac{1}{2} & 0 & \frac{1}{2} & 0 \\ \hline & 0 & 0 & 1 \end{array} \end{split}\]

• stages: 3

• order: 2

• stages order: 1

• stability function:

  \[ \frac{z^{3}}{4} + \frac{z^{2}}{2} + z + 1 \]

• bibliography: Crouseilles, Nicolas & Einkemmer, Lukas & Massot, Josselin, , 2020, Journal of Computational Physics

Template Parameters:

value_t – type of coefficient (doubleby default)

template<typename value_t>
using ponio::runge_kutta::rk_33_t = explicit_runge_kutta::explicit_runge_kutta<butcher_rk_33<value_t>>

RK (3,3) method.

see more on ponio

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|ccc} 0 & 0 & 0 & 0 \\ \frac{1}{2} & \frac{1}{2} & 0 & 0 \\ 1 & -1 & 2 & 0 \\ \hline & \frac{1}{6} & \frac{2}{3} & \frac{1}{6} \end{array} \end{split}\]

• stages: 3

• order: 3

• stages order: 1

• stability function:

  \[ \frac{z^{3}}{6} + \frac{z^{2}}{2} + z + 1 \]

Template Parameters:

value_t – type of coefficient (doubleby default)

template<typename value_t>
using ponio::runge_kutta::rk_33_233e_t = explicit_runge_kutta::explicit_runge_kutta<butcher_rk_33_233e<value_t>>

RK (3,3) 233e method.

see more on ponio

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|ccc} 0 & 0 & 0 & 0 \\ \frac{2}{3} & \frac{2}{3} & 0 & 0 \\ \frac{2}{3} & \frac{1}{3} & \frac{1}{3} & 0 \\ \hline & \frac{1}{4} & 0 & \frac{3}{4} \end{array} \end{split}\]

• stages: 3

• order: 3

• stages order: 1

• stability function:

  \[ \frac{z^{3}}{6} + \frac{z^{2}}{2} + z + 1 \]

• bibliography: Butcher, J. C., , 2008, Wiley

Template Parameters:

value_t – type of coefficient (doubleby default)

template<typename value_t>
using ponio::runge_kutta::rk_33_bogackishampine_t = explicit_runge_kutta::explicit_runge_kutta<butcher_rk_33_bogackishampine<value_t>>

RK (3,3) Bogacki-Shampine method.

see more on ponio

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|ccc} 0 & 0 & 0 & 0 \\ \frac{1}{2} & \frac{1}{2} & 0 & 0 \\ \frac{3}{4} & 0 & \frac{3}{4} & 0 \\ \hline & \frac{2}{9} & \frac{1}{3} & \frac{4}{9} \end{array} \end{split}\]

• stages: 3

• order: 3

• stages order: 1

• stability function:

  \[ \frac{z^{3}}{6} + \frac{z^{2}}{2} + z + 1 \]

Template Parameters:

value_t – type of coefficient (doubleby default)

template<typename value_t>
using ponio::runge_kutta::rk_33_heun_t = explicit_runge_kutta::explicit_runge_kutta<butcher_rk_33_heun<value_t>>

RK (3,3) Heun method.

see more on ponio

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|ccc} 0 & 0 & 0 & 0 \\ \frac{1}{3} & \frac{1}{3} & 0 & 0 \\ \frac{2}{3} & 0 & \frac{2}{3} & 0 \\ \hline & \frac{1}{4} & 0 & \frac{3}{4} \end{array} \end{split}\]

• stages: 3

• order: 3

• stages order: 1

• stability function:

  \[ \frac{z^{3}}{6} + \frac{z^{2}}{2} + z + 1 \]

Template Parameters:

value_t – type of coefficient (doubleby default)

template<typename value_t>
using ponio::runge_kutta::rk_33_ralston_t = explicit_runge_kutta::explicit_runge_kutta<butcher_rk_33_ralston<value_t>>

RK (3,3) Ralston method.

see more on ponio

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|ccc} 0 & 0 & 0 & 0 \\ \frac{1}{2} & \frac{1}{2} & 0 & 0 \\ \frac{3}{4} & 0 & \frac{3}{4} & 0 \\ \hline & \frac{2}{9} & \frac{1}{3} & \frac{4}{9} \end{array} \end{split}\]

• stages: 3

• order: 3

• stages order: 1

• stability function:

  \[ \frac{z^{3}}{6} + \frac{z^{2}}{2} + z + 1 \]

Template Parameters:

value_t – type of coefficient (doubleby default)

template<typename value_t>
using ponio::runge_kutta::rk_33_van_der_houwen_t = explicit_runge_kutta::explicit_runge_kutta<butcher_rk_33_van_der_houwen<value_t>>

RK (3,3) Van der Houwen method.

see more on ponio

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|ccc} 0 & 0 & 0 & 0 \\ \frac{8}{15} & \frac{8}{15} & 0 & 0 \\ \frac{2}{3} & \frac{1}{4} & \frac{5}{12} & 0 \\ \hline & \frac{1}{4} & 0 & \frac{3}{4} \end{array} \end{split}\]

• stages: 3

• order: 3

• stages order: 1

• stability function:

  \[ \frac{z^{3}}{6} + \frac{z^{2}}{2} + z + 1 \]

Template Parameters:

value_t – type of coefficient (doubleby default)

template<typename value_t>
using ponio::runge_kutta::rk_44_t = explicit_runge_kutta::explicit_runge_kutta<butcher_rk_44<value_t>>

RK (4,4) method.

see more on ponio

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|cccc} 0 & 0 & 0 & 0 & 0 \\ \frac{1}{2} & \frac{1}{2} & 0 & 0 & 0 \\ \frac{1}{2} & 0 & \frac{1}{2} & 0 & 0 \\ 1 & 0 & 0 & 1 & 0 \\ \hline & \frac{1}{6} & \frac{1}{3} & \frac{1}{3} & \frac{1}{6} \end{array} \end{split}\]

• stages: 4

• order: 4

• stages order: 1

• stability function:

  \[ \frac{z^{4}}{24} + \frac{z^{3}}{6} + \frac{z^{2}}{2} + z + 1 \]

Template Parameters:

value_t – type of coefficient (doubleby default)

template<typename value_t>
using ponio::runge_kutta::rk_44_235j_t = explicit_runge_kutta::explicit_runge_kutta<butcher_rk_44_235j<value_t>>

RK (4,4) 235j method.

see more on ponio

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|cccc} 0 & 0 & 0 & 0 & 0 \\ \frac{1}{4} & \frac{1}{4} & 0 & 0 & 0 \\ \frac{1}{2} & 0 & \frac{1}{2} & 0 & 0 \\ 1 & 1 & -2 & 2 & 0 \\ \hline & \frac{1}{6} & 0 & \frac{2}{3} & \frac{1}{6} \end{array} \end{split}\]

• stages: 4

• order: 4

• stages order: 1

• stability function:

  \[ \frac{z^{4}}{24} + \frac{z^{3}}{6} + \frac{z^{2}}{2} + z + 1 \]

• bibliography: Butcher, J. C., , 2008, Wiley

Template Parameters:

value_t – type of coefficient (doubleby default)

template<typename value_t>
using ponio::runge_kutta::rk_44_38_t = explicit_runge_kutta::explicit_runge_kutta<butcher_rk_44_38<value_t>>

RK (4,4) 3/8 method.

see more on ponio

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|cccc} 0 & 0 & 0 & 0 & 0 \\ \frac{1}{3} & \frac{1}{3} & 0 & 0 & 0 \\ \frac{2}{3} & - \frac{1}{3} & 1 & 0 & 0 \\ 1 & 1 & -1 & 1 & 0 \\ \hline & \frac{1}{8} & \frac{3}{8} & \frac{3}{8} & \frac{1}{8} \end{array} \end{split}\]

• stages: 4

• order: 4

• stages order: 1

• stability function:

  \[ \frac{z^{4}}{24} + \frac{z^{3}}{6} + \frac{z^{2}}{2} + z + 1 \]

Template Parameters:

value_t – type of coefficient (doubleby default)

template<typename value_t>
using ponio::runge_kutta::rk_44_ralston_t = explicit_runge_kutta::explicit_runge_kutta<butcher_rk_44_ralston<value_t>>

RK (4,4) Ralston method.

see more on ponio

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|cccc} 0 & 0 & 0 & 0 & 0 \\ \frac{2}{5} & \frac{2}{5} & 0 & 0 & 0 \\ - \frac{1822949}{4000000} & \frac{29697761}{100000000} & \frac{3968991}{25000000} & 0 & 0 \\ 1 & \frac{545251}{2500000} & - \frac{76274129}{25000000} & \frac{95821619}{25000000} & 0 \\ \hline & \frac{4369007}{25000000} & - \frac{27574033}{50000000} & \frac{3013839}{2500000} & \frac{8559239}{50000000} \end{array} \end{split}\]

• stages: 4

• order: 4

• stages order: 0

• stability function:

  \[ \frac{3255208207820492337531 z^{4}}{78125000000000000000000} + \frac{833333300435997058009 z^{3}}{5000000000000000000000} + \frac{4999999951211 z^{2}}{10000000000000} + z + 1 \]

Template Parameters:

value_t – type of coefficient (doubleby default)

template<typename value_t>
using ponio::runge_kutta::rk_65_t = explicit_runge_kutta::explicit_runge_kutta<butcher_rk_65<value_t>>

RK (6,5) method.

see more on ponio

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|cccccc} 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ \frac{1}{5} & \frac{1}{5} & 0 & 0 & 0 & 0 & 0 \\ \frac{3}{10} & \frac{3}{40} & \frac{9}{40} & 0 & 0 & 0 & 0 \\ \frac{4}{5} & \frac{44}{45} & - \frac{56}{15} & \frac{32}{9} & 0 & 0 & 0 \\ \frac{8}{9} & \frac{19372}{6561} & - \frac{25360}{2187} & \frac{64448}{6561} & - \frac{212}{729} & 0 & 0 \\ 1 & \frac{9017}{3168} & - \frac{355}{33} & \frac{46732}{5247} & \frac{49}{176} & - \frac{5103}{18656} & 0 \\ \hline & \frac{35}{384} & 0 & \frac{500}{1113} & \frac{125}{192} & - \frac{2187}{6784} & \frac{11}{84} \end{array} \end{split}\]

• stages: 6

• order: 5

• stages order: 1

• stability function:

  \[ \frac{z^{6}}{600} + \frac{z^{5}}{120} + \frac{z^{4}}{24} + \frac{z^{3}}{6} + \frac{z^{2}}{2} + z + 1 \]

Template Parameters:

value_t – type of coefficient (doubleby default)

template<typename value_t>
using ponio::runge_kutta::rk_65_236a_t = explicit_runge_kutta::explicit_runge_kutta<butcher_rk_65_236a<value_t>>

RK (6,5) 236a method.

see more on ponio

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|cccccc} 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ \frac{1}{4} & \frac{1}{4} & 0 & 0 & 0 & 0 & 0 \\ \frac{1}{4} & \frac{1}{8} & \frac{1}{8} & 0 & 0 & 0 & 0 \\ \frac{1}{2} & 0 & 0 & \frac{1}{2} & 0 & 0 & 0 \\ \frac{3}{4} & \frac{3}{16} & - \frac{3}{8} & \frac{3}{8} & \frac{9}{16} & 0 & 0 \\ 1 & - \frac{3}{7} & \frac{8}{7} & \frac{6}{7} & - \frac{12}{7} & \frac{8}{7} & 0 \\ \hline & \frac{7}{90} & 0 & \frac{16}{45} & \frac{2}{15} & \frac{16}{45} & \frac{7}{90} \end{array} \end{split}\]

• stages: 6

• order: 5

• stages order: 1

• stability function:

  \[ \frac{z^{6}}{1280} + \frac{z^{5}}{120} + \frac{z^{4}}{24} + \frac{z^{3}}{6} + \frac{z^{2}}{2} + z + 1 \]

• bibliography: Butcher, J. C., , 2008, Wiley

Template Parameters:

value_t – type of coefficient (doubleby default)

template<typename value_t>
using ponio::runge_kutta::rk_76_t = explicit_runge_kutta::explicit_runge_kutta<butcher_rk_76<value_t>>

RK (7,6) method.

see more on ponio

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|ccccccc} 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ 1 & 1 & 0 & 0 & 0 & 0 & 0 & 0 \\ \frac{1}{2} & \frac{3}{8} & \frac{1}{8} & 0 & 0 & 0 & 0 & 0 \\ \frac{2}{3} & \frac{8}{27} & \frac{2}{27} & \frac{8}{27} & 0 & 0 & 0 & 0 \\ \frac{1}{2} - \frac{\sqrt{21}}{14} & - \frac{3}{56} + \frac{9 \sqrt{21}}{392} & - \frac{1}{7} + \frac{\sqrt{21}}{49} & \frac{6}{7} - \frac{6 \sqrt{21}}{49} & - \frac{9}{56} + \frac{3 \sqrt{21}}{392} & 0 & 0 & 0 \\ \frac{\sqrt{21}}{14} + \frac{1}{2} & - \frac{51 \sqrt{21}}{392} - \frac{33}{56} & - \frac{1}{7} - \frac{\sqrt{21}}{49} & - \frac{8 \sqrt{21}}{49} & \frac{9}{280} + \frac{363 \sqrt{21}}{1960} & \frac{\sqrt{21}}{5} + \frac{6}{5} & 0 & 0 \\ 1 & \frac{11}{6} + \frac{7 \sqrt{21}}{12} & \frac{2}{3} & - \frac{10}{9} + \frac{14 \sqrt{21}}{9} & \frac{7}{10} - \frac{21 \sqrt{21}}{20} & - \frac{343}{90} - \frac{7 \sqrt{21}}{10} & \frac{49}{18} - \frac{7 \sqrt{21}}{18} & 0 \\ \hline & \frac{1}{20} & 0 & \frac{16}{45} & 0 & \frac{49}{180} & \frac{49}{180} & \frac{1}{20} \end{array} \end{split}\]

• stages: 7

• order: 6

• stages order: 1

• stability function:

  \[ - \frac{z^{7}}{2160} + \frac{z^{6}}{720} + \frac{z^{5}}{120} + \frac{z^{4}}{24} + \frac{z^{3}}{6} + \frac{z^{2}}{2} + z + 1 \]

• bibliography: Luther, H. A., , 1968, Math. Comp.

Template Parameters:

value_t – type of coefficient (doubleby default)

template<typename value_t>
using ponio::runge_kutta::rk_86_t = explicit_runge_kutta::explicit_runge_kutta<butcher_rk_86<value_t>>

RK (8,6) method.

see more on ponio

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|cccccccc} 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ \frac{1}{9} & \frac{1}{9} & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ \frac{1}{6} & \frac{1}{24} & \frac{1}{8} & 0 & 0 & 0 & 0 & 0 & 0 \\ \frac{1}{3} & \frac{1}{6} & - \frac{1}{2} & \frac{2}{3} & 0 & 0 & 0 & 0 & 0 \\ \frac{1}{2} & \frac{935}{2536} & - \frac{2781}{2536} & \frac{309}{317} & \frac{321}{1268} & 0 & 0 & 0 & 0 \\ \frac{2}{3} & - \frac{12710}{951} & \frac{8287}{317} & - \frac{40}{317} & - \frac{6335}{317} & 8 & 0 & 0 & 0 \\ \frac{5}{6} & \frac{5840285}{3104064} & - \frac{7019}{2536} & - \frac{52213}{86224} & \frac{1278709}{517344} & - \frac{433}{2448} & \frac{33}{1088} & 0 & 0 \\ 1 & - \frac{5101675}{1767592} & \frac{112077}{25994} & \frac{334875}{441898} & - \frac{973617}{883796} & - \frac{1421}{1394} & \frac{333}{5576} & \frac{36}{41} & 0 \\ \hline & \frac{41}{840} & 0 & \frac{9}{35} & \frac{9}{280} & \frac{34}{105} & \frac{9}{280} & \frac{9}{35} & \frac{41}{840} \end{array} \end{split}\]

• stages: 8

• order: 6

• stages order: 1

• stability function:

  \[ \frac{1177 z^{8}}{48285440} + \frac{18713 z^{7}}{81481680} + \frac{z^{6}}{720} + \frac{z^{5}}{120} + \frac{z^{4}}{24} + \frac{z^{3}}{6} + \frac{z^{2}}{2} + z + 1 \]

Template Parameters:

value_t – type of coefficient (doubleby default)

template<typename value_t>
using ponio::runge_kutta::rk_nssp_21_t = explicit_runge_kutta::explicit_runge_kutta<butcher_rk_nssp_21<value_t>>

RK NSSP (2,1) method.

see more on ponio

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|cc} 0 & 0 & 0 \\ \frac{3}{4} & \frac{3}{4} & 0 \\ \hline & 0 & 1 \end{array} \end{split}\]

• stages: 2

• order: 1

• stages order: 1

• stability function:

  \[ \frac{3 z^{2}}{4} + z + 1 \]

Template Parameters:

value_t – type of coefficient (doubleby default)

template<typename value_t>
using ponio::runge_kutta::rk_nssp_32_t = explicit_runge_kutta::explicit_runge_kutta<butcher_rk_nssp_32<value_t>>

RK NSSP (3,2) method.

see more on ponio

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|ccc} 0 & 0 & 0 & 0 \\ \frac{1}{3} & \frac{1}{3} & 0 & 0 \\ 1 & 0 & 1 & 0 \\ \hline & \frac{1}{2} & 0 & \frac{1}{2} \end{array} \end{split}\]

• stages: 3

• order: 3

• stages order: 1

• stability function:

  \[ \frac{z^{3}}{6} + \frac{z^{2}}{2} + z + 1 \]

Template Parameters:

value_t – type of coefficient (doubleby default)

template<typename value_t>
using ponio::runge_kutta::rk_nssp_33_t = explicit_runge_kutta::explicit_runge_kutta<butcher_rk_nssp_33<value_t>>

RK NSSP (3,3) method.

see more on ponio

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|ccc} 0 & 0 & 0 & 0 \\ - \frac{4}{9} & - \frac{4}{9} & 0 & 0 \\ \frac{2}{3} & \frac{7}{6} & - \frac{1}{2} & 0 \\ \hline & \frac{1}{4} & 0 & \frac{3}{4} \end{array} \end{split}\]

• stages: 3

• order: 3

• stages order: 1

• stability function:

  \[ \frac{z^{3}}{6} + \frac{z^{2}}{2} + z + 1 \]

Template Parameters:

value_t – type of coefficient (doubleby default)

template<typename value_t>
using ponio::runge_kutta::rk_nssp_53_t = explicit_runge_kutta::explicit_runge_kutta<butcher_rk_nssp_53<value_t>>

RK NSSP (5,3) method.

see more on ponio

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|ccccc} 0 & 0 & 0 & 0 & 0 & 0 \\ \frac{1}{7} & \frac{1}{7} & 0 & 0 & 0 & 0 \\ \frac{3}{16} & 0 & \frac{3}{16} & 0 & 0 & 0 \\ \frac{1}{3} & 0 & 0 & \frac{1}{3} & 0 & 0 \\ \frac{2}{3} & 0 & 0 & 0 & \frac{2}{3} & 0 \\ \hline & \frac{1}{4} & 0 & 0 & 0 & \frac{3}{4} \end{array} \end{split}\]

• stages: 5

• order: 3

• stages order: 1

• stability function:

  \[ \frac{z^{5}}{224} + \frac{z^{4}}{32} + \frac{z^{3}}{6} + \frac{z^{2}}{2} + z + 1 \]

Template Parameters:

value_t – type of coefficient (doubleby default)

template<typename value_t>
using ponio::runge_kutta::rk_spp_43_t = explicit_runge_kutta::explicit_runge_kutta<butcher_rk_spp_43<value_t>>

RK SPP (4,3) method.

see more on ponio

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|cccc} 0 & 0 & 0 & 0 & 0 \\ \frac{1}{2} & \frac{1}{2} & 0 & 0 & 0 \\ 1 & \frac{1}{2} & \frac{1}{2} & 0 & 0 \\ \frac{1}{2} & \frac{1}{6} & \frac{1}{6} & \frac{1}{6} & 0 \\ \hline & \frac{1}{6} & \frac{1}{6} & \frac{1}{6} & \frac{1}{2} \end{array} \end{split}\]

• stages: 4

• order: 3

• stages order: 1

• stability function:

  \[ \frac{z^{4}}{48} + \frac{z^{3}}{6} + \frac{z^{2}}{2} + z + 1 \]

• bibliography: Spiteri, Raymond J. & Ruuth, Steven J., , 2002, SIAM J. Numer. Anal.

Template Parameters:

value_t – type of coefficient (doubleby default)

template<typename value_t>
using ponio::runge_kutta::rk_ssp_22_heun_t = explicit_runge_kutta::explicit_runge_kutta<butcher_rk_ssp_22_heun<value_t>>

RK SSP (2,2) Heun method.

see more on ponio

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|cc} 0 & 0 & 0 \\ 1 & 1 & 0 \\ \hline & \frac{1}{2} & \frac{1}{2} \end{array} \end{split}\]

• stages: 2

• order: 2

• stages order: 1

• stability function:

  \[ \frac{z^{2}}{2} + z + 1 \]

Template Parameters:

value_t – type of coefficient (doubleby default)

template<typename value_t>
using ponio::runge_kutta::rk_ssp_32_t = explicit_runge_kutta::explicit_runge_kutta<butcher_rk_ssp_32<value_t>>

RK SSP (3,2) method.

see more on ponio

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|ccc} 0 & 0 & 0 & 0 \\ \frac{1}{2} & \frac{1}{2} & 0 & 0 \\ 1 & \frac{1}{2} & \frac{1}{2} & 0 \\ \hline & \frac{1}{3} & \frac{1}{3} & \frac{1}{3} \end{array} \end{split}\]

• stages: 3

• order: 2

• stages order: 1

• stability function:

  \[ \frac{z^{3}}{12} + \frac{z^{2}}{2} + z + 1 \]

• bibliography: Spiteri, Raymond J. & Ruuth, Steven J., , 2002, SIAM J. Numer. Anal.

Template Parameters:

value_t – type of coefficient (doubleby default)

template<typename value_t>
using ponio::runge_kutta::rk_ssp_33_t = explicit_runge_kutta::explicit_runge_kutta<butcher_rk_ssp_33<value_t>>

RK SSP (3,3) method.

see more on ponio

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|ccc} 0 & 0 & 0 & 0 \\ 1 & 1 & 0 & 0 \\ \frac{1}{2} & \frac{1}{4} & \frac{1}{4} & 0 \\ \hline & \frac{1}{6} & \frac{1}{6} & \frac{2}{3} \end{array} \end{split}\]

• stages: 3

• order: 3

• stages order: 1

• stability function:

  \[ \frac{z^{3}}{6} + \frac{z^{2}}{2} + z + 1 \]

• bibliography: Shu, Chi-Wang & Osher, Stanley, , 1988, Journal of Computational Physics

Template Parameters:

value_t – type of coefficient (doubleby default)

template<typename value_t>
using ponio::runge_kutta::rk_ssp_42_t = explicit_runge_kutta::explicit_runge_kutta<butcher_rk_ssp_42<value_t>>

RK SSP (4,2) method.

see more on ponio

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|cccc} 0 & 0 & 0 & 0 & 0 \\ \frac{1}{3} & \frac{1}{3} & 0 & 0 & 0 \\ \frac{2}{3} & \frac{1}{3} & \frac{1}{3} & 0 & 0 \\ 1 & \frac{1}{3} & \frac{1}{3} & \frac{1}{3} & 0 \\ \hline & \frac{1}{4} & \frac{1}{4} & \frac{1}{4} & \frac{1}{4} \end{array} \end{split}\]

• stages: 4

• order: 2

• stages order: 1

• stability function:

  \[ \frac{z^{4}}{108} + \frac{z^{3}}{9} + \frac{z^{2}}{2} + z + 1 \]

• bibliography: Spiteri, Raymond J. & Ruuth, Steven J., , 2002, SIAM J. Numer. Anal.

Template Parameters:

value_t – type of coefficient (doubleby default)

template<typename value_t>
using ponio::runge_kutta::rk_ssp_53_t = explicit_runge_kutta::explicit_runge_kutta<butcher_rk_ssp_53<value_t>>

RK SSP (5,3) method.

see more on ponio

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|ccccc} 0 & 0 & 0 & 0 & 0 & 0 \\ \frac{3772689151171}{10000000000000} & \frac{3772689151171}{10000000000000} & 0 & 0 & 0 & 0 \\ \frac{75453783023419}{100000000000000} & \frac{3772689151171}{10000000000000} & \frac{3772689151171}{10000000000000} & 0 & 0 & 0 \\ \frac{24528441134657}{50000000000000} & \frac{16352294089771}{100000000000000} & \frac{16352294089771}{100000000000000} & \frac{16352294089771}{100000000000000} & 0 & 0 \\ \frac{78784303014311}{100000000000000} & \frac{1863007424357}{12500000000000} & \frac{3707818346181}{25000000000000} & \frac{3707818346181}{25000000000000} & \frac{4277212106251}{12500000000000} & 0 \\ \hline & \frac{19707596384481}{100000000000000} & \frac{2356063301953}{20000000000000} & \frac{2927431298443}{25000000000000} & \frac{27015874934251}{100000000000000} & \frac{3723310876263}{12500000000000} \end{array} \end{split}\]

• stages: 5

• order: 3

• stages order: 0

• stability function:

  \[ \frac{3706558184075319324316398834852923260111749972090780583449705543 z^{5}}{1562500000000000000000000000000000000000000000000000000000000000000} + \frac{157195381248937277398992604182081645724718375717682053 z^{4}}{5000000000000000000000000000000000000000000000000000000} + \frac{83333333306904819032998927912645085112413 z^{3}}{500000000000000000000000000000000000000000} + \frac{5000000001360194746491536401 z^{2}}{10000000000000000000000000000} + \frac{100000000032373 z}{100000000000000} + 1 \]

• bibliography: Spiteri, Raymond J. & Ruuth, Steven J., , 2002, SIAM J. Numer. Anal.

Template Parameters:

value_t – type of coefficient (doubleby default)

template<typename value_t>
using ponio::runge_kutta::rk_ssp_54_t = explicit_runge_kutta::explicit_runge_kutta<butcher_rk_ssp_54<value_t>>

RK SSP (5,4) method.

see more on ponio

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|ccccc} 0 & 0 & 0 & 0 & 0 & 0 \\ \frac{4896902837549}{12500000000000} & \frac{4896902837549}{12500000000000} & 0 & 0 & 0 & 0 \\ \frac{58607968896779}{100000000000000} & \frac{21766909633821}{100000000000000} & \frac{36841059262959}{100000000000000} & 0 & 0 & 0 \\ \frac{47454236302687}{100000000000000} & \frac{165384173419}{2000000000000} & \frac{13995850206999}{100000000000000} & \frac{12594588712369}{50000000000000} & 0 & 0 \\ \frac{23375265775231}{25000000000000} & \frac{84957854629}{1250000000000} & \frac{5751734922219}{50000000000000} & \frac{20703489864929}{100000000000000} & \frac{54497475021237}{100000000000000} & 0 \\ \hline & \frac{14681187618661}{100000000000000} & \frac{6212072731139}{25000000000000} & \frac{208517660733}{2000000000000} & \frac{686097252299}{2500000000000} & \frac{4520149663879}{20000000000000} \end{array} \end{split}\]

• stages: 5

• order: 4

• stages order: 0

• stability function:

  \[ \frac{559714787685486221089408690222017843875490474417403023535939367817 z^{5}}{125000000000000000000000000000000000000000000000000000000000000000000} + \frac{20833333324892802302249888799459495814640610965580531 z^{4}}{500000000000000000000000000000000000000000000000000000} + \frac{33333333323124540472487090140319832509121 z^{3}}{200000000000000000000000000000000000000000} + \frac{1249999999913938597141832613 z^{2}}{2500000000000000000000000000} + \frac{49999999995611 z}{50000000000000} + 1 \]

• bibliography: Spiteri, Raymond J. & Ruuth, Steven J., , 2002, SIAM J. Numer. Anal.

Template Parameters:

value_t – type of coefficient (doubleby default)

template<typename value_t>
using ponio::runge_kutta::rkc_202_t = explicit_runge_kutta::explicit_runge_kutta<butcher_rkc_202<value_t>>

RKC (20,2) method.

see more on ponio

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|cccccccccccccccccccc} 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ \frac{1}{133} & \frac{1}{133} & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ \frac{1}{133} & \frac{1}{266} & \frac{1}{266} & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ \frac{8}{399} & - \frac{8}{1197} & \frac{32}{3591} & \frac{64}{3591} & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ \frac{5}{133} & - \frac{255}{8512} & \frac{15}{1064} & \frac{5}{133} & \frac{135}{8512} & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ \frac{8}{133} & - \frac{1076}{16625} & \frac{64}{3325} & \frac{192}{3325} & \frac{108}{3325} & \frac{256}{16625} & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ \frac{5}{57} & - \frac{151}{1368} & \frac{25}{1026} & \frac{40}{513} & \frac{15}{304} & \frac{16}{513} & \frac{125}{8208} & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ \frac{16}{133} & - \frac{38008}{228095} & \frac{192}{6517} & \frac{640}{6517} & \frac{432}{6517} & \frac{1536}{32585} & \frac{200}{6517} & \frac{3456}{228095} & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ \frac{3}{19} & - \frac{79437}{340480} & \frac{21}{608} & \frac{9}{76} & \frac{405}{4864} & \frac{6}{95} & \frac{225}{4864} & \frac{81}{2660} & \frac{147}{9728} & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ \frac{80}{399} & - \frac{210604}{678699} & \frac{1280}{32319} & \frac{640}{4617} & \frac{40}{399} & \frac{2560}{32319} & \frac{2000}{32319} & \frac{128}{2793} & \frac{140}{4617} & \frac{10240}{678699} & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ \frac{33}{133} & - \frac{148049}{372400} & \frac{297}{6650} & \frac{528}{3325} & \frac{891}{7600} & \frac{1584}{16625} & \frac{165}{2128} & \frac{7128}{116375} & \frac{693}{15200} & \frac{704}{23275} & \frac{8019}{532000} & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ \frac{40}{133} & - \frac{1840178}{3717483} & \frac{800}{16093} & \frac{2880}{16093} & \frac{2160}{16093} & \frac{256}{2299} & \frac{1500}{16093} & \frac{8640}{112651} & \frac{140}{2299} & \frac{5120}{112651} & \frac{486}{16093} & \frac{8000}{531069} & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ \frac{143}{399} & - \frac{191347}{317520} & \frac{1573}{28728} & \frac{715}{3591} & \frac{1287}{8512} & \frac{2288}{17955} & \frac{3575}{32832} & \frac{429}{4655} & \frac{5005}{65664} & \frac{4576}{75411} & \frac{3861}{85120} & \frac{325}{10773} & \frac{17303}{1149120} & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ \frac{8}{19} & - \frac{34732694}{48213165} & \frac{192}{3211} & \frac{704}{3211} & \frac{540}{3211} & \frac{2304}{16055} & \frac{400}{3211} & \frac{1728}{16055} & \frac{294}{3211} & \frac{5120}{67431} & \frac{972}{16055} & \frac{1600}{35321} & \frac{484}{16055} & \frac{6912}{459173} & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ \frac{65}{133} & - \frac{20432707}{24086832} & \frac{845}{13034} & \frac{1560}{6517} & \frac{19305}{104272} & \frac{1040}{6517} & \frac{14625}{104272} & \frac{5616}{45619} & \frac{65}{608} & \frac{4160}{45619} & \frac{15795}{208544} & \frac{13000}{215061} & \frac{4719}{104272} & \frac{2160}{71687} & \frac{10985}{729904} & 0 & 0 & 0 & 0 & 0 & 0 \\ \frac{32}{57} & - \frac{37985848}{38513475} & \frac{896}{12825} & \frac{3328}{12825} & \frac{96}{475} & \frac{11264}{64125} & \frac{80}{513} & \frac{2304}{16625} & \frac{1568}{12825} & \frac{4096}{38475} & \frac{216}{2375} & \frac{1280}{16929} & \frac{3872}{64125} & \frac{3072}{67925} & \frac{2704}{89775} & \frac{12544}{833625} & 0 & 0 & 0 & 0 & 0 \\ \frac{85}{133} & - \frac{1855818719}{1635938304} & \frac{1275}{17024} & \frac{85}{304} & \frac{29835}{136192} & \frac{51}{266} & \frac{23375}{136192} & \frac{2295}{14896} & \frac{5355}{38912} & \frac{340}{2793} & \frac{4131}{38912} & \frac{2125}{23408} & \frac{10285}{136192} & \frac{2295}{38038} & \frac{43095}{953344} & \frac{119}{3952} & \frac{57375}{3813376} & 0 & 0 & 0 & 0 \\ \frac{96}{133} & - \frac{4227342628}{3270412145} & \frac{3072}{38437} & \frac{11520}{38437} & \frac{1296}{5491} & \frac{39936}{192185} & \frac{7200}{38437} & \frac{228096}{1345295} & \frac{840}{5491} & \frac{36864}{269059} & \frac{23328}{192185} & \frac{6400}{60401} & \frac{17424}{192185} & \frac{414720}{5496491} & \frac{16224}{269059} & \frac{16128}{356915} & \frac{8100}{269059} & \frac{49152}{3267145} & 0 & 0 & 0 \\ \frac{17}{21} & - \frac{1193972657}{817296480} & \frac{289}{3402} & \frac{544}{1701} & \frac{85}{336} & \frac{272}{1215} & \frac{5525}{27216} & \frac{136}{735} & \frac{1309}{7776} & \frac{5440}{35721} & \frac{153}{1120} & \frac{6800}{56133} & \frac{2057}{19440} & \frac{272}{3003} & \frac{14365}{190512} & \frac{952}{15795} & \frac{425}{9408} & \frac{256}{8505} & \frac{4913}{326592} & 0 & 0 \\ \frac{120}{133} & - \frac{25446999243}{15523707199} & \frac{4320}{48013} & \frac{16320}{48013} & \frac{12960}{48013} & \frac{11520}{48013} & \frac{1500}{6859} & \frac{67392}{336091} & \frac{1260}{6859} & \frac{56320}{336091} & \frac{7290}{48013} & \frac{72000}{528143} & \frac{5808}{48013} & \frac{103680}{980837} & \frac{30420}{336091} & \frac{6720}{89167} & \frac{20250}{336091} & \frac{36864}{816221} & \frac{1445}{48013} & \frac{233280}{15508199} & 0 \\ \hline & - \frac{9454786563}{5173168000} & \frac{19}{200} & \frac{9}{25} & \frac{459}{1600} & \frac{32}{125} & \frac{15}{64} & \frac{27}{125} & \frac{637}{3200} & \frac{32}{175} & \frac{2673}{16000} & \frac{5}{33} & \frac{1089}{8000} & \frac{432}{3575} & \frac{169}{1600} & \frac{147}{1625} & \frac{135}{1792} & \frac{128}{2125} & \frac{289}{6400} & \frac{243}{8075} & \frac{361}{24000} \end{array} \end{split}\]

• stages: 20

• order: 2

• stages order: 1

• stability function:

  \[ \frac{32768 z^{20}}{563794957188848500305686224618190097399925} + \frac{131072 z^{19}}{847811965697516541813061991907052778045} + \frac{1212416 z^{18}}{6374526057876064224158360841406411865} + \frac{196608 z^{17}}{1369393352927188877370217151752183} + \frac{765952 z^{16}}{10296190623512698326091858283851} + \frac{379912192 z^{15}}{13547619241464076744857708268225} + \frac{162021376 z^{14}}{20372359761600115405801065065} + \frac{265125888 z^{13}}{153175637305264025607526805} + \frac{67317120 z^{12}}{230339304218442143770717} + \frac{66593280 z^{11}}{1731874467807835667449} + \frac{36626304 z^{10}}{9301151814220384895} + \frac{229632 z^{9}}{736141813551277} + \frac{731952 z^{8}}{38744305976383} + \frac{250240 z^{7}}{291310571251} + \frac{62560 z^{6}}{2190305047} + \frac{275264 z^{5}}{411711475} + \frac{12903 z^{4}}{1238230} + \frac{66 z^{3}}{665} + \frac{z^{2}}{2} + z + 1 \]

Template Parameters:

value_t – type of coefficient (doubleby default)

template<typename value_t>
using ponio::runge_kutta::rkc_51_t = explicit_runge_kutta::explicit_runge_kutta<butcher_rkc_51<value_t>>

RKC (5,1) method.

see more on ponio

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|ccccc} 0 & 0 & 0 & 0 & 0 & 0 \\ \frac{1}{25} & \frac{1}{25} & 0 & 0 & 0 & 0 \\ \frac{4}{25} & \frac{2}{25} & \frac{2}{25} & 0 & 0 & 0 \\ \frac{9}{25} & \frac{3}{25} & \frac{4}{25} & \frac{2}{25} & 0 & 0 \\ \frac{16}{25} & \frac{4}{25} & \frac{6}{25} & \frac{4}{25} & \frac{2}{25} & 0 \\ \hline & \frac{1}{5} & \frac{8}{25} & \frac{6}{25} & \frac{4}{25} & \frac{2}{25} \end{array} \end{split}\]

• stages: 5

• order: 1

• stages order: 1

• stability function:

  \[ \frac{16 z^{5}}{9765625} + \frac{16 z^{4}}{78125} + \frac{28 z^{3}}{3125} + \frac{4 z^{2}}{25} + z + 1 \]

Template Parameters:

value_t – type of coefficient (doubleby default)

template<typename value_t>
using ponio::runge_kutta::rkc_52_t = explicit_runge_kutta::explicit_runge_kutta<butcher_rkc_52<value_t>>

RKC (5,2) method.

see more on ponio

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|ccccc} 0 & 0 & 0 & 0 & 0 & 0 \\ \frac{1}{8} & \frac{1}{8} & 0 & 0 & 0 & 0 \\ \frac{1}{8} & \frac{1}{16} & \frac{1}{16} & 0 & 0 & 0 \\ \frac{1}{3} & - \frac{1}{9} & \frac{4}{27} & \frac{8}{27} & 0 & 0 \\ \frac{5}{8} & - \frac{255}{512} & \frac{15}{64} & \frac{5}{8} & \frac{135}{512} & 0 \\ \hline & - \frac{269}{250} & \frac{8}{25} & \frac{24}{25} & \frac{27}{50} & \frac{32}{125} \end{array} \end{split}\]

• stages: 5

• order: 2

• stages order: 1

• stability function:

  \[ \frac{z^{5}}{6400} + \frac{z^{4}}{160} + \frac{7 z^{3}}{80} + \frac{z^{2}}{2} + z + 1 \]

Template Parameters:

value_t – type of coefficient (doubleby default)

Embedded methods

The ponio library provides also adaptive time step methods.

template<typename value_t>
using ponio::runge_kutta::rk54_6m_t = explicit_runge_kutta::explicit_runge_kutta<butcher_rk54_6m<value_t>>

RK5(4) 6M method.

see more on ponio

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|cccccc} 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ \frac{1}{5} & \frac{1}{5} & 0 & 0 & 0 & 0 & 0 \\ \frac{3}{10} & \frac{3}{40} & \frac{9}{40} & 0 & 0 & 0 & 0 \\ \frac{3}{5} & \frac{3}{10} & - \frac{9}{10} & \frac{6}{5} & 0 & 0 & 0 \\ \frac{2}{3} & \frac{226}{729} & - \frac{25}{27} & \frac{880}{729} & \frac{55}{729} & 0 & 0 \\ 1 & - \frac{181}{270} & \frac{5}{2} & - \frac{266}{297} & - \frac{91}{27} & \frac{189}{55} & 0 \\ \hline & \frac{19}{216} & 0 & \frac{1000}{2079} & - \frac{125}{216} & \frac{81}{88} & \frac{5}{56} \\ \hline & \frac{31}{540} & 0 & \frac{190}{297} & - \frac{145}{108} & \frac{351}{220} & \frac{1}{20} \end{array} \end{split}\]

• stages: 6

• order: 5

• stages order: 1

• stability function:

  \[ \frac{z^{6}}{800} + \frac{z^{5}}{120} + \frac{z^{4}}{24} + \frac{z^{3}}{6} + \frac{z^{2}}{2} + z + 1 \]

• bibliography: Dormand, J.R. & Prince, P.J., , 1980, Journal of Computational and Applied Mathematics

Template Parameters:

value_t – type of coefficient (doubleby default)

template<typename value_t>
using ponio::runge_kutta::rk54_7m_t = explicit_runge_kutta::explicit_runge_kutta<butcher_rk54_7m<value_t>>

RK5(4) 7M method.

see more on ponio

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|ccccccc} 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ \frac{1}{5} & \frac{1}{5} & 0 & 0 & 0 & 0 & 0 & 0 \\ \frac{3}{10} & \frac{3}{40} & \frac{9}{40} & 0 & 0 & 0 & 0 & 0 \\ \frac{4}{5} & \frac{44}{45} & - \frac{56}{15} & \frac{32}{9} & 0 & 0 & 0 & 0 \\ \frac{8}{9} & \frac{19372}{6561} & - \frac{25360}{2187} & \frac{64448}{6561} & - \frac{212}{729} & 0 & 0 & 0 \\ 1 & \frac{9017}{3168} & - \frac{355}{33} & \frac{46732}{5247} & \frac{49}{176} & - \frac{5103}{18656} & 0 & 0 \\ 1 & \frac{35}{384} & 0 & \frac{500}{1113} & \frac{125}{192} & - \frac{2187}{6784} & \frac{11}{84} & 0 \\ \hline & \frac{35}{384} & 0 & \frac{500}{1113} & \frac{125}{192} & - \frac{2187}{6784} & \frac{11}{84} & 0 \\ \hline & \frac{5179}{57600} & 0 & \frac{7571}{16695} & \frac{393}{640} & - \frac{92097}{339200} & \frac{187}{2100} & \frac{1}{40} \end{array} \end{split}\]

• stages: 7

• order: 5

• stages order: 1

• stability function:

  \[ \frac{z^{6}}{600} + \frac{z^{5}}{120} + \frac{z^{4}}{24} + \frac{z^{3}}{6} + \frac{z^{2}}{2} + z + 1 \]

• bibliography: Dormand, J.R. & Prince, P.J., , 1980, Journal of Computational and Applied Mathematics

Template Parameters:

value_t – type of coefficient (doubleby default)

template<typename value_t>
using ponio::runge_kutta::rk54_7s_t = explicit_runge_kutta::explicit_runge_kutta<butcher_rk54_7s<value_t>>

RK5(4) 7S method.

see more on ponio

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|ccccccc} 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ \frac{2}{9} & \frac{2}{9} & 0 & 0 & 0 & 0 & 0 & 0 \\ \frac{1}{3} & \frac{1}{12} & \frac{1}{4} & 0 & 0 & 0 & 0 & 0 \\ \frac{5}{9} & \frac{55}{324} & - \frac{25}{108} & \frac{50}{81} & 0 & 0 & 0 & 0 \\ \frac{2}{3} & \frac{83}{330} & - \frac{13}{22} & \frac{61}{66} & \frac{9}{110} & 0 & 0 & 0 \\ 1 & - \frac{19}{28} & \frac{9}{4} & \frac{1}{7} & - \frac{27}{7} & \frac{22}{7} & 0 & 0 \\ 1 & \frac{19}{200} & 0 & \frac{3}{5} & - \frac{243}{400} & \frac{33}{40} & \frac{7}{80} & 0 \\ \hline & \frac{19}{200} & 0 & \frac{3}{5} & - \frac{243}{400} & \frac{33}{40} & \frac{7}{80} & 0 \\ \hline & \frac{431}{5000} & 0 & \frac{333}{500} & - \frac{7857}{10000} & \frac{957}{1000} & \frac{193}{2000} & - \frac{1}{50} \end{array} \end{split}\]

• stages: 7

• order: 5

• stages order: 1

• stability function:

  \[ \frac{z^{6}}{1296} + \frac{z^{5}}{120} + \frac{z^{4}}{24} + \frac{z^{3}}{6} + \frac{z^{2}}{2} + z + 1 \]

• bibliography: Dormand, J.R. & Prince, P.J., , 1980, Journal of Computational and Applied Mathematics

Template Parameters:

value_t – type of coefficient (doubleby default)

template<typename value_t>
using ponio::runge_kutta::rk65_8m_t = explicit_runge_kutta::explicit_runge_kutta<butcher_rk65_8m<value_t>>

RK6(5) 8M method.

see more on ponio

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|cccccccc} 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ \frac{1}{10} & \frac{1}{10} & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ \frac{2}{9} & - \frac{2}{81} & \frac{20}{81} & 0 & 0 & 0 & 0 & 0 & 0 \\ \frac{3}{7} & \frac{615}{1372} & - \frac{270}{343} & \frac{1053}{1372} & 0 & 0 & 0 & 0 & 0 \\ \frac{3}{5} & \frac{3243}{5500} & - \frac{54}{55} & \frac{50949}{71500} & \frac{4998}{17875} & 0 & 0 & 0 & 0 \\ \frac{4}{5} & - \frac{26492}{37125} & \frac{72}{55} & \frac{2808}{23375} & - \frac{24206}{37125} & \frac{338}{459} & 0 & 0 & 0 \\ 1 & \frac{5561}{2376} & - \frac{35}{11} & - \frac{24117}{31603} & \frac{899983}{200772} & - \frac{5225}{1836} & \frac{3925}{4056} & 0 & 0 \\ 1 & \frac{465467}{266112} & - \frac{2945}{1232} & - \frac{5610201}{14158144} & \frac{10513573}{3212352} & - \frac{424325}{205632} & \frac{376225}{454272} & 0 & 0 \\ \hline & \frac{821}{10800} & 0 & \frac{19683}{71825} & \frac{175273}{912600} & \frac{395}{3672} & \frac{785}{2704} & \frac{3}{50} & 0 \\ \hline & \frac{61}{864} & 0 & \frac{98415}{321776} & \frac{16807}{146016} & \frac{1375}{7344} & \frac{1375}{5408} & - \frac{37}{1120} & \frac{1}{10} \end{array} \end{split}\]

• stages: 8

• order: 5

• stages order: 1

• stability function:

  \[ \frac{157 z^{7}}{693000} + \frac{943 z^{6}}{693000} + \frac{z^{5}}{120} + \frac{z^{4}}{24} + \frac{z^{3}}{6} + \frac{z^{2}}{2} + z + 1 \]

• bibliography: Prince, P.J. & Dormand, J.R., , 1981, Journal of Computational and Applied Mathematics

Template Parameters:

value_t – type of coefficient (doubleby default)

template<typename value_t>
using ponio::runge_kutta::rk87_13m_t = explicit_runge_kutta::explicit_runge_kutta<butcher_rk87_13m<value_t>>

RK8(7) 13M method.

see more on ponio

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|ccccccccccccc} 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ \frac{1}{18} & \frac{1}{18} & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ \frac{1}{12} & \frac{1}{48} & \frac{1}{16} & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ \frac{1}{8} & \frac{1}{32} & 0 & \frac{3}{32} & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ \frac{5}{16} & \frac{5}{16} & 0 & - \frac{75}{64} & \frac{75}{64} & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ \frac{3}{8} & \frac{3}{80} & 0 & 0 & \frac{3}{16} & \frac{3}{20} & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ \frac{59}{400} & \frac{29443841}{614563906} & 0 & 0 & \frac{77736538}{692538347} & - \frac{28693883}{1125000000} & \frac{23124283}{1800000000} & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ \frac{93}{200} & \frac{16016141}{946692911} & 0 & 0 & \frac{61564180}{158732637} & \frac{22789713}{633445777} & \frac{545815736}{2771057229} & - \frac{180193667}{1043307555} & 0 & 0 & 0 & 0 & 0 & 0 \\ \frac{5490023248}{9719169821} & \frac{39632708}{573591083} & 0 & 0 & - \frac{433636366}{683701615} & - \frac{421739975}{2616292301} & \frac{100302831}{723423059} & \frac{790204164}{839813087} & \frac{800635310}{3783071287} & 0 & 0 & 0 & 0 & 0 \\ \frac{13}{20} & \frac{246121993}{1340847787} & 0 & 0 & - \frac{37695042795}{15268766246} & - \frac{309121744}{1061227803} & - \frac{12992083}{490766935} & \frac{6005943493}{2108947869} & \frac{393006217}{1396673457} & \frac{123872331}{1001029789} & 0 & 0 & 0 & 0 \\ \frac{1201146811}{1299019798} & - \frac{1028468189}{846180014} & 0 & 0 & \frac{8478235783}{508512852} & \frac{1311729495}{1432422823} & - \frac{10304129995}{1701304382} & - \frac{48777925059}{3047939560} & \frac{15336726248}{1032824649} & - \frac{45442868181}{3398467696} & \frac{3065993473}{597172653} & 0 & 0 & 0 \\ 1 & \frac{185892177}{718116043} & 0 & 0 & - \frac{3185094517}{667107341} & - \frac{477755414}{1098053517} & - \frac{703635378}{230739211} & \frac{5731566787}{1027545527} & \frac{5232866602}{850066563} & - \frac{4093664535}{808688257} & \frac{3962137247}{1805957418} & \frac{65686358}{487910083} & 0 & 0 \\ 1 & \frac{403863854}{491063109} & 0 & 0 & - \frac{5068492393}{434740067} & - \frac{411421997}{543043805} & \frac{652783627}{914296604} & \frac{11173962825}{925320556} & - \frac{13158990841}{6184727034} & \frac{3936647629}{1978049680} & - \frac{160528059}{685178525} & \frac{248638103}{1413531060} & 0 & 0 \\ \hline & \frac{13451932}{455176623} & 0 & 0 & 0 & 0 & - \frac{808719846}{976000145} & \frac{1757004468}{5645159321} & \frac{656045339}{265891186} & - \frac{3867574721}{1518517206} & \frac{465885868}{322736535} & \frac{53011238}{667516719} & \frac{2}{45} & 0 \\ \hline & \frac{14005451}{335480064} & 0 & 0 & 0 & 0 & - \frac{59238493}{1068277825} & \frac{181606767}{758867731} & \frac{561292985}{797845732} & - \frac{1041891430}{1371343529} & \frac{760417239}{1151165299} & \frac{118820643}{751138087} & - \frac{528747749}{2220607170} & \frac{1}{4} \end{array} \end{split}\]

• stages: 13

• order: 5

• stages order: 0

• stability function:

  \[ - \frac{18014461410191061851044337424226964161393267397 z^{12}}{176362113233915443784632481296954662860466153062400000000} + \frac{579610672882944040965459831292043121517203402453769233950141611115653229058190630611561896792980377451 z^{11}}{52930149274503007565740611232642107727595980255919940377981010519187847281299356065113509881736331264000000000} + \frac{1311690294666324576022540737002853396744479072059268864101155037520235267958573549649709719841929880743320656239063971762630287492634974480323740498115819303 z^{10}}{4688900342046827221233515138559599920429581165305922099071888609522466458955325926837138275611106507974272813904493125923254971797262611160871875006955520000000000} + \frac{173202026928930674956862643792583472430250030082292305387802401908069242973872705211790880025189048547738958626947936068664047772105017256485381526909165350626899349196961435582733119709433309706069823297 z^{9}}{67105009555446125023068824548786436707601134466887950725124522889091129532097475245188150926338518916219849168487030705367949449358686050694711115285704900502867563299072552908079406625899017839509504000000000} + \frac{7715232892371949774979784590658427020425594166350073790724222667294167202182718396183370280501115600063120008672418123668744560041229464535772623192406109198524961399341137876727533206440531469493852666690918433774379091701003201795555337533344723 z^{8}}{308063994448581211308820363469127573445794487885051155265053442800141273324921936449975003043547483729993905998094764633305974019177609230957044120104890909146869771787154457113674382372738353117527574872488176875678515767874036429961250734080000000000} + \frac{1115410772446920399165305019719214725765371472850245379183763357941915347092294926234456635437733656099186960427101598549137987457933329442390392644569310693844705793492448010408280781272564584162428916553203011074648454741226478889523182041000669641484874761652737883698112559927732662571 z^{7}}{5621670293132478897734681633585456635782010950377593787681816513737078370374455354028115066931825468927808340615671300736801048788883707461876157917746358391571722654115519102224487269020296777426775904906708003430717505307068899653397091192232294293767812034172984650357204859617280000000000} + \frac{2334184210892632413432406632760118538337965198713072315996282819882547044427650573815526110550168509804990463654159134528362052293639226380264133101339961105460479644073451098753650308153072999048686027975471718977265895400336367229277347300063482314939584551590892019533087498395484565799320354149668843956037 z^{6}}{1680612631842695365419782967854313884096448524191177364353254872846657435518138538104378993769154867373255541490261239744797413240336864205184799981310002908344334745670117778934423186016602134782127099076328585467665019379670118441448699823611353002416723647861324475407438047381724399955314181034188800000000000} + \frac{101938791853692372195797298931657797035446315257702626996539029237170275359971335430427362704404014005955811485110901810821115594835218928148422697586247351886417158308973747725084896276348838296766115956166200620137082882467258262588569033910706668118005115419413632925651871658299708165056330850217590331728453024615360637675899 z^{5}}{12232655022443084867186848515744314278060261084231049838318011538894205145990936296390374381515621795156936523017935899423481812590503754314178625220100007212342174547711683064105274245610997023365929547084040244234165614715135707927688581300109712475395341921292807529093266366872602285862829615310416554653646160845529600000000000} + \frac{138820097270640751120649850690533109535436253518333598268111135955940008885783475477190314584157308544073572073500149585900701977685573813193247885201454753963496045709198563626005911211296652210100215694168128651787914370606904284332938258325794036365064828248882400532748197533758274819705163326831259697511341082278673947921917179437873 z^{4}}{3331682334495378028266699558154999686693825757134914203767818382252623002541210392436705070441207772961670211088202513385121752112616157038057116570905551150962337027779221405580658896786293713771395640894876869351627596338589008467714916771851280949811277658571238421987834003438977586204595534034553650306587157519688552641233400000000000} + \frac{3216737384368291044698829700110298764369664287360835592258622015286755792022434243736588442844033580530423099367532361370615245170560337176132882900424940350802558077478589846160796918495618861855441103784076715411009578045816743811248521650498723955189031730729100471151548928419603476636427196447952583568864143090650768677307983523126179 z^{3}}{19300424306209746327953771365225335142148564894753303271242687236259452113118809358360232773977835790060757284835517651354864672562460022508282823989970241522919638048352573623033815278023599290575127094633593342804859458902451341482245913263407288839790392247464745045036013070971092101256021450813070321717626259367871188699629220000000000} + \frac{54736314961905438367044618304357346820367656411141526630641353828331284274482492355150836021591733858246532000654016442241006910581380001863585052317000771310994173682808059050685863646266334772248540751685793756206483101990512959326950535438526941201674760403119908113929728847761616112749381005775214952273840810844959153843281976167508849 z^{2}}{109472629923810877303151802180443887829663773452582787281691697382229123444715213354774726603313486065343651849174748730456837655237271123163198252479200918706478090870824248798633506615516456844589767154259107323781280356251508456322869018913086115015505821109371806986065156372896859837999835823419158261983104914773298644450774170750000000} + \frac{370952526618963512497086801444530779494344988945804575254 z}{370952526618963512211046321800629321994711078958870761915} + 1 \]

• bibliography: Prince, P.J. & Dormand, J.R., , 1981, Journal of Computational and Applied Mathematics

Template Parameters:

value_t – type of coefficient (doubleby default)

Diagonal implicit methods

When the matrix \(A\) is lower triangular with a diagonal, the Runge-Kutta method is called diagonal implicit (or DIRK). You have to provide a Jacobian function that returns the Jacobian matrix in point \((t, u)\) (see ponio::implicit_problem). You can also provide an operator base definition (see ponio::implicit_operator_problem).

template<typename value_t, typename linear_algebra_t = void, typename ...Args>
auto ponio::runge_kutta::backward_euler_t(Args... args)

backward Euler method

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|c} 1 & 1 \\ \hline & 1 \end{array} \end{split}\]

• stages: 1

• order: 1

• stages order: 1

• stability function:

  \[ - \frac{1}{z - 1} \]

Template Parameters:

• value_t – type of coefficient (doubleby default)

• linear_algebra_t – type that provides linear algebra if it is undefined for state_t (void by default)

• Args – optional arguments to build linear_algebra_t object

template<typename value_t, typename linear_algebra_t = void, typename ...Args>
auto ponio::runge_kutta::crancknicolson_t(Args... args)

Cranck-Nicolson method.

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|cc} 0 & 0 & 0 \\ 1 & \frac{1}{2} & \frac{1}{2} \\ \hline & \frac{1}{2} & \frac{1}{2} \end{array} \end{split}\]

• stages: 2

• order: 2

• stages order: 2

• stability function:

  \[ \frac{- z - 2}{z - 2} \]

Template Parameters:

• value_t – type of coefficient (doubleby default)

• linear_algebra_t – type that provides linear algebra if it is undefined for state_t (void by default)

• Args – optional arguments to build linear_algebra_t object

template<typename value_t, typename linear_algebra_t = void, typename ...Args>
auto ponio::runge_kutta::crancknicolson_2_t(Args... args)

Cranck-Nicolson 2 method.

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|cc} 0 & 0 & 0 \\ 1 & 0 & 1 \\ \hline & \frac{1}{2} & \frac{1}{2} \end{array} \end{split}\]

• stages: 2

• order: 2

• stages order: 1

• stability function:

  \[ \frac{z^{2} - 2}{2 \left(z - 1\right)} \]

Template Parameters:

• value_t – type of coefficient (doubleby default)

• linear_algebra_t – type that provides linear algebra if it is undefined for state_t (void by default)

• Args – optional arguments to build linear_algebra_t object

template<typename value_t, typename linear_algebra_t = void, typename ...Args>
auto ponio::runge_kutta::dirk23_t(Args... args)

DIRK(2,3) method.

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|cc} \frac{\sqrt{3}}{6} + \frac{1}{2} & \frac{\sqrt{3}}{6} + \frac{1}{2} & 0 \\ \frac{1}{2} - \frac{\sqrt{3}}{6} & - \frac{\sqrt{3}}{3} & \frac{\sqrt{3}}{6} + \frac{1}{2} \\ \hline & \frac{1}{2} & \frac{1}{2} \end{array} \end{split}\]

• stages: 2

• order: 3

• stages order: 1

• stability function:

  \[ \frac{z^{2} \left(- \sqrt{3} - 1\right) - 2 \sqrt{3} z + 6}{z^{2} \left(\sqrt{3} + 2\right) + z \left(-6 - 2 \sqrt{3}\right) + 6} \]

Template Parameters:

• value_t – type of coefficient (doubleby default)

• linear_algebra_t – type that provides linear algebra if it is undefined for state_t (void by default)

• Args – optional arguments to build linear_algebra_t object

template<typename value_t, typename linear_algebra_t = void, typename ...Args>
auto ponio::runge_kutta::dirk23_crouzeix_t(Args... args)

DIRK(2,3) Crouzeix method.

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|cc} \frac{\sqrt{3}}{6} + \frac{1}{2} & \frac{\sqrt{3}}{6} + \frac{1}{2} & 0 \\ \frac{1}{2} - \frac{\sqrt{3}}{6} & - \frac{\sqrt{3}}{3} & \frac{\sqrt{3}}{6} + \frac{1}{2} \\ \hline & \frac{1}{2} & \frac{1}{2} \end{array} \end{split}\]

• stages: 2

• order: 3

• stages order: 1

• stability function:

  \[ \frac{z^{2} \left(- \sqrt{3} - 1\right) - 2 \sqrt{3} z + 6}{z^{2} \left(\sqrt{3} + 2\right) + z \left(-6 - 2 \sqrt{3}\right) + 6} \]

Template Parameters:

• value_t – type of coefficient (doubleby default)

• linear_algebra_t – type that provides linear algebra if it is undefined for state_t (void by default)

• Args – optional arguments to build linear_algebra_t object

template<typename value_t, typename linear_algebra_t = void, typename ...Args>
auto ponio::runge_kutta::dirk34_t(Args... args)

DIRK(3,4) method.

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|ccc} \frac{1}{2} + \frac{\sqrt{3} \cos{\left(\frac{\pi}{18} \right)}}{3} & \frac{1}{2} + \frac{\sqrt{3} \cos{\left(\frac{\pi}{18} \right)}}{3} & 0 & 0 \\ \frac{1}{2} & - \frac{\sqrt{3} \cos{\left(\frac{\pi}{18} \right)}}{3} & \frac{1}{2} + \frac{\sqrt{3} \cos{\left(\frac{\pi}{18} \right)}}{3} & 0 \\ - \frac{\sqrt{3} \cos{\left(\frac{\pi}{18} \right)}}{3} + \frac{1}{2} & 1 + \frac{2 \sqrt{3} \cos{\left(\frac{\pi}{18} \right)}}{3} & - \frac{4 \sqrt{3} \cos{\left(\frac{\pi}{18} \right)}}{3} - 1 & \frac{1}{2} + \frac{\sqrt{3} \cos{\left(\frac{\pi}{18} \right)}}{3} \\ \hline & \frac{1}{8 \cos^{2}{\left(\frac{\pi}{18} \right)}} & 1 - \frac{1}{4 \cos^{2}{\left(\frac{\pi}{18} \right)}} & \frac{1}{8 \cos^{2}{\left(\frac{\pi}{18} \right)}} \end{array} \end{split}\]

• stages: 3

• order: 4

• stages order: 1

• stability function:

  \[ \frac{z^{3} \left(- 30 \cos{\left(\frac{\pi}{18} \right)} - \frac{27 \sqrt{3}}{2} - \sqrt{3} \cos{\left(\frac{\pi}{9} \right)} + 2 \sqrt{3} \left(1 - \cos{\left(\frac{\pi}{9} \right)}\right)^{2}\right) + z^{2} \left(- 36 \cos{\left(\frac{\pi}{18} \right)} - 9 \sqrt{3}\right) + z \left(36 \cos{\left(\frac{\pi}{18} \right)} + 36 \sqrt{3} \cos{\left(\frac{\pi}{9} \right)} + 36 \sqrt{3}\right) - 72 \cos{\left(\frac{\pi}{18} \right)}}{\left(z^{3} \left(18 \cos{\left(\frac{\pi}{9} \right)} + 30 + 24 \sqrt{3} \cos{\left(\frac{\pi}{18} \right)}\right) + z^{2} \left(- 72 \sqrt{3} \cos{\left(\frac{\pi}{18} \right)} - 90 - 36 \cos{\left(\frac{\pi}{9} \right)}\right) + z \left(108 + 72 \sqrt{3} \cos{\left(\frac{\pi}{18} \right)}\right) - 72\right) \cos{\left(\frac{\pi}{18} \right)}} \]

Template Parameters:

• value_t – type of coefficient (doubleby default)

• linear_algebra_t – type that provides linear algebra if it is undefined for state_t (void by default)

• Args – optional arguments to build linear_algebra_t object

template<typename value_t, typename linear_algebra_t = void, typename ...Args>
auto ponio::runge_kutta::dirk_qin_zhang_t(Args... args)

DIRK Qin Zhang method.

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|cc} \frac{1}{4} & \frac{1}{4} & 0 \\ \frac{3}{4} & \frac{1}{2} & \frac{1}{4} \\ \hline & \frac{1}{2} & \frac{1}{2} \end{array} \end{split}\]

• stages: 2

• order: 2

• stages order: 1

• stability function:

  \[ \frac{z^{2} + 8 z + 16}{z^{2} - 8 z + 16} \]

Template Parameters:

• value_t – type of coefficient (doubleby default)

• linear_algebra_t – type that provides linear algebra if it is undefined for state_t (void by default)

• Args – optional arguments to build linear_algebra_t object

template<typename value_t, typename linear_algebra_t = void, typename ...Args>
auto ponio::runge_kutta::esdirk_54_a_t(Args... args)

ESDIRK (5,4) a method.

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|ccccccc} 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ \frac{13}{25} & \frac{13}{50} & \frac{13}{50} & 0 & 0 & 0 & 0 & 0 \\ \frac{123033320996791}{100000000000000} & \frac{13}{100} & \frac{210083302491977}{250000000000000} & \frac{13}{50} & 0 & 0 & 0 & 0 \\ \frac{223941496087519}{250000000000000} & \frac{44743922956641}{200000000000000} & \frac{476755323197997}{1000000000000000} & - \frac{323544768155631}{5000000000000000} & \frac{13}{50} & 0 & 0 & 0 \\ \frac{54549201232331}{125000000000000} & \frac{166485643232483}{1000000000000000} & \frac{104500188415917}{1000000000000000} & \frac{363148227209871}{10000000000000000} & - \frac{6545352225537}{50000000000000} & \frac{13}{50} & 0 & 0 \\ 1 & \frac{69278201156341}{500000000000000} & 0 & - \frac{106133430043801}{2500000000000000} & \frac{122332894900157}{5000000000000000} & \frac{619430390724807}{1000000000000000} & \frac{13}{50} & 0 \\ 1 & \frac{136597511776403}{1000000000000000} & 0 & - \frac{274845439826919}{5000000000000000} & - \frac{82372534566421}{2000000000000000} & \frac{157483262247541}{250000000000000} & \frac{696247944820273}{10000000000000000} & \frac{13}{50} \\ \hline & \frac{136597511776403}{1000000000000000} & 0 & - \frac{274845439826919}{5000000000000000} & - \frac{82372534566421}{2000000000000000} & \frac{157483262247541}{250000000000000} & \frac{696247944820273}{10000000000000000} & \frac{13}{50} \\ \hline & \frac{69278201156341}{500000000000000} & 0 & - \frac{106133430043801}{2500000000000000} & \frac{122332894900157}{5000000000000000} & \frac{619430390724807}{1000000000000000} & \frac{13}{50} & 0 \end{array} \end{split}\]

• stages: 7

• order: 5

• stages order: 0

• stability function:

  \[ \frac{- 16398976328036028934688055215985721008498119004548655437683 z^{6} - 62784133333333052867655311368879985464447647030969785505495979417905529550 z^{5} + 177908333333331334958183544805291543462556056513069659813105250000000000000 z^{4} + 1535833333333341158027214497462647329913534272687500000000000000000000000000 z^{3} - 1437500000000011502293449947559375000000000000000000000000000000000000000000 z^{2} - 17500000000000000000000000000000000000000000000000000000000000000000000000000 z + 31250000000000000000000000000000000000000000000000000000000000000000000000000}{2000000000000000000000000000000000000000000000000000000000000000000 \left(4826809 z^{6} - 111387900 z^{5} + 1071037500 z^{4} - 5492500000 z^{3} + 15843750000 z^{2} - 24375000000 z + 15625000000\right)} \]

• bibliography: Kværnø, A., , 2004, BIT Numerical Mathematics

Template Parameters:

• value_t – type of coefficient (doubleby default)

• linear_algebra_t – type that provides linear algebra if it is undefined for state_t (void by default)

• Args – optional arguments to build linear_algebra_t object

template<typename value_t, typename linear_algebra_t = void, typename ...Args>
auto ponio::runge_kutta::implicit_midpoint_t(Args... args)

implicit midpoint method

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|c} \frac{1}{2} & \frac{1}{2} \\ \hline & 1 \end{array} \end{split}\]

• stages: 1

• order: 2

• stages order: 1

• stability function:

  \[ \frac{- z - 2}{z - 2} \]

Template Parameters:

• value_t – type of coefficient (doubleby default)

• linear_algebra_t – type that provides linear algebra if it is undefined for state_t (void by default)

• Args – optional arguments to build linear_algebra_t object

template<typename value_t, typename linear_algebra_t = void, typename ...Args>
auto ponio::runge_kutta::lsdirk22qz_t(Args... args)

L-SDIRK-22-QZ method.

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|cc} \frac{\sqrt{2}}{2} + 1 & \frac{\sqrt{2}}{2} + 1 & 0 \\ 1 & - \frac{\sqrt{2}}{2} & \frac{\sqrt{2}}{2} + 1 \\ \hline & \frac{1}{2} & \frac{1}{2} \end{array} \end{split}\]

• stages: 2

• order: 1

• stages order: 1

• stability function:

  \[ \frac{z^{2} \left(\sqrt{2} + 2\right) + z \left(- 4 \sqrt{2} - 4\right) + 4}{2 \left(z^{2} \left(2 \sqrt{2} + 3\right) + z \left(-4 - 2 \sqrt{2}\right) + 2\right)} \]

Template Parameters:

• value_t – type of coefficient (doubleby default)

• linear_algebra_t – type that provides linear algebra if it is undefined for state_t (void by default)

• Args – optional arguments to build linear_algebra_t object

template<typename value_t, typename linear_algebra_t = void, typename ...Args>
auto ponio::runge_kutta::lsdirk33_t(Args... args)

L-SDIRK-33 method.

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|ccc} \frac{871733043}{2000000000} & \frac{871733043}{2000000000} & 0 & 0 \\ \frac{2871733043}{4000000000} & \frac{1128266957}{4000000000} & \frac{871733043}{2000000000} & 0 \\ 1 & \frac{120849664915323}{100000000000000} & - \frac{128872634130647}{200000000000000} & \frac{871733043}{2000000000} \\ \hline & \frac{120849664915323}{100000000000000} & - \frac{128872634130647}{200000000000000} & \frac{871733043}{2000000000} \end{array} \end{split}\]

• stages: 3

• order: 3

• stages order: 0

• stability function:

  \[ \frac{10000 \left(190128552647941748065093 z^{2} + 246079651600004000000000 z - 800000000000000000000000\right)}{662446064918471276784030507 z^{3} - 4559510989548239094000000000 z^{2} + 10460796516000000000000000000 z - 8000000000000000000000000000} \]

Template Parameters:

• value_t – type of coefficient (doubleby default)

• linear_algebra_t – type that provides linear algebra if it is undefined for state_t (void by default)

• Args – optional arguments to build linear_algebra_t object

template<typename value_t, typename linear_algebra_t = void, typename ...Args>
auto ponio::runge_kutta::lsdirk43_t(Args... args)

L-SDIRK-43 method.

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|cccc} \frac{1}{2} & \frac{1}{2} & 0 & 0 & 0 \\ \frac{2}{3} & \frac{1}{6} & \frac{1}{2} & 0 & 0 \\ \frac{1}{2} & - \frac{1}{2} & \frac{1}{2} & \frac{1}{2} & 0 \\ 1 & \frac{3}{2} & - \frac{3}{2} & \frac{1}{2} & \frac{1}{2} \\ \hline & \frac{3}{2} & - \frac{3}{2} & \frac{1}{2} & \frac{1}{2} \end{array} \end{split}\]

• stages: 4

• order: 3

• stages order: 1

• stability function:

  \[ \frac{8 \left(z^{3} - 6 z + 6\right)}{3 \left(z^{4} - 8 z^{3} + 24 z^{2} - 32 z + 16\right)} \]

Template Parameters:

• value_t – type of coefficient (doubleby default)

• linear_algebra_t – type that provides linear algebra if it is undefined for state_t (void by default)

• Args – optional arguments to build linear_algebra_t object

template<typename value_t, typename linear_algebra_t = void, typename ...Args>
auto ponio::runge_kutta::sdirk_23_t(Args... args)

SDIRK (2,3) method.

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|cc} \frac{\sqrt{3}}{6} + \frac{1}{2} & \frac{\sqrt{3}}{6} + \frac{1}{2} & 0 \\ \frac{1}{2} - \frac{\sqrt{3}}{6} & - \frac{\sqrt{3}}{3} & \frac{\sqrt{3}}{6} + \frac{1}{2} \\ \hline & \frac{1}{2} & \frac{1}{2} \end{array} \end{split}\]

• stages: 2

• order: 3

• stages order: 1

• stability function:

  \[ \frac{z^{2} \left(- \sqrt{3} - 1\right) - 2 \sqrt{3} z + 6}{z^{2} \left(\sqrt{3} + 2\right) + z \left(-6 - 2 \sqrt{3}\right) + 6} \]

Template Parameters:

• value_t – type of coefficient (doubleby default)

• linear_algebra_t – type that provides linear algebra if it is undefined for state_t (void by default)

• Args – optional arguments to build linear_algebra_t object

template<typename value_t, typename linear_algebra_t = void, typename ...Args>
auto ponio::runge_kutta::sdirk_34_t(Args... args)

SDIRK (3,4) method.

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|ccc} \frac{1}{2} + \frac{\sqrt{3} \cos{\left(\frac{\pi}{18} \right)}}{3} & \frac{1}{2} + \frac{\sqrt{3} \cos{\left(\frac{\pi}{18} \right)}}{3} & 0 & 0 \\ \frac{1}{2} & - \frac{\sqrt{3} \cos{\left(\frac{\pi}{18} \right)}}{3} & \frac{1}{2} + \frac{\sqrt{3} \cos{\left(\frac{\pi}{18} \right)}}{3} & 0 \\ - \frac{\sqrt{3} \cos{\left(\frac{\pi}{18} \right)}}{3} + \frac{1}{2} & 1 + \frac{2 \sqrt{3} \cos{\left(\frac{\pi}{18} \right)}}{3} & - \frac{4 \sqrt{3} \cos{\left(\frac{\pi}{18} \right)}}{3} - 1 & \frac{1}{2} + \frac{\sqrt{3} \cos{\left(\frac{\pi}{18} \right)}}{3} \\ \hline & \frac{1}{8 \cos^{2}{\left(\frac{\pi}{18} \right)}} & 1 - \frac{1}{4 \cos^{2}{\left(\frac{\pi}{18} \right)}} & \frac{1}{8 \cos^{2}{\left(\frac{\pi}{18} \right)}} \end{array} \end{split}\]

• stages: 3

• order: 4

• stages order: 1

• stability function:

  \[ \frac{z^{3} \left(- 30 \cos{\left(\frac{\pi}{18} \right)} - \frac{27 \sqrt{3}}{2} - \sqrt{3} \cos{\left(\frac{\pi}{9} \right)} + 2 \sqrt{3} \left(1 - \cos{\left(\frac{\pi}{9} \right)}\right)^{2}\right) + z^{2} \left(- 36 \cos{\left(\frac{\pi}{18} \right)} - 9 \sqrt{3}\right) + z \left(36 \cos{\left(\frac{\pi}{18} \right)} + 36 \sqrt{3} \cos{\left(\frac{\pi}{9} \right)} + 36 \sqrt{3}\right) - 72 \cos{\left(\frac{\pi}{18} \right)}}{\left(z^{3} \left(18 \cos{\left(\frac{\pi}{9} \right)} + 30 + 24 \sqrt{3} \cos{\left(\frac{\pi}{18} \right)}\right) + z^{2} \left(- 72 \sqrt{3} \cos{\left(\frac{\pi}{18} \right)} - 90 - 36 \cos{\left(\frac{\pi}{9} \right)}\right) + z \left(108 + 72 \sqrt{3} \cos{\left(\frac{\pi}{18} \right)}\right) - 72\right) \cos{\left(\frac{\pi}{18} \right)}} \]

Template Parameters:

• value_t – type of coefficient (doubleby default)

• linear_algebra_t – type that provides linear algebra if it is undefined for state_t (void by default)

• Args – optional arguments to build linear_algebra_t object

template<typename value_t, typename linear_algebra_t = void, typename ...Args>
auto ponio::runge_kutta::sdirk_54_t(Args... args)

SDIRK (5,4) method.

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|ccccc} \frac{1}{4} & \frac{1}{4} & 0 & 0 & 0 & 0 \\ \frac{3}{4} & \frac{1}{2} & \frac{1}{4} & 0 & 0 & 0 \\ \frac{11}{20} & \frac{17}{50} & - \frac{1}{25} & \frac{1}{4} & 0 & 0 \\ \frac{1}{2} & \frac{371}{1360} & - \frac{137}{2720} & \frac{15}{544} & \frac{1}{4} & 0 \\ 1 & \frac{25}{24} & - \frac{49}{48} & \frac{125}{16} & - \frac{85}{12} & \frac{1}{4} \\ \hline & \frac{25}{24} & - \frac{49}{48} & \frac{125}{16} & - \frac{85}{12} & \frac{1}{4} \\ \hline & \frac{59}{48} & - \frac{17}{96} & \frac{225}{32} & - \frac{85}{12} & 0 \end{array} \end{split}\]

• stages: 5

• order: 4

• stages order: 1

• stability function:

  \[ \frac{4 \left(- 7 z^{4} - 8 z^{3} + 96 z^{2} + 192 z - 768\right)}{3 \left(z^{5} - 20 z^{4} + 160 z^{3} - 640 z^{2} + 1280 z - 1024\right)} \]

• bibliography: Hairer, Ernst & Wanner, Gerhard, , 1996, Springer Berlin Heidelberg

Template Parameters:

• value_t – type of coefficient (doubleby default)

• linear_algebra_t – type that provides linear algebra if it is undefined for state_t (void by default)

• Args – optional arguments to build linear_algebra_t object

template<typename value_t, typename linear_algebra_t = void, typename ...Args>
auto ponio::runge_kutta::sspirk33_t(Args... args)

SSPIRK33 method.

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|ccc} \frac{1}{2} - \frac{\sqrt{2}}{4} & \frac{1}{2} - \frac{\sqrt{2}}{4} & 0 & 0 \\ \frac{1}{2} & \frac{\sqrt{2}}{4} & \frac{1}{2} - \frac{\sqrt{2}}{4} & 0 \\ \frac{\sqrt{2}}{4} + \frac{1}{2} & \frac{\sqrt{2}}{4} & \frac{\sqrt{2}}{4} & \frac{1}{2} - \frac{\sqrt{2}}{4} \\ \hline & \frac{1}{3} & \frac{1}{3} & \frac{1}{3} \end{array} \end{split}\]

• stages: 3

• order: 3

• stages order: 1

• stability function:

  \[ \frac{z^{3} \left(22 - 15 \sqrt{2}\right) + 12 z^{2} + z \left(-48 + 72 \sqrt{2}\right) + 96}{3 \left(z^{3} \left(-10 + 7 \sqrt{2}\right) + z^{2} \left(36 - 24 \sqrt{2}\right) + z \left(-48 + 24 \sqrt{2}\right) + 32\right)} \]

Template Parameters:

• value_t – type of coefficient (doubleby default)

• linear_algebra_t – type that provides linear algebra if it is undefined for state_t (void by default)

• Args – optional arguments to build linear_algebra_t object

Lawson methods

Lawson methods was initialy propose into [Law67]. Lawson methods are build to solve a problem with a linear and nonlinear part, to solve exactly the problem when the nonlinear part goes to zero. This class of problem can be write as

\[\dot{u} = L u + N(t, u)\]

First, we introduce the change of variable

\[v(t) = e^{-Lt}u(t)\]

which yields the equation

\[\dot{v}(t) = -Le^{-Lt}u(t) + e^{-Lt}\dot{u}(t)\]

which can be rewrite in therm of \(v\) as

\[\dot{v}(t) = e^{-Lt}N(t, e^{-Lt}v) = \tilde{N}(t, v)\]

We solve this equation with a classical Runge-Kutta method RK(\(s\), \(n\)) with \(s\) stages and of order \(n\). We rewrite the scheme in terms of the variable \(u\), which yields the following scheme

\[\begin{split}\begin{aligned} u^{(i)} &= u^n + \Delta t \sum_j a_{ij} k_j, \quad i = 1, \dots, s \\ k_i &= e^{-c_i\Delta t L} N(t^n + c_i\Delta t, e^{c_i\Delta t L} u^{(i)}) \\ u^{n+1} &= e^{\Delta t L}\left( u^n + \Delta t \sum_i b_i k_i \right) \end{aligned}\end{split}\]

In ponio, Lawson methods have the same name of the underlying Runge-Kutta method prefixed by l.

Explicit methods

template<typename value_t, typename Exp_t>
auto ponio::runge_kutta::leuler_t = []( Exp_t exp_ , double tol=ponio::default_config::tol ){returnlawson_runge_kutta::make_lawson<butcher_euler<value_t>,Exp_t>(exp_,tol);}

Lawson Euler method.

see more on ponio

• stages: 1

• order: 1

• underlying Runge-Kutta method: Euler (for more information see ponio::runge_kutta::euler_t)

Template Parameters:

value_t – type of coefficient (doubleby default)

Param exp_:

exponential function for the Lawson method

Param tol:

tolerance (only for adaptive time step methods)

template<typename value_t, typename Exp_t>
auto ponio::runge_kutta::lexplicit_euler_sub4_t = []( Exp_t exp_ , double tol=ponio::default_config::tol ){returnlawson_runge_kutta::make_lawson<butcher_explicit_euler_sub4<value_t>,Exp_t>(exp_,tol);}

Lawson Explicit Euler sub4 method.

see more on ponio

• stages: 4

• order: 1

• underlying Runge-Kutta method: Explicit Euler sub4 (for more information see ponio::runge_kutta::explicit_euler_sub4_t)

Template Parameters:

value_t – type of coefficient (doubleby default)

Param exp_:

exponential function for the Lawson method

Param tol:

tolerance (only for adaptive time step methods)

template<typename value_t, typename Exp_t>
auto ponio::runge_kutta::lrk56_t = []( Exp_t exp_ , double tol=ponio::default_config::tol ){returnlawson_runge_kutta::make_lawson<butcher_rk56<value_t>,Exp_t>(exp_,tol);}

Lawson RK56 method.

see more on ponio

• stages: 6

• order: 5

• underlying Runge-Kutta method: RK56 (for more information see ponio::runge_kutta::rk56_t)

Template Parameters:

value_t – type of coefficient (doubleby default)

Param exp_:

exponential function for the Lawson method

Param tol:

tolerance (only for adaptive time step methods)

template<typename value_t, typename Exp_t>
auto ponio::runge_kutta::lrk56a_t = []( Exp_t exp_ , double tol=ponio::default_config::tol ){returnlawson_runge_kutta::make_lawson<butcher_rk56a<value_t>,Exp_t>(exp_,tol);}

Lawson RK56a method.

see more on ponio

• stages: 6

• order: 5

• underlying Runge-Kutta method: RK56a (for more information see ponio::runge_kutta::rk56a_t)

Template Parameters:

value_t – type of coefficient (doubleby default)

Param exp_:

exponential function for the Lawson method

Param tol:

tolerance (only for adaptive time step methods)

template<typename value_t, typename Exp_t>
auto ponio::runge_kutta::lrk6es_t = []( Exp_t exp_ , double tol=ponio::default_config::tol ){returnlawson_runge_kutta::make_lawson<butcher_rk6es<value_t>,Exp_t>(exp_,tol);}

Lawson RK6ES method.

see more on ponio

• stages: 7

• order: 6

• underlying Runge-Kutta method: RK6ES (for more information see ponio::runge_kutta::rk6es_t)

Template Parameters:

value_t – type of coefficient (doubleby default)

Param exp_:

exponential function for the Lawson method

Param tol:

tolerance (only for adaptive time step methods)

template<typename value_t, typename Exp_t>
auto ponio::runge_kutta::lrk_118_t = []( Exp_t exp_ , double tol=ponio::default_config::tol ){returnlawson_runge_kutta::make_lawson<butcher_rk_118<value_t>,Exp_t>(exp_,tol);}

Lawson RK (11,8) method.

see more on ponio

• stages: 11

• order: 8

• underlying Runge-Kutta method: RK (11,8) (for more information see ponio::runge_kutta::rk_118_t)

Template Parameters:

value_t – type of coefficient (doubleby default)

Param exp_:

exponential function for the Lawson method

Param tol:

tolerance (only for adaptive time step methods)

template<typename value_t, typename Exp_t>
auto ponio::runge_kutta::lrk_21_t = []( Exp_t exp_ , double tol=ponio::default_config::tol ){returnlawson_runge_kutta::make_lawson<butcher_rk_21<value_t>,Exp_t>(exp_,tol);}

Lawson RK (2,1) method.

see more on ponio

• stages: 2

• order: 2

• underlying Runge-Kutta method: RK (2,1) (for more information see ponio::runge_kutta::rk_21_t)

Template Parameters:

value_t – type of coefficient (doubleby default)

Param exp_:

exponential function for the Lawson method

Param tol:

tolerance (only for adaptive time step methods)

template<typename value_t, typename Exp_t>
auto ponio::runge_kutta::lrk_22_midpoint_t = []( Exp_t exp_ , double tol=ponio::default_config::tol ){returnlawson_runge_kutta::make_lawson<butcher_rk_22_midpoint<value_t>,Exp_t>(exp_,tol);}

Lawson RK (2,2) mid-point method.

see more on ponio

• stages: 2

• order: 2

• underlying Runge-Kutta method: RK (2,2) mid-point (for more information see ponio::runge_kutta::rk_22_midpoint_t)

Template Parameters:

value_t – type of coefficient (doubleby default)

Param exp_:

exponential function for the Lawson method

Param tol:

tolerance (only for adaptive time step methods)

template<typename value_t, typename Exp_t>
auto ponio::runge_kutta::lrk_22_ralston_t = []( Exp_t exp_ , double tol=ponio::default_config::tol ){returnlawson_runge_kutta::make_lawson<butcher_rk_22_ralston<value_t>,Exp_t>(exp_,tol);}

Lawson RK (2,2) Ralston method.

see more on ponio

• stages: 2

• order: 2

• underlying Runge-Kutta method: RK (2,2) Ralston (for more information see ponio::runge_kutta::rk_22_ralston_t)

Template Parameters:

value_t – type of coefficient (doubleby default)

Param exp_:

exponential function for the Lawson method

Param tol:

tolerance (only for adaptive time step methods)

template<typename value_t, typename Exp_t>
auto ponio::runge_kutta::lrk_32_best_t = []( Exp_t exp_ , double tol=ponio::default_config::tol ){returnlawson_runge_kutta::make_lawson<butcher_rk_32_best<value_t>,Exp_t>(exp_,tol);}

Lawson RK (3,2) best method.

see more on ponio

• stages: 3

• order: 2

• underlying Runge-Kutta method: RK (3,2) best (for more information see ponio::runge_kutta::rk_32_best_t)

Template Parameters:

value_t – type of coefficient (doubleby default)

Param exp_:

exponential function for the Lawson method

Param tol:

tolerance (only for adaptive time step methods)

template<typename value_t, typename Exp_t>
auto ponio::runge_kutta::lrk_33_t = []( Exp_t exp_ , double tol=ponio::default_config::tol ){returnlawson_runge_kutta::make_lawson<butcher_rk_33<value_t>,Exp_t>(exp_,tol);}

Lawson RK (3,3) method.

see more on ponio

• stages: 3

• order: 3

• underlying Runge-Kutta method: RK (3,3) (for more information see ponio::runge_kutta::rk_33_t)

Template Parameters:

value_t – type of coefficient (doubleby default)

Param exp_:

exponential function for the Lawson method

Param tol:

tolerance (only for adaptive time step methods)

template<typename value_t, typename Exp_t>
auto ponio::runge_kutta::lrk_33_233e_t = []( Exp_t exp_ , double tol=ponio::default_config::tol ){returnlawson_runge_kutta::make_lawson<butcher_rk_33_233e<value_t>,Exp_t>(exp_,tol);}

Lawson RK (3,3) 233e method.

see more on ponio

• stages: 3

• order: 3

• underlying Runge-Kutta method: RK (3,3) 233e (for more information see ponio::runge_kutta::rk_33_233e_t)

Template Parameters:

value_t – type of coefficient (doubleby default)

Param exp_:

exponential function for the Lawson method

Param tol:

tolerance (only for adaptive time step methods)

template<typename value_t, typename Exp_t>
auto ponio::runge_kutta::lrk_33_bogackishampine_t = []( Exp_t exp_ , double tol=ponio::default_config::tol ){returnlawson_runge_kutta::make_lawson<butcher_rk_33_bogackishampine<value_t>,Exp_t>(exp_,tol);}

Lawson RK (3,3) Bogacki-Shampine method.

see more on ponio

• stages: 3

• order: 3

• underlying Runge-Kutta method: RK (3,3) Bogacki-Shampine (for more information see ponio::runge_kutta::rk_33_bogackishampine_t)

Template Parameters:

value_t – type of coefficient (doubleby default)

Param exp_:

exponential function for the Lawson method

Param tol:

tolerance (only for adaptive time step methods)

template<typename value_t, typename Exp_t>
auto ponio::runge_kutta::lrk_33_heun_t = []( Exp_t exp_ , double tol=ponio::default_config::tol ){returnlawson_runge_kutta::make_lawson<butcher_rk_33_heun<value_t>,Exp_t>(exp_,tol);}

Lawson RK (3,3) Heun method.

see more on ponio

• stages: 3

• order: 3

• underlying Runge-Kutta method: RK (3,3) Heun (for more information see ponio::runge_kutta::rk_33_heun_t)

Template Parameters:

value_t – type of coefficient (doubleby default)

Param exp_:

exponential function for the Lawson method

Param tol:

tolerance (only for adaptive time step methods)

template<typename value_t, typename Exp_t>
auto ponio::runge_kutta::lrk_33_ralston_t = []( Exp_t exp_ , double tol=ponio::default_config::tol ){returnlawson_runge_kutta::make_lawson<butcher_rk_33_ralston<value_t>,Exp_t>(exp_,tol);}

Lawson RK (3,3) Ralston method.

see more on ponio

• stages: 3

• order: 3

• underlying Runge-Kutta method: RK (3,3) Ralston (for more information see ponio::runge_kutta::rk_33_ralston_t)

Template Parameters:

value_t – type of coefficient (doubleby default)

Param exp_:

exponential function for the Lawson method

Param tol:

tolerance (only for adaptive time step methods)

template<typename value_t, typename Exp_t>
auto ponio::runge_kutta::lrk_33_van_der_houwen_t = []( Exp_t exp_ , double tol=ponio::default_config::tol ){returnlawson_runge_kutta::make_lawson<butcher_rk_33_van_der_houwen<value_t>,Exp_t>(exp_,tol);}

Lawson RK (3,3) Van der Houwen method.

see more on ponio

• stages: 3

• order: 3

• underlying Runge-Kutta method: RK (3,3) Van der Houwen (for more information see ponio::runge_kutta::rk_33_van_der_houwen_t)

Template Parameters:

value_t – type of coefficient (doubleby default)

Param exp_:

exponential function for the Lawson method

Param tol:

tolerance (only for adaptive time step methods)

template<typename value_t, typename Exp_t>
auto ponio::runge_kutta::lrk_44_t = []( Exp_t exp_ , double tol=ponio::default_config::tol ){returnlawson_runge_kutta::make_lawson<butcher_rk_44<value_t>,Exp_t>(exp_,tol);}

Lawson RK (4,4) method.

see more on ponio

• stages: 4

• order: 4

• underlying Runge-Kutta method: RK (4,4) (for more information see ponio::runge_kutta::rk_44_t)

Template Parameters:

value_t – type of coefficient (doubleby default)

Param exp_:

exponential function for the Lawson method

Param tol:

tolerance (only for adaptive time step methods)

template<typename value_t, typename Exp_t>
auto ponio::runge_kutta::lrk_44_235j_t = []( Exp_t exp_ , double tol=ponio::default_config::tol ){returnlawson_runge_kutta::make_lawson<butcher_rk_44_235j<value_t>,Exp_t>(exp_,tol);}

Lawson RK (4,4) 235j method.

see more on ponio

• stages: 4

• order: 4

• underlying Runge-Kutta method: RK (4,4) 235j (for more information see ponio::runge_kutta::rk_44_235j_t)

Template Parameters:

value_t – type of coefficient (doubleby default)

Param exp_:

exponential function for the Lawson method

Param tol:

tolerance (only for adaptive time step methods)

template<typename value_t, typename Exp_t>
auto ponio::runge_kutta::lrk_44_38_t = []( Exp_t exp_ , double tol=ponio::default_config::tol ){returnlawson_runge_kutta::make_lawson<butcher_rk_44_38<value_t>,Exp_t>(exp_,tol);}

Lawson RK (4,4) 3/8 method.

see more on ponio

• stages: 4

• order: 4

• underlying Runge-Kutta method: RK (4,4) 3/8 (for more information see ponio::runge_kutta::rk_44_38_t)

Template Parameters:

value_t – type of coefficient (doubleby default)

Param exp_:

exponential function for the Lawson method

Param tol:

tolerance (only for adaptive time step methods)

template<typename value_t, typename Exp_t>
auto ponio::runge_kutta::lrk_44_ralston_t = []( Exp_t exp_ , double tol=ponio::default_config::tol ){returnlawson_runge_kutta::make_lawson<butcher_rk_44_ralston<value_t>,Exp_t>(exp_,tol);}

Lawson RK (4,4) Ralston method.

see more on ponio

• stages: 4

• order: 4

• underlying Runge-Kutta method: RK (4,4) Ralston (for more information see ponio::runge_kutta::rk_44_ralston_t)

Template Parameters:

value_t – type of coefficient (doubleby default)

Param exp_:

exponential function for the Lawson method

Param tol:

tolerance (only for adaptive time step methods)

template<typename value_t, typename Exp_t>
auto ponio::runge_kutta::lrk_65_t = []( Exp_t exp_ , double tol=ponio::default_config::tol ){returnlawson_runge_kutta::make_lawson<butcher_rk_65<value_t>,Exp_t>(exp_,tol);}

Lawson RK (6,5) method.

see more on ponio

• stages: 6

• order: 5

• underlying Runge-Kutta method: RK (6,5) (for more information see ponio::runge_kutta::rk_65_t)

Template Parameters:

value_t – type of coefficient (doubleby default)

Param exp_:

exponential function for the Lawson method

Param tol:

tolerance (only for adaptive time step methods)

template<typename value_t, typename Exp_t>
auto ponio::runge_kutta::lrk_65_236a_t = []( Exp_t exp_ , double tol=ponio::default_config::tol ){returnlawson_runge_kutta::make_lawson<butcher_rk_65_236a<value_t>,Exp_t>(exp_,tol);}

Lawson RK (6,5) 236a method.

see more on ponio

• stages: 6

• order: 5

• underlying Runge-Kutta method: RK (6,5) 236a (for more information see ponio::runge_kutta::rk_65_236a_t)

Template Parameters:

value_t – type of coefficient (doubleby default)

Param exp_:

exponential function for the Lawson method

Param tol:

tolerance (only for adaptive time step methods)

template<typename value_t, typename Exp_t>
auto ponio::runge_kutta::lrk_76_t = []( Exp_t exp_ , double tol=ponio::default_config::tol ){returnlawson_runge_kutta::make_lawson<butcher_rk_76<value_t>,Exp_t>(exp_,tol);}

Lawson RK (7,6) method.

see more on ponio

• stages: 7

• order: 6

• underlying Runge-Kutta method: RK (7,6) (for more information see ponio::runge_kutta::rk_76_t)

Template Parameters:

value_t – type of coefficient (doubleby default)

Param exp_:

exponential function for the Lawson method

Param tol:

tolerance (only for adaptive time step methods)

template<typename value_t, typename Exp_t>
auto ponio::runge_kutta::lrk_86_t = []( Exp_t exp_ , double tol=ponio::default_config::tol ){returnlawson_runge_kutta::make_lawson<butcher_rk_86<value_t>,Exp_t>(exp_,tol);}

Lawson RK (8,6) method.

see more on ponio

• stages: 8

• order: 6

• underlying Runge-Kutta method: RK (8,6) (for more information see ponio::runge_kutta::rk_86_t)

Template Parameters:

value_t – type of coefficient (doubleby default)

Param exp_:

exponential function for the Lawson method

Param tol:

tolerance (only for adaptive time step methods)

template<typename value_t, typename Exp_t>
auto ponio::runge_kutta::lrk_nssp_21_t = []( Exp_t exp_ , double tol=ponio::default_config::tol ){returnlawson_runge_kutta::make_lawson<butcher_rk_nssp_21<value_t>,Exp_t>(exp_,tol);}

Lawson RK NSSP (2,1) method.

see more on ponio

• stages: 2

• order: 1

• underlying Runge-Kutta method: RK NSSP (2,1) (for more information see ponio::runge_kutta::rk_nssp_21_t)

Template Parameters:

value_t – type of coefficient (doubleby default)

Param exp_:

exponential function for the Lawson method

Param tol:

tolerance (only for adaptive time step methods)

template<typename value_t, typename Exp_t>
auto ponio::runge_kutta::lrk_nssp_32_t = []( Exp_t exp_ , double tol=ponio::default_config::tol ){returnlawson_runge_kutta::make_lawson<butcher_rk_nssp_32<value_t>,Exp_t>(exp_,tol);}

Lawson RK NSSP (3,2) method.

see more on ponio

• stages: 3

• order: 3

• underlying Runge-Kutta method: RK NSSP (3,2) (for more information see ponio::runge_kutta::rk_nssp_32_t)

Template Parameters:

value_t – type of coefficient (doubleby default)

Param exp_:

exponential function for the Lawson method

Param tol:

tolerance (only for adaptive time step methods)

template<typename value_t, typename Exp_t>
auto ponio::runge_kutta::lrk_nssp_33_t = []( Exp_t exp_ , double tol=ponio::default_config::tol ){returnlawson_runge_kutta::make_lawson<butcher_rk_nssp_33<value_t>,Exp_t>(exp_,tol);}

Lawson RK NSSP (3,3) method.

see more on ponio

• stages: 3

• order: 3

• underlying Runge-Kutta method: RK NSSP (3,3) (for more information see ponio::runge_kutta::rk_nssp_33_t)

Template Parameters:

value_t – type of coefficient (doubleby default)

Param exp_:

exponential function for the Lawson method

Param tol:

tolerance (only for adaptive time step methods)

template<typename value_t, typename Exp_t>
auto ponio::runge_kutta::lrk_nssp_53_t = []( Exp_t exp_ , double tol=ponio::default_config::tol ){returnlawson_runge_kutta::make_lawson<butcher_rk_nssp_53<value_t>,Exp_t>(exp_,tol);}

Lawson RK NSSP (5,3) method.

see more on ponio

• stages: 5

• order: 3

• underlying Runge-Kutta method: RK NSSP (5,3) (for more information see ponio::runge_kutta::rk_nssp_53_t)

Template Parameters:

value_t – type of coefficient (doubleby default)

Param exp_:

exponential function for the Lawson method

Param tol:

tolerance (only for adaptive time step methods)

template<typename value_t, typename Exp_t>
auto ponio::runge_kutta::lrk_spp_43_t = []( Exp_t exp_ , double tol=ponio::default_config::tol ){returnlawson_runge_kutta::make_lawson<butcher_rk_spp_43<value_t>,Exp_t>(exp_,tol);}

Lawson RK SPP (4,3) method.

see more on ponio

• stages: 4

• order: 3

• underlying Runge-Kutta method: RK SPP (4,3) (for more information see ponio::runge_kutta::rk_spp_43_t)

Template Parameters:

value_t – type of coefficient (doubleby default)

Param exp_:

exponential function for the Lawson method

Param tol:

tolerance (only for adaptive time step methods)

template<typename value_t, typename Exp_t>
auto ponio::runge_kutta::lrk_ssp_22_heun_t = []( Exp_t exp_ , double tol=ponio::default_config::tol ){returnlawson_runge_kutta::make_lawson<butcher_rk_ssp_22_heun<value_t>,Exp_t>(exp_,tol);}

Lawson RK SSP (2,2) Heun method.

see more on ponio

• stages: 2

• order: 2

• underlying Runge-Kutta method: RK SSP (2,2) Heun (for more information see ponio::runge_kutta::rk_ssp_22_heun_t)

Template Parameters:

value_t – type of coefficient (doubleby default)

Param exp_:

exponential function for the Lawson method

Param tol:

tolerance (only for adaptive time step methods)

template<typename value_t, typename Exp_t>
auto ponio::runge_kutta::lrk_ssp_32_t = []( Exp_t exp_ , double tol=ponio::default_config::tol ){returnlawson_runge_kutta::make_lawson<butcher_rk_ssp_32<value_t>,Exp_t>(exp_,tol);}

Lawson RK SSP (3,2) method.

see more on ponio

• stages: 3

• order: 2

• underlying Runge-Kutta method: RK SSP (3,2) (for more information see ponio::runge_kutta::rk_ssp_32_t)

Template Parameters:

value_t – type of coefficient (doubleby default)

Param exp_:

exponential function for the Lawson method

Param tol:

tolerance (only for adaptive time step methods)

template<typename value_t, typename Exp_t>
auto ponio::runge_kutta::lrk_ssp_33_t = []( Exp_t exp_ , double tol=ponio::default_config::tol ){returnlawson_runge_kutta::make_lawson<butcher_rk_ssp_33<value_t>,Exp_t>(exp_,tol);}

Lawson RK SSP (3,3) method.

see more on ponio

• stages: 3

• order: 3

• underlying Runge-Kutta method: RK SSP (3,3) (for more information see ponio::runge_kutta::rk_ssp_33_t)

Template Parameters:

value_t – type of coefficient (doubleby default)

Param exp_:

exponential function for the Lawson method

Param tol:

tolerance (only for adaptive time step methods)

template<typename value_t, typename Exp_t>
auto ponio::runge_kutta::lrk_ssp_42_t = []( Exp_t exp_ , double tol=ponio::default_config::tol ){returnlawson_runge_kutta::make_lawson<butcher_rk_ssp_42<value_t>,Exp_t>(exp_,tol);}

Lawson RK SSP (4,2) method.

see more on ponio

• stages: 4

• order: 2

• underlying Runge-Kutta method: RK SSP (4,2) (for more information see ponio::runge_kutta::rk_ssp_42_t)

Template Parameters:

value_t – type of coefficient (doubleby default)

Param exp_:

exponential function for the Lawson method

Param tol:

tolerance (only for adaptive time step methods)

template<typename value_t, typename Exp_t>
auto ponio::runge_kutta::lrk_ssp_53_t = []( Exp_t exp_ , double tol=ponio::default_config::tol ){returnlawson_runge_kutta::make_lawson<butcher_rk_ssp_53<value_t>,Exp_t>(exp_,tol);}

Lawson RK SSP (5,3) method.

see more on ponio

• stages: 5

• order: 3

• underlying Runge-Kutta method: RK SSP (5,3) (for more information see ponio::runge_kutta::rk_ssp_53_t)

Template Parameters:

value_t – type of coefficient (doubleby default)

Param exp_:

exponential function for the Lawson method

Param tol:

tolerance (only for adaptive time step methods)

template<typename value_t, typename Exp_t>
auto ponio::runge_kutta::lrk_ssp_54_t = []( Exp_t exp_ , double tol=ponio::default_config::tol ){returnlawson_runge_kutta::make_lawson<butcher_rk_ssp_54<value_t>,Exp_t>(exp_,tol);}

Lawson RK SSP (5,4) method.

see more on ponio

• stages: 5

• order: 4

• underlying Runge-Kutta method: RK SSP (5,4) (for more information see ponio::runge_kutta::rk_ssp_54_t)

Template Parameters:

value_t – type of coefficient (doubleby default)

Param exp_:

exponential function for the Lawson method

Param tol:

tolerance (only for adaptive time step methods)

template<typename value_t, typename Exp_t>
auto ponio::runge_kutta::lrkc_202_t = []( Exp_t exp_ , double tol=ponio::default_config::tol ){returnlawson_runge_kutta::make_lawson<butcher_rkc_202<value_t>,Exp_t>(exp_,tol);}

Lawson RKC (20,2) method.

see more on ponio

• stages: 20

• order: 2

• underlying Runge-Kutta method: RKC (20,2) (for more information see ponio::runge_kutta::rkc_202_t)

Template Parameters:

value_t – type of coefficient (doubleby default)

Param exp_:

exponential function for the Lawson method

Param tol:

tolerance (only for adaptive time step methods)

template<typename value_t, typename Exp_t>
auto ponio::runge_kutta::lrkc_51_t = []( Exp_t exp_ , double tol=ponio::default_config::tol ){returnlawson_runge_kutta::make_lawson<butcher_rkc_51<value_t>,Exp_t>(exp_,tol);}

Lawson RKC (5,1) method.

see more on ponio

• stages: 5

• order: 1

• underlying Runge-Kutta method: RKC (5,1) (for more information see ponio::runge_kutta::rkc_51_t)

Template Parameters:

value_t – type of coefficient (doubleby default)

Param exp_:

exponential function for the Lawson method

Param tol:

tolerance (only for adaptive time step methods)

template<typename value_t, typename Exp_t>
auto ponio::runge_kutta::lrkc_52_t = []( Exp_t exp_ , double tol=ponio::default_config::tol ){returnlawson_runge_kutta::make_lawson<butcher_rkc_52<value_t>,Exp_t>(exp_,tol);}

Lawson RKC (5,2) method.

see more on ponio

• stages: 5

• order: 2

• underlying Runge-Kutta method: RKC (5,2) (for more information see ponio::runge_kutta::rkc_52_t)

Template Parameters:

value_t – type of coefficient (doubleby default)

Param exp_:

exponential function for the Lawson method

Param tol:

tolerance (only for adaptive time step methods)

Embedded methods

template<typename value_t, typename Exp_t>
auto ponio::runge_kutta::lrk54_6m_t = []( Exp_t exp_ , double tol=ponio::default_config::tol ){returnlawson_runge_kutta::make_lawson<butcher_rk54_6m<value_t>,Exp_t>(exp_,tol);}

Lawson RK5(4) 6M method.

see more on ponio

• stages: 6

• order: 5

• underlying Runge-Kutta method: RK5(4) 6M (for more information see ponio::runge_kutta::rk54_6m_t)

Template Parameters:

value_t – type of coefficient (doubleby default)

Param exp_:

exponential function for the Lawson method

Param tol:

tolerance (only for adaptive time step methods)

template<typename value_t, typename Exp_t>
auto ponio::runge_kutta::lrk54_7m_t = []( Exp_t exp_ , double tol=ponio::default_config::tol ){returnlawson_runge_kutta::make_lawson<butcher_rk54_7m<value_t>,Exp_t>(exp_,tol);}

Lawson RK5(4) 7M method.

see more on ponio

• stages: 7

• order: 5

• underlying Runge-Kutta method: RK5(4) 7M (for more information see ponio::runge_kutta::rk54_7m_t)

Template Parameters:

value_t – type of coefficient (doubleby default)

Param exp_:

exponential function for the Lawson method

Param tol:

tolerance (only for adaptive time step methods)

template<typename value_t, typename Exp_t>
auto ponio::runge_kutta::lrk54_7s_t = []( Exp_t exp_ , double tol=ponio::default_config::tol ){returnlawson_runge_kutta::make_lawson<butcher_rk54_7s<value_t>,Exp_t>(exp_,tol);}

Lawson RK5(4) 7S method.

see more on ponio

• stages: 7

• order: 5

• underlying Runge-Kutta method: RK5(4) 7S (for more information see ponio::runge_kutta::rk54_7s_t)

Template Parameters:

value_t – type of coefficient (doubleby default)

Param exp_:

exponential function for the Lawson method

Param tol:

tolerance (only for adaptive time step methods)

template<typename value_t, typename Exp_t>
auto ponio::runge_kutta::lrk65_8m_t = []( Exp_t exp_ , double tol=ponio::default_config::tol ){returnlawson_runge_kutta::make_lawson<butcher_rk65_8m<value_t>,Exp_t>(exp_,tol);}

Lawson RK6(5) 8M method.

see more on ponio

• stages: 8

• order: 5

• underlying Runge-Kutta method: RK6(5) 8M (for more information see ponio::runge_kutta::rk65_8m_t)

Template Parameters:

value_t – type of coefficient (doubleby default)

Param exp_:

exponential function for the Lawson method

Param tol:

tolerance (only for adaptive time step methods)

template<typename value_t, typename Exp_t>
auto ponio::runge_kutta::lrk87_13m_t = []( Exp_t exp_ , double tol=ponio::default_config::tol ){returnlawson_runge_kutta::make_lawson<butcher_rk87_13m<value_t>,Exp_t>(exp_,tol);}

Lawson RK8(7) 13M method.

see more on ponio

• stages: 13

• order: 5

• underlying Runge-Kutta method: RK8(7) 13M (for more information see ponio::runge_kutta::rk87_13m_t)

Template Parameters:

value_t – type of coefficient (doubleby default)

Param exp_:

exponential function for the Lawson method

Param tol:

tolerance (only for adaptive time step methods)

Exponential Runge-Kutta methods

Like Lawson methods, exponential Runge-Kutta methods are build to solve a problem with a linear and non-linear part, to solve exactly the problem when the nonlinear part goes to zero. This class of problem can be write as

\[\dot{u} = L u + N(t, u)\]

We solve this on a time step between \(0\) and \(\Delta t\)

\[u(t^n + \Delta t) = e^{\Delta t L} U + \int_0^{\Delta t} e^{(\Delta t - s)L}N(t^n + s, u(t^n+s)) \,\mathrm{d}s\]

Interpolation of the integral yields to build a custom Runge-Kutta method which reads as

\[\begin{split}\begin{aligned} u^{(i)} &= u^n + \Delta t\sum_j a_{ij}(\Delta t L)\cdot ( k_j + Lu^n ), \quad i = 1, \dots, s \\ k_i &= N(t^n + c_i\Delta t, u^{(i)}) \\ u^{n+1} &= u^n + \Delta t \sum_i b_i(\Delta t L)\cdot( k_i + Lu^n) \end{aligned}\end{split}\]

Note

Matrix \(A\) and vector \(b\) could contain some functions defined by

\[\varphi_\ell(z) = \frac{e^z - \sum_{k=0}^{\ell-1} \frac{1}{k!}z^k }{z^\ell}\]

and we use the notations :math:varphi_ell = varphi_ell(Delta t L) and \(\varphi_{\ell,j} = \varphi_\ell(c_j \Delta t L)\).

template<typename value_t, typename linear_t>
using ponio::runge_kutta::etd2cf3_t = exponential_runge_kutta::explicit_exp_rk_butcher<butcher_etd2cf3<value_t, linear_t>>

ETD2CF3 method.

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|ccc} 0 & 0 & 0 & 0 \\ \frac{1}{3} & \frac{\varphi_{1 2}}{3} & 0 & 0 \\ \frac{2}{3} & \frac{2 \varphi_{1 3}}{3} - \frac{4 \varphi_{2 3}}{3} & \frac{4 \varphi_{2 3}}{3} & 0 \\ \hline & \varphi_{1} - \frac{9 \varphi_{2}}{2} + 9 \varphi_{3} & 6 \varphi_{2} - 18 \varphi_{3} & - \frac{3 \varphi_{2}}{2} + 9 \varphi_{3} \end{array} \end{split}\]

Template Parameters:

• value_t – type of coefficient (doubleby default)

• linear_t – type of coefficient (doubleby default)

template<typename value_t, typename linear_t>
using ponio::runge_kutta::etd3rk_t = exponential_runge_kutta::explicit_exp_rk_butcher<butcher_etd3rk<value_t, linear_t>>

ETD3RK method.

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|ccc} 0 & 0 & 0 & 0 \\ \frac{1}{2} & \frac{\varphi_{1 2}}{2} & 0 & 0 \\ 1 & - \varphi_{1 3} & 2 \varphi_{1 3} & 0 \\ \hline & \varphi_{1} - 3 \varphi_{2} + 4 \varphi_{3} & 4 \varphi_{2} - 8 \varphi_{3} & - \varphi_{2} + 4 \varphi_{3} \end{array} \end{split}\]

Template Parameters:

• value_t – type of coefficient (doubleby default)

• linear_t – type of coefficient (doubleby default)

template<typename value_t, typename linear_t>
using ponio::runge_kutta::cox_matthews_t = exponential_runge_kutta::explicit_exp_rk_butcher<butcher_cox_matthews<value_t, linear_t>>

Cox Matthews method.

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|cccc} 0 & 0 & 0 & 0 & 0 \\ \frac{1}{2} & \frac{\varphi_{1 2}}{2} & 0 & 0 & 0 \\ \frac{1}{2} & 0 & \frac{\varphi_{1 3}}{2} & 0 & 0 \\ 1 & \frac{\varphi_{1 3} \left(\varphi_{0 3} - 1\right)}{2} & 0 & \varphi_{1 3} & 0 \\ \hline & \varphi_{1} - 3 \varphi_{2} + 4 \varphi_{3} & 2 \varphi_{2} - 4 \varphi_{3} & 2 \varphi_{2} - 4 \varphi_{3} & - \varphi_{2} + 4 \varphi_{3} \end{array} \end{split}\]

Template Parameters:

• value_t – type of coefficient (doubleby default)

• linear_t – type of coefficient (doubleby default)

template<typename value_t, typename linear_t>
using ponio::runge_kutta::exprk22_t = exponential_runge_kutta::explicit_exp_rk_butcher<butcher_exprk22<value_t, linear_t>>

expRK(2,2) method

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|cc} 0 & 0 & 0 \\ 1 & \varphi_{1 2} & 0 \\ \hline & \varphi_{1} - \varphi_{2} & \varphi_{2} \end{array} \end{split}\]

Template Parameters:

• value_t – type of coefficient (doubleby default)

• linear_t – type of coefficient (doubleby default)

template<typename value_t, typename linear_t>
using ponio::runge_kutta::hochbruck_ostermann_t = exponential_runge_kutta::explicit_exp_rk_butcher<butcher_hochbruck_ostermann<value_t, linear_t>>

Hochbruck Ostermann method.

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|ccccc} 0 & 0 & 0 & 0 & 0 & 0 \\ \frac{1}{2} & \frac{\varphi_{1 2}}{2} & 0 & 0 & 0 & 0 \\ \frac{1}{2} & \frac{\varphi_{1 3}}{2} - \varphi_{2 3} & \varphi_{2 3} & 0 & 0 & 0 \\ 1 & \varphi_{1 4} - 2 \varphi_{2 4} & \varphi_{2 4} & \varphi_{2 4} & 0 & 0 \\ \frac{1}{2} & \frac{\varphi_{1 5}}{2} - \frac{\varphi_{2 4}}{4} - \frac{3 \varphi_{2 5}}{4} + \varphi_{3 4} + \frac{\varphi_{3 5}}{2} & \frac{\varphi_{2 4}}{4} + \frac{\varphi_{2 5}}{2} - \varphi_{3 4} - \frac{\varphi_{3 5}}{2} & \frac{\varphi_{2 4}}{4} + \frac{\varphi_{2 5}}{2} - \varphi_{3 4} - \frac{\varphi_{3 5}}{2} & - \frac{\varphi_{2 4}}{4} - \frac{\varphi_{2 5}}{4} + \varphi_{3 4} + \frac{\varphi_{3 5}}{2} & 0 \\ \hline & \varphi_{1} - 3 \varphi_{2} + 4 \varphi_{3} & 0 & 0 & - \varphi_{2} + 4 \varphi_{3} & 4 \varphi_{2} - 8 \varphi_{3} \end{array} \end{split}\]

Template Parameters:

• value_t – type of coefficient (doubleby default)

• linear_t – type of coefficient (doubleby default)

template<typename value_t, typename linear_t>
using ponio::runge_kutta::krogstad_t = exponential_runge_kutta::explicit_exp_rk_butcher<butcher_krogstad<value_t, linear_t>>

Krogstad method.

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|cccc} 0 & 0 & 0 & 0 & 0 \\ \frac{1}{2} & \frac{\varphi_{1 2}}{2} & 0 & 0 & 0 \\ \frac{1}{2} & \frac{\varphi_{1 3}}{2} - \varphi_{2 3} & \varphi_{2 3} & 0 & 0 \\ 1 & \varphi_{1 4} - 2 \varphi_{2 4} & 0 & 2 \varphi_{2 4} & 0 \\ \hline & \varphi_{1} - 3 \varphi_{2} + 4 \varphi_{3} & 2 \varphi_{2} - 4 \varphi_{3} & 2 \varphi_{2} - 4 \varphi_{3} & - \varphi_{2} + 4 \varphi_{3} \end{array} \end{split}\]

Template Parameters:

• value_t – type of coefficient (doubleby default)

• linear_t – type of coefficient (doubleby default)

template<typename value_t, typename linear_t>
using ponio::runge_kutta::strehmel_weiner_t = exponential_runge_kutta::explicit_exp_rk_butcher<butcher_strehmel_weiner<value_t, linear_t>>

Strehmel Weiner method.

This method is based on the following Butcher tableau:

\[\begin{split} \begin{array}{c|cccc} 0 & 0 & 0 & 0 & 0 \\ \frac{1}{2} & \frac{\varphi_{1 2}}{2} & 0 & 0 & 0 \\ \frac{1}{2} & \frac{\varphi_{1 3}}{2} - \frac{\varphi_{2 3}}{2} & \frac{\varphi_{2 3}}{2} & 0 & 0 \\ 1 & \varphi_{1 4} - 2 \varphi_{2 4} & - 2 \varphi_{2 4} & 4 \varphi_{2 4} & 0 \\ \hline & \varphi_{1} - 3 \varphi_{2} + 4 \varphi_{3} & 0 & 4 \varphi_{2} - 8 \varphi_{3} & - \varphi_{2} + 4 \varphi_{3} \end{array} \end{split}\]

Template Parameters:

• value_t – type of coefficient (doubleby default)

• linear_t – type of coefficient (doubleby default)

Explicit stabilized Runge-Kutta methods

Some problems, like heat equation, require methods stabilized on the negative real axis. The ponio library provides a Runge-Kutta Chebyshev method of order 2, ROCK2 method (of order 2), ROCK4 method (of order 4) and a Runge-Kutta Legendre method of order 1 or 2.

Runge-Kutta Chebyshev method

The algorithm of RKC2 is the following:

\[\begin{split}\begin{aligned} y_0 &= f(t^n, y^n) \\ y_1 &= y_0 + \tilde{\mu}_1\Delta t f(t^n, y^n) \\ y_j &= (1 - \mu_j - \nu_j)y^n + \mu_j y_{j-1} + \nu_j y_{j-2} \\ &~~~~~~ + \tilde{\mu}_j \Delta t f(t^n + c_j\Delta t, y_{j-1}) + \tilde{\gamma}_j\Delta t f(t^n, y^n), \quad j=2,\dots, s \end{aligned}\end{split}\]

with coefficients given by

\[\tilde{\mu}_1 = b_1 \omega_1\]

\[\mu_j = \frac{2b_j}{b_{j-1}}\omega_0, \quad \nu_j = -\frac{b_j}{b_{j-2}}, \quad \tilde{\mu}_j = \frac{2b_j}{b_{j-1}}\omega_1\]

\[\tilde{\gamma}_j = -(1 - b_{j-1}T_{j-1}(\omega_0))\]

where

\[b_0 = b_2, \quad b_1 = \frac{1}{\omega_0}, \quad b_j = \frac{T_j''(\omega_0)}{(T_j'(\omega_0))^2},\ j=2,\dots, s\]

\[c_0 = 0, \quad c_1 = c_2, \quad c_j = \frac{T_s'(\omega_0)}{T_s''(\omega_0)}\frac{T_j''(\omega_0)}{T_j'(\omega_0)}, \quad c_s = 1\]

and

\[\omega_0 = 1 + \frac{\epsilon}{s^2},\quad \omega_1 = \frac{T_s'(\omega_0)}{T_s''(\omega_0)},\quad \epsilon \approx \frac{2}{13}\]

and where \(T_j(x)\) is the Chebyshev polynomial.

template<std::size_t N_stages_, typename _value_t = double>
class explicit_rkc2

define RKC2 with user defined number of stages

Template Parameters:

• N_stages_ – number of stages

• value_t – type of coefficients

Public Functions

inline explicit_rkc2(value_t eps = 2. / 13.)

Construct a new explicit RKC2 algorithm.

Parameters:

eps – value of relaxation

template<typename problem_t, typename state_t, typename array_ki_t, std::size_t j>
inline void stage(Stage<j>, problem_t &f, value_t tn, state_t &yn, array_ki_t const &Yj, value_t dt, state_t &yi)

computes stage \(j\) of RKC2 method

\(y_j = (1 - \mu_j - \nu_j)y^n + \mu_j y_{j-1} + \nu_j y_{j-2} + \tilde{\mu}_j\Delta tf(t^n + c_{j-1}\Delta t, y_{j-1}) + \tilde{\gamma}_j\Delta t f(t^n, y^n)\)

Template Parameters:

• problem_t – type of operator \(f\)

• state_t – type of current state

• array_ki_t – type of array of temporary stages

• j – integer of stage

Parameters:

• f – operator \(f\)

• tn – current time

• yn – current state

• Yj – array of temporary stages

• dt – current time step

• yi – computed output step

template<typename problem_t, typename state_t, typename array_ki_t>
inline void stage(Stage<0>, problem_t &f, value_t tn, state_t &yn, array_ki_t const&, value_t, state_t &yi)

computes pseudo first stage of RKC2 method

\(y_0 = f(t^n, y^n)\)

Template Parameters:

• problem_t – type of operator \(f\)

• state_t – type of current state

• array_ki_t – type of array of temporary stages

• j – integer of stage

Parameters:

• f – operator \(f\)

• tn – current time

• yn – current state

• yi – computed output step

template<typename problem_t, typename state_t, typename array_ki_t>
inline void stage(Stage<1>, problem_t&, value_t, state_t &yn, array_ki_t const &Yj, value_t dt, state_t &yi)

computes first stage of RKC2 method

\(y_1 = y^n + \tilde{\mu}_1 \Delta t f(t^n, y^n)\)

Template Parameters:

• problem_t – type of operator \(f\)

• state_t – type of current state

• array_ki_t – type of array of temporary stages

• j – integer of stage

Parameters:

• yn – current state

• Yj – array of temporary stages

• dt – current time step

• yi – computed output step

template<typename problem_t, typename state_t, typename array_ki_t>
inline void stage(Stage<2>, problem_t &f, value_t tn, state_t &yn, array_ki_t const &Yj, value_t dt, state_t &yi)

computes second stage of RKC2 method

\(y_2 = (1-\mu_2 -\nu_2)y^n + \mu_2y_1 + \nu_2y^n + \tilde{\mu}_2\Delta t f(t^n + c_1\Delta t, y_1) + \tilde{\gamma}_2 \Delta t f(t^n, y^n)\)

Template Parameters:

• problem_t – type of operator \(f\)

• state_t – type of current state

• array_ki_t – type of array of temporary stages

• j – integer of stage

Parameters:

• f – operator \(f\)

• tn – current time

• yn – current state

• Yj – array of temporary stages

• dt – current time step

• yi – computed output step

inline auto &info()

gets iteration_info object

inline auto const &info() const

gets iteration_info object (constant version)

Public Static Functions

template<std::size_t J>
static inline value_t b(value_t const x)

computes \(b_j = \) coefficients with following formula: \(b_0=b_2\), \(b_1=\frac{1}{\omega_0}\), \(b_j = \frac{T_j''(\omega_0)}{(T_j'(\omega_0))^2}\)

Template Parameters:

J – index \(j\)

Parameters:

x – value of \(\omega_0\)

ROCK2 method

We write the method ROCK2 presented in [AM01]. The algorithm of ROCK2 method is the following:

\[\begin{split}\begin{aligned} y_0 &= y^n \\ y_1 &= y^n + \Delta t \mu_1 f(y^n) \\ y_j &= \Delta t \mu_j f(y_{j-1}) - \nu_j y_{j-1} - \kappa_j y_{j-2}, \quad j=2,\dots, s-2\\ y_{s-1} &= y_{s-2} + \Delta t \sigma f(u_{s-2}) \\ y_{s}^\star &= y_{s-1} + \Delta t \sigma f(y_{s-1}) \\ y^{n+1} &= y_s^\star - \Delta t \sigma (1 - \frac{\tau}{\sigma})( f(y_{s-1}) - f(y_{s-2}) ) \end{aligned}\end{split}\]

where \(\mu_j\), \(\nu_j\) and \(\kappa_j\) coefficients coming from a minimization problem.

template<typename eig_computer_t, bool _is_embedded = false, typename _value_t = double>
class rock2_impl

define ROCK2 method

Template Parameters:

• eig_computer_t – type of computer of maximal eigenvalue

• _is_embedded – define if method is used as adaptive or constant time step method [default is false]

• _value_t – type of coefficients

Public Functions

inline rock2_impl(eig_computer_t &&_eig_computer, value_t _a_tol = 1e-4, value_t _r_tol = 1e-4)

Construct a new ROCK2 algorithm object.

Parameters:

• _eig_computer – estimator of spectral radius

• _a_tol – absolute tolerance (for adaptive time step method)

• _r_tol – relative tolerance (for adaptive time step method)

template<typename state_t>
inline auto error(state_t &&unp1, state_t &&un, state_t &&tmp)

computes an error for adaptive time step method

Template Parameters:

state_t – type of current state

Parameters:

• unp1 – estimation of solution at time \(t^{n+1} = t^n+\Delta t\)

• un – solution at time \(t^n\)

• tmp – \(-\Delta t \sigma(1 - \frac{\tau}{\sigma^2})(f(u_{s-1}) - f(u_{s-2}))\)

Returns:

auto estimation of error to compare to 1

template<typename problem_t, typename state_t, typename array_ki_t>
inline void operator()(problem_t &f, value_t &tn, state_t &un, array_ki_t &G, value_t &dt, state_t &unp1)

iteration of ROCK2 method

Template Parameters:

• problem_t – type of \(f\)

• state_t – type of current state

• array_ki_t – type of temporary stages (only 3 needed for ROCK2)

Parameters:

• f – operator \(f\)

• tn – current time

• un – current state

• G – array of temporary stages

• dt – current time step

template<bool is_embedded = false, typename value_t = double, typename eig_computer_t>
auto ponio::runge_kutta::rock::rock2(eig_computer_t &&eig_computer)

helper to build a rock2_impl object

Template Parameters:

• is_embedded – define if method is used as adaptive or constant time step method [default is false]

• value_t – type of coefficients

• eig_computer_t – type of computer of maximal eigenvalue

Parameters:

eig_computer – computer of maximal eigenvalues

template<bool is_embedded = false, typename value_t = double>
auto ponio::runge_kutta::rock::rock2()

helper to build a rock2_impl object with power method to compute maximal eigen value

Template Parameters:

• is_embedded – define if method is used as adaptive or constant time step method [default is false]

• value_t – type of coefficients

ROCK4 method

We write the method ROCK2 presented in [Abd02]. The algorithm of ROCK4 method is the following:

\[\begin{split}\begin{aligned} y_0 &= y^n \\ y_1 &= y^n + \Delta t \mu_1 f(y^n) \\ y_j &= \Delta t \mu_j f(y_{j-1}) - \nu_j y_{j-1} - \kappa_j y_{j-2}, \quad j=2,\dots, s-4 \\ y_{s-3} &= y_{s-4} + a_{21} \Delta t f(u_{s-4}) \\ y_{s-2} &= y_{s-4} + a_{31} \Delta t f(u_{s-4}) + a_{32} \Delta t f(y_{s-3}) \\ y_{s-1} &= y_{s-4} + a_{41} \Delta t f(u_{s-4}) + a_{42} \Delta t f(y_{s-3}) + a_{43} \Delta t f(y_{s-2}) \\ y^{n+1} &= y_{s-4} + b_1 \Delta t f(y_{s-4}) + b_2 \Delta t f(y_{s-3}) + b_3 \Delta t f(y_{s-2}) + b_4 \Delta t f(y_{s-1}) \end{aligned}\end{split}\]

where \(\mu_j\), \(\nu_j\) and \(\kappa_j\) coefficients coming from a minimization problem, and \(a_{ij}\), \(b_i\) coming from an order 4 method build with \(y_{s-4}\) as initial condition.

template<typename eig_computer_t, bool _is_embedded = false, typename _value_t = double>
class rock4_impl

define ROCK4 method

Template Parameters:

• eig_computer_t – type of computer of maximal eigenvalue

• _is_embedded – define if method is used as adaptive or constant time step method [default is false]

• _value_t – type of coefficients

Public Functions

inline rock4_impl(eig_computer_t &&_eig_computer, value_t _a_tol = 1e-4, value_t _r_tol = 1e-4)

Construct a new ROCK4 algorithm object.

Parameters:

• _eig_computer – estimator of spectral radius

• _a_tol – absolute tolerance (for adaptive time step method)

• _r_tol – relative tolerance (for adaptive time step method)

template<typename state_t>
inline auto error(state_t const &unp1, state_t const &tmp)

computes an error for adaptive time step method

Template Parameters:

state_t – type of current state

Parameters:

• unp1 – estimation of solution at time \(t^{n+1} = t^n+\Delta t\)

• tmp – other estimation of solution at time \(t^{n+1} = t^n+\Delta t\)

Returns:

auto estimation of error to compare to 1

template<typename problem_t, typename state_t, typename array_ki_t>
inline void operator()(problem_t &f, value_t &tn, state_t &un, array_ki_t &G, value_t &dt, state_t &unp1)

iteration of ROCK4 method

Template Parameters:

• problem_t – type of \(f\)

• state_t – type of current state

• array_ki_t – type of temporary stages (only 6 needed for ROCK4)

Parameters:

• f – operator \(f\)

• tn – current time

• un – current state

• G – array of temporary stages

• dt – current time step

• u_np1 – solution \(u^{n+1}\) at time \(t^{n+1} = t^n + \Delta t\)

inline auto &info()

gets iteration_info object

inline auto const &info() const

gets iteration_info object (constant version)

template<bool is_embedded = false, typename value_t = double, typename eig_computer_t>
auto ponio::runge_kutta::rock::rock4(eig_computer_t &&eig_computer)

helper to build a rock4_impl object

Template Parameters:

• is_embedded – define if method is used as adaptive or constant time step method [default is false]

• value_t – type of coefficients

• eig_computer_t – type of computer of maximal eigenvalue

Parameters:

eig_computer – computer of maximal eigenvalues

template<bool is_embedded = false, typename value_t = double>
auto ponio::runge_kutta::rock::rock4()

helper to build a rock4_impl object with power method to compute maximal eigen value

Template Parameters:

• is_embedded – define if method is used as adaptive or constant time step method [default is false]

• value_t – type of coefficients

Runge-Kutta Legendre method

An other way to get a stabilized Runge-Kutta method is to use Legendre polynomials, we follow presentation in [MBA14].

The algorithm of RKL1 is the following:

\[\begin{split}\begin{aligned} y^{(0)} &= y^n \\ y^{(1)} &= y^n + \tilde{\mu}_1\Delta t f(t^n, y^n) \\ y^{(j)} &= \mu_jy^{(j-1)} + \nu_j y^{(j-2)} + \tilde{\mu}_j\Delta t f(t^n, y^{(j-1)}), \quad j=2,\dots, s \\ y^{(n+1)} &= y^{(s)} \end{aligned}\end{split}\]

with coefficients given by

\[\mu_j = \frac{2j-1}{j}, \qquad \nu_j = \frac{1-j}{j}\]

\[\tilde{\mu}_j = \frac{2j-1}{j}\frac{2}{s^2 + s}\]

where \(s\) is the number of stages of the method.

template<std::size_t N_stages_, typename _value_t = double>
class explicit_rkl1

define RKL1 (Runge-Kutta-Legendre method of ordre 1) with user defined number of stages

Warning

the method is only presented with an autonomous problem, ie \(\dot{y} = f(y)\)

Template Parameters:

• N_stages_ – number of stages

• _value_t – type of coefficients

Public Functions

template<typename problem_t, typename state_t, typename array_ki_t, std::size_t j>
inline void stage(Stage<j>, problem_t &f, value_t tn, state_t&, array_ki_t const &Yj, value_t dt, state_t &yi)

compute stage \(j\) of RKL1 method

\(y^{(j)} = \mu_j y^{(j-1)} + \nu_j y^{(j-2)} + \tilde{\mu}_j \Delta t f(t^n, y^{(j-1)})\)

Template Parameters:

• problem_t – type of operator \(f\)

• state_t – type of current state

• array_ki_t – type of array of temporary stages

• j – integer of stage

Parameters:

• f – operator \(f\)

• tn – current time

• yn – current state

• Yj – array of temporary stages

• dt – current time step

• yi – computed output step

template<typename problem_t, typename state_t, typename array_ki_t>
inline void stage(Stage<0>, problem_t&, value_t, state_t &yn, array_ki_t const&, value_t, state_t &yi)

compute stage \(0\) of RKL1 method

\(y^{(0)} = y^n\)

Template Parameters:

• problem_t – type of operator \(f\)

• state_t – type of current state

• array_ki_t – type of array of temporary stages

Parameters:

• yn – current state

• yi – computed output step

template<typename problem_t, typename state_t, typename array_ki_t>
inline void stage(Stage<1>, problem_t &f, value_t tn, state_t &yn, array_ki_t const&, value_t dt, state_t &yi)

compute stage \(1\) of RKL1 method

\(y^{(1)} = y^n + \tilde{\mu}_1 \Delta t f(t^n, y^n)\)

Template Parameters:

• problem_t – type of operator \(f\)

• state_t – type of current state

• array_ki_t – type of array of temporary stages

Parameters:

• f – operator \(f\)

• tn – current time

• yn – current state

• dt – current time step

• yi – computed output step

inline auto &info()

gets iteration_info object

inline auto const &info() const

gets iteration_info object (constant version)

Public Static Functions

template<std::size_t j>
static inline value_t mu()

compute \(\mu_j = \frac{2j-1}{j}\)

Template Parameters:

j – index \(j\)

template<std::size_t j>
static inline value_t nu()

compute \(\nu_j = \frac{1-j}{j}\)

Template Parameters:

j – index \(j\)

template<std::size_t j>
static inline value_t mu_t()

compute \(\tilde{\mu}_j = \frac{2j-1}{j}\frac{2}{s^2 + s}\) with \(s\) the number of stages

Template Parameters:

j – index \(j\)

The algorithm of RKL2 is the following:

\[\begin{split}\begin{aligned} y^{(0)} &= y^n \\ y^{(1)} &= y^{(0)} + \tilde{\mu}_1\Delta t f(t^n, y^{(0)}) \\ y^{(j)} &= \mu_jy^{(j-1)} + \nu_j y^{(j-2)} + (1-\mu_j-\nu_j)y^{(0)} + \tilde{\mu}_j\Delta t f(t^n, y^{(j-1)}) + \tilde{\gamma}_j\Delta tf(t^n, y^{(0)}), \quad j=2,\dots, s \\ y^{(n+1)} &= y^{(s)} \end{aligned}\end{split}\]

with coefficients given by

\[\mu_j = \frac{2j-1}{j}\frac{b_j}{b_{j-1}}, \qquad \nu_j = -\frac{j-1}{j}\frac{b_j}{b_{j-2}}\]

\[\tilde{\mu}_j = \mu_j w_1,\ 1<j \qquad \tilde{\mu}_1 = b_1 w_1\]

\[\tilde{\gamma}_j = -a_{j-1}\tilde{\mu}_j\]

where

\[b_0 = b_1 = b_2 = \frac{1}{3},\qquad b_j = \frac{j^2 + j - 2}{2j(j+1)},\ 2\leq j\]

and

\[a_j = 1 - b_j, \qquad w_1 = \frac{4}{s^2 + s - 2}\]

where \(s\) is the number of stages of the method

template<std::size_t N_stages_, typename _value_t = double>
class explicit_rkl2

define RKL2 (Runge-Kutta-Legendre method of ordre 2) with user defined number of stages

Warning

the method is only presented with an autonomous problem, ie \(\dot{y} = f(y)\)

Template Parameters:

• N_stages_ – number of stages

• _value_t – type of coefficients

Public Functions

template<typename problem_t, typename state_t, typename array_ki_t, std::size_t j>
inline void stage(Stage<j>, problem_t &f, value_t tn, state_t &yn, array_ki_t const &Yj, value_t dt, state_t &yi)

compute stage \(j\) of RKL2 method

\(y^{(j)} = \mu_j y^{(j-1)} + \nu_j y^{(j-2)} + (1-\mu_j-\nu_j)y^{(0)} + \tilde{\mu}_j \Delta t f(t^n, y^{(j-1)}) + \gamma_j\Delta t f(t^n, y^{(0)})\)

Template Parameters:

• problem_t – type of operator \(f\)

• state_t – type of current state

• array_ki_t – type of array of temporary stages

• j – integer of stage

Parameters:

• f – operator \(f\)

• tn – current time

• yn – current state

• Yj – array of temporary stages

• dt – current time step

• yi – computed output step

template<typename problem_t, typename state_t, typename array_ki_t>
inline void stage(Stage<0>, problem_t &f, value_t tn, state_t &yn, array_ki_t const&, value_t dt, state_t &yi)

compute \(\Delta t f(t^n, y^n)\) term (needs for all other stages)

Warning

first stage of RKL2 method is \(y^n\) but here we compute \(\Delta t f(t^n, y^n)\) term

Template Parameters:

• problem_t – type of operator \(f\)

• state_t – type of current state

• array_ki_t – type of array of temporary stages

Parameters:

• f – operator \(f\)

• tn – current time

• yn – current state

• dt – current time step

• yi – computed output step

template<typename problem_t, typename state_t, typename array_ki_t>
inline void stage(Stage<1>, problem_t&, value_t, state_t &yn, array_ki_t const &Yj, value_t, state_t &yi)

compute stage \(1\) of RKL2 method

\(y^{(1)} = y^n + \tilde{\mu}_1 \Delta t f(t^n, y^n)\)

Template Parameters:

• problem_t – type of operator \(f\)

• state_t – type of current state

• array_ki_t – type of array of temporary stages

Parameters:

• yn – current state

• Yj – array of temporary stages

• yi – computed output step

template<typename problem_t, typename state_t, typename array_ki_t>
inline void stage(Stage<2>, problem_t &f, value_t tn, state_t &yn, array_ki_t const &Yj, value_t dt, state_t &yi)

compute stage \(2\) of RKL2 method

\(y^{(2)} = \mu_2 y^{(1)} + \nu_2y^n + (1-\mu_2-\nu_2)y^n + \tilde{\mu}_2\Delta tf(t^n, y^{(1)}) + \gamma_2\Delta t f(t^n, y^n)\)

Template Parameters:

• problem_t – type of operator \(f\)

• state_t – type of current state

• array_ki_t – type of array of temporary stages

Parameters:

• f – operator \(f\)

• tn – current time

• yn – current state

• Yj – array of temporary stages

• dt – current time step

• yi – computed output step

inline auto &info()

gets iteration_info object

inline auto const &info() const

gets iteration_info object (constant version)

Public Static Functions

static inline value_t w1()

compute \(w_1\) coefficient

\(w_1 = \frac{4}{s^2 + s -2}\) with \(s\) the number of stages of the method

template<std::size_t j>
static inline value_t mu()

compute \(\mu_j\) coefficient

\(b_1 = \frac{1}{3}\)

Template Parameters:

j – index \(j\)

template<std::size_t j>
static inline value_t nu()

compute \(\mu_j\) coefficient

\(mu_j = \frac{2j-1}{j}\frac{b_j}{b_{j-1}}\)

Template Parameters:

j – index \(j\)

template<std::size_t j>
static inline value_t mu_t()

compute \(\tilde{\mu}_j\) coefficient

\(\tilde{\mu}_j = \mu_j w_1\) for \(1<j\), \(\tilde{\mu}_1 - b_1w_1\)

Template Parameters:

j – index \(j\)

template<std::size_t j>
static inline value_t gamma_t()

compute \(\gamma_j\) coefficient

\(gamma_j = -a_{j-1}\tilde{\mu}_j\)

Template Parameters:

j – index \(j\)

IMEX stabilized methods

The PIROCK method is introduce in [AV13], the complete scheme is a IMEX scheme that allows for solving an equation of the following form (with 3 operators):

\[\dot{u} = F_R(u) + F_D(u) + F_A(u)\]

with \(F_R\) a reaction operator (solved implicitly), \(F_D\) a diffusion operator (solved explicitly by a modified ROCK2 method) and \(F_A\) an advection operator (solved with an explicit RK3 method). The first implementation of PIROCK scheme in ponio makes to solve only reaction-diffusion problem (e.g. \(F_A = 0\)), the second one to solve a complet reaction-diffusion-advection problem.

The method is divided into 5 steps:

1. Diffusion stabilization procedure (modified ROCK2 method)

   \[\begin{split}\begin{aligned} u^{(0)} &= u^n \\ u^{(1)} &= u^n + \alpha \mu_1 \Delta t F_D(u^n) \\ u^{(j)} &= \alpha \mu_j \Delta t F_D(u^{(j-1)}) - \nu_j u^{(j-1)} - \kappa_j u^{(j-2)},\quad j=2,\dots, s-2+\ell \end{aligned}\end{split}\]

2. Finishing procedure for diffusion

   \[\begin{split}\begin{aligned} u^{*(s-1)} &= u^{(s-2)} + \sigma_\alpha \Delta t F_D(u^{(s-2)}) \\ u^{*(s)} &= u^{*(s-1)} + \sigma_\alpha \Delta t F_D(u^{*(s-1)}) \end{aligned}\end{split}\]

3. Starting value for advection-reaction

   \[U = u^{(s-2+\ell)}\]

4. Finishing procedure for advection-reaction coupling

   \[\begin{split}\begin{aligned} u^{(s+1)} &= u^{(s-2+\ell)} + \gamma \Delta t F_R(u^{(s+1)}) \\ u^{(s+2)} &= u^{(s-2+\ell)} + \beta \Delta t F_D(u^{(s+1)}) + \Delta t F_A(u^{(s+1)}) + (1-2\gamma)\Delta t F_R(u^{(s+1)}) + \gamma \Delta t F_R(u^{(s+2)}) \\ u^{(s+3)} &= u^{(s-2+\ell)} + (1-2\gamma)\Delta t F_A(u^{(s+1)}) + (1-\gamma)\Delta t F_R(u^{(s+1)}) \\ u^{(s+4)} &= u^{(s-2+\ell)} + \frac{1}{3}\Delta t F_A(u^{(s+1)}) \\ u^{(s+5)} &= u^{(s-2+\ell)} + \frac{2\beta}{3} \Delta t F_D(u^{(s+1)}) + \frac{2}{3}\Delta t J_R^{-1} F_A(u^{(s+4)}) + \left(\frac{2}{3} - \gamma\right) \Delta t F_R(u^{(s+1)}) + \frac{2\gamma}{3}\Delta t F_R(u^{(s+2)}) \end{aligned}\end{split}\]

5. Computation of the integration step

   \[\begin{split}\begin{aligned} u^{n+1} =& u^{*(s)} \\ & - \sigma_\alpha\left( 1 - \frac{\tau_\alpha}{\sigma_\alpha^2}\right)\Delta t (F_D(u^{*(s-1)}) - F_D(u^{(s-2)})) \\ & + \frac{1}{4}\Delta t F_A(u^{(s+1)}) + \frac{3}{4}\Delta t F_A(u^{(s+5)}) \\ & + \frac{1}{2}\Delta t F_R(u^{(s+1)}) + \frac{1}{2}\Delta t F_R(u^{(s+2)}) \\ & + \frac{J_R^{-\ell}}{2-4\gamma}\Delta t (F_D(u^{(s+3)}) - F_D(u^{(s+1)})) \end{aligned}\end{split}\]

   where \(\mu_j\), \(\nu_j\), \(\kappa_j\) are the same coefficients as for the standard ROCK2 method and:

   • \(\gamma = 1-\frac{\sqrt{2}}{2}\)

   • \(\beta = 1-2\alpha P'_{s-2+\ell}(0)\) where \(P'_{s-2+\ell}\) is the stability polynomial of the underlying ROCK2 method

   • \(J_R = I - \gamma \Delta t \frac{\partial F_R}{\partial u}(u^{(s-2+\ell)})\)

   • \(\sigma_\alpha = \frac{1-\alpha}{2} + \alpha \sigma\), where \(\sigma\) coming from the underlying ROCK2 method

   • \(\tau_\alpha = \frac{(\alpha-1)^2}{2} + 2\alpha(1-\alpha)\sigma + \alpha^2\tau\), where \(\sigma\) and \(\tau\) coming from the underlying ROCK2 method

For adaptive time step method, need to compute one error by operator and take the maximum value

\[\begin{split}\begin{aligned} err_D &= \sigma_\alpha \Delta t \left( 1 - \frac{\tau_\alpha}{\sigma_\alpha^2} \right) ( F_D( u^{*(s-1)} - u^{(s-2)} )) \\ err_R &= \Delta t J_R^{-1}\left( \frac{1}{6}F_R(u^{(s+1)}) - \frac{1}{6}F_R(u^{(s+2)}) \right) \\ err_A &= -\frac{3}{20}\Delta t F_A(u^{(s+1)}) + \frac{3}{10}\Delta t F_A(u^{(s+4)}) - \frac{3}{20}\Delta t F_A(u^{(s+5)}) \end{aligned}\end{split}\]

The error is computed with the maximum:

\[err = \max\left( \|err_D\|, \|err_R\|, \|err_A\|^{2/3} \right)\]

with the following error norm:

\[\|\cdot\| : err \mapsto \sum_i \left( \frac{err_i}{a_{tol} + \max(y^{n+1}_i, y^n_i) } \right)^2\]

where \(y^n\) and \(y^{n+1}\) are estimation of the solution respectively at time \(t^n\) and after an iteration \(t^{n+1}=t^n+\Delta t\), and \(a_{tol}\) and \(r_{tol}\) respectively absolute and relative tolerances. This norm give a normalized error, so its compare to \(1\) and new time step is computed as:

\[\Delta t_{new} = \sqrt{\frac{1}{err}} \Delta t\]

Note

As in a lot of adaptive time step method, the new time step stays in \(\Delta t_{new} \in [0.5\Delta t, 2 \Delta t]\) set and there is a safety factor, so the new time step is multiply by \(0.8\).

Warning

Adaptive time step method is not yet implemented for complet PIROCK ponio::runge_kutta::pirock::pirock_RDA_impl.

Still keep two free parameters \(\ell\) and \(\alpha\) given free to user, ponio provides two choices for this parameters:

• \(\ell=2\) and \(\alpha=1\), in this case, if \(F_A=0\) and \(F_R=0\) we have the standard ROCK2 method;

• \(\ell=1\) and \(\alpha = \frac{1}{2P'_{s-2+\ell}(0)}\), so \(\beta=0\) to minimized computation cost.

ponio provides two computers for \(\alpha\) and \(\beta\) values.

template<typename value_t = double>
class alpha_fixed

Computer of \(\alpha\) and \(\beta\) parameters with fixed value of \(\alpha\).

Template Parameters:

value_t – Type of coefficients

Public Functions

inline value_t alpha(std::size_t, std::size_t) const

Return \(\alpha\) (given value in constructor)

inline value_t beta(std::size_t s, std::size_t l) const

Return \(\beta = 1 - 2\alpha P'_{s-2+\ell}(0)\).

Parameters:

• s – number of stages

• l – choosen parameter \(\ell\)

template<typename value_t = double>
class beta_0

Computer of \(\alpha\) and \(\beta\) parameters with fixed value of \(\beta\) to 0.

Template Parameters:

value_t – Type of coefficients

Public Functions

inline value_t alpha(std::size_t s, std::size_t l) const

Return \(\alpha = \frac{1}{2P'_{s-2+\ell}(0)}\) such \(\beta = 0\).

inline value_t beta(std::size_t, std::size_t) const

Return \(\beta = 0\).

PIROCK for reaction-diffusion problem

In this section we present the implementation of PIROCK method for only reaction-diffusion problem (i.e. \(F_A = 0\)).

template<std::size_t _l, typename alpha_beta_computer_t, typename eig_computer_t, typename _shampine_trick_caller_t = void, bool _is_embedded = false, typename _value_t = double>
class pirock_impl

implementation of PIROCK method

Template Parameters:

• _l – Number of augmented stages (l = 1, 2)

• alpha_beta_computer_t – Choice of computing of \(\alpha\) and \(\beta\) parameters

• eig_computer_t – Computing of eigenvalues of explicit part (diffusion)

• _shampine_trick_caller_t – Computing of Shampine’s trick

• _is_embedded – Set adaptive time step method (default: false)

• value_t – Type of coefficients

Public Functions

pirock_impl() = default

Construct a new pirock impl object.

inline pirock_impl(alpha_beta_computer_t &&_alpha_beta_computer, eig_computer_t &&_eig_computer)

Construct a new pirock impl object.

Parameters:

• _alpha_beta_computer – computer of parameters alpha and beta object that has two member functions which take number of stages (s) l parameter

• _eig_computer – eigenvalue computer functor (that take as argument the function of explicit part, the current time, the current state and the current time step)

template<typename _shampine_trick_caller_t_>
inline pirock_impl(alpha_beta_computer_t &&_alpha_beta_computer, eig_computer_t &&_eig_computer, _shampine_trick_caller_t_ &&_shampine_trick_caller, value_t a_tol = default_config::tol, value_t r_tol = default_config::tol)

Construct a new pirock impl object with Shampine’s trick.

Parameters:

• _alpha_beta_computer – alpha and beta computer object

• _eig_computer – eigenvalue computer functor

• _shampine_trick_caller – Shampine’s trick functor

• a_tol – absolute tolerance

• r_tol – relative tolerance

template<typename problem_t, typename state_t, typename array_ki_t>
inline void operator()(problem_t &pb, value_t &tn, state_t &un, array_ki_t &U, value_t &dt, state_t &u_np1)

iteration of PIROCK method

Template Parameters:

• problem_t – type of \(f\)

• state_t – type of current state

• array_ki_t – type of temporary stages (only 3 needed for ROCK2)

Parameters:

• pb – problem \((F_D, F_R)\) and the Jacibian of reaction part \(\frac{\partial F_R}{\partial u}\)

• tn – current time

• un – current state

• U – array of temporary stages

• dt – current time step

• u_np1 – solution \(u^{n+1}\) at time \(t^{n+1} = t^n + \Delta t\)

inline auto &info()

gets iteration_info object

inline auto const &info() const

gets iteration_info object (constant version)

Helper functions

User interface functions to build a PIROCK method.

template<std::size_t l = 1, bool is_embedded = false, typename value_t = double, typename alpha_beta_computer_t, typename eig_computer_t, typename shampine_trick_caller_t>
auto ponio::runge_kutta::pirock::pirock(alpha_beta_computer_t &&alpha_beta_computer, eig_computer_t &&eig_computer, shampine_trick_caller_t &&shampine_trick_caller, value_t absolute_tolerance = default_config::tol, value_t relative_tolerance = default_config::tol)

helper function to build PIROCK algorithm

Template Parameters:

• l – Number of augmented stages (l = 1, 2)

• is_embedded – Set adaptive time step method (default: false)

• value_t – Type of coefficients

• alpha_beta_computer_t – Choice of computing of \(\alpha\) and \(\beta\) parameters

• eig_computer_t – Computing of eigenvalues of explicit part (diffusion)

• shampine_trick_caller_t – Computing of Shampine’s trick

Parameters:

• alpha_beta_computer – \(\alpha\) and \(\beta\) computer object

• eig_computer – Eigenvalue computer of explicit part of the problem

• shampine_trick_caller – Shampine’s trick computer

• absolute_tolerance – Absolute tolerance for adaptive time step method (default: ponio::default_config::tol)

• relative_tolerance – Relative tolerance for adaptive time step method (default: ponio::default_config::tol)

template<std::size_t l = 1, typename value_t = double, typename alpha_beta_computer_t, typename eig_computer_t>
auto ponio::runge_kutta::pirock::pirock(alpha_beta_computer_t &&alpha_beta_computer, eig_computer_t &&eig_computer)

helper function to build PIROCK algorithm

Note

Without a Shampine’s trick caller, this method is fixed time step.

Template Parameters:

• l – Number of augmented stages (l = 1, 2)

• value_t – Type of coefficients

• alpha_beta_computer_t – Choice of computing of \(\alpha\) and \(\beta\) parameters

• eig_computer_t – Computing of eigenvalues of explicit part (diffusion)

Parameters:

• alpha_beta_computer – \(\alpha\) and \(\beta\) computer object

• eig_computer – Eigenvalue computer of explicit part of the problem

template<std::size_t l = 1, typename value_t = double, typename eig_computer_t>
auto ponio::runge_kutta::pirock::pirock(eig_computer_t &&eig_computer)

helper function to build PIROCK algorithm

Note

Without a Shampine’s trick caller, this method is fixed time step, and without a \(\alpha\) and \(\beta\) computer, parameters are fixed to \(\beta = 0\).

Template Parameters:

• l – Number of augmented stages (l = 1, 2)

• value_t – Type of coefficients

• eig_computer_t – Computing of eigenvalues of explicit part (diffusion)

Parameters:

eig_computer – Eigenvalue computer of explicit part of the problem

template<std::size_t l = 1, typename value_t = double>
auto ponio::runge_kutta::pirock::pirock()

helper function to build PIROCK algorithm

Note

Without a Shampine’s trick caller, this method is fixed time step, without a \(\alpha\) and \(\beta\) computer, parameters are fixed to \(\beta = 0\) and without a eigenvalues computer this method use power method to estimate spectral radius of explicit part (diffusion operator).

Template Parameters:

• l – Number of augmented stages (l = 1, 2)

• value_t – Type of coefficients

Following functions are useful for to build a PIROCK method with \(\ell=2\) and \(\alpha = 1\) (with ponio::runge_kutta::pirock::alpha_fixed computer).

template<typename value_t = double, typename eig_computer_t>
auto ponio::runge_kutta::pirock::pirock_a1(eig_computer_t &&eig_computer)

helper function to build PIROCK algorithm

Note

This function builds a PIROCK algorithm with parameters \(\ell=2\) and \(\alpha = 1\) and specific spectral radius estimator of explicit part (diffusion operator).

Template Parameters:

• value_t – Type of coefficients

• eig_computer_t – Computing of eigenvalues of explicit part (diffusion)

Parameters:

eig_computer – Eigenvalue computer of explicit part of the problem

template<typename value_t = double>
auto ponio::runge_kutta::pirock::pirock_a1()

helper function to build PIROCK algorithm

Note

This function builds a PIROCK algorithm with parameters \(\ell=2\) and \(\alpha = 1\) and power method to estimate spectral radius.

Template Parameters:

value_t – Type of coefficients

Following functions are useful for to build a PIROCK method with \(\ell=1\) and \(\beta = 0\) (with ponio::runge_kutta::pirock::beta_0 computer).

template<typename value_t = double, typename eig_computer_t>
auto ponio::runge_kutta::pirock::pirock_b0(eig_computer_t &&eig_computer)

helper function to build PIROCK algorithm

Note

This function builds a PIROCK algorithm with parameters \(\ell=1\) and \(\beta = 0\) and specific spectral radius estimator of explicit part (diffusion operator).

Template Parameters:

• value_t – Type of coefficients

• eig_computer_t – Computing of eigenvalues of explicit part (diffusion)

Parameters:

eig_computer – Eigenvalue computer of explicit part of the problem

template<typename value_t = double>
auto ponio::runge_kutta::pirock::pirock_b0()

helper function to build PIROCK algorithm

Note

This function builds a PIROCK algorithm with parameters \(\ell=1\) and \(\beta = 0\) and power method to estimate spectral radius.

Template Parameters:

value_t – Type of coefficients

PIROCK for reaction-diffusion-advection problem

In this section we present the implementation of complet PIROCK method (i.e. for reaction-diffusion-advection problem).

template<std::size_t _l, typename alpha_beta_computer_t, typename eig_computer_t, typename _shampine_trick_caller_t = void, bool _is_embedded = false, typename _value_t = double>
class pirock_RDA_impl

implementation of PIROCK method for reaction-diffusion-advection problem

Template Parameters:

• _l – Number of augmented stages (l = 1, 2)

• alpha_beta_computer_t – Choice of computing of \(\alpha\) and \(\beta\) parameters

• eig_computer_t – Computing of eigenvalues of explicit part (diffusion)

• _shampine_trick_caller_t – Computing of Shampine’s trick

• _is_embedded – Set adaptive time step method (default: false)

• value_t – Type of coefficients

Public Functions

pirock_RDA_impl() = default

Construct a new pirock_RDA_impl object.

inline pirock_RDA_impl(alpha_beta_computer_t &&_alpha_beta_computer, eig_computer_t &&_eig_computer)

Construct a new pirock_RDA_impl object.

Parameters:

• _alpha_beta_computer – computer of parameters alpha and beta object that has two member functions which take number of stages (s) l parameter

• _eig_computer – eigenvalue computer functor (that take as argument the function of explicit part, the current time, the current state and the current time step)

template<typename _shampine_trick_caller_t_>
inline pirock_RDA_impl(alpha_beta_computer_t &&_alpha_beta_computer, eig_computer_t &&_eig_computer, _shampine_trick_caller_t_ &&_shampine_trick_caller, value_t a_tol = default_config::tol, value_t r_tol = default_config::tol)

Construct a new pirock_RDA_impl object with Shampine’s trick.

Parameters:

• _alpha_beta_computer – alpha and beta computer object

• _eig_computer – eigenvalue computer functor

• _shampine_trick_caller – Shampine’s trick functor

• a_tol – absolute tolerance

• r_tol – relative tolerance

template<typename problem_t, typename state_t, typename array_ki_t>
inline void operator()(problem_t &pb, value_t &tn, state_t &un, array_ki_t &U, value_t &dt, state_t &u_np1)

iteration of PIROCK_RDA method

Template Parameters:

• problem_t – type of \(f\)

• state_t – type of current state

• array_ki_t – type of temporary stages (13 or 18 stages needed)

Parameters:

• pb – problem with 3 operators: \(F_R\), \(F_D\) and \(F_A\)

• tn – current time

• un – current state

• U – array of temporary stages

• dt – current time step

• u_np1 – solution \(u^{n+1}\) at time \(t^{n+1} = t^n + \Delta t\)

inline auto &info()

gets iteration_info object

inline auto const &info() const

gets iteration_info object (constant version)

Helper functions

User interface functions to build a PIROCK method.

template<std::size_t l = 1, bool is_embedded = false, typename value_t = double, typename alpha_beta_computer_t, typename eig_computer_t, typename shampine_trick_caller_t>
auto ponio::runge_kutta::pirock::pirock_RDA(alpha_beta_computer_t &&alpha_beta_computer, eig_computer_t &&eig_computer, shampine_trick_caller_t &&shampine_trick_caller, value_t absolute_tolerance = default_config::tol, value_t relative_tolerance = default_config::tol)

helper function to build PIROCK algorithm

Template Parameters:

• l – Number of augmented stages (l = 1, 2)

• is_embedded – Set adaptive time step method (default: false)

• value_t – Type of coefficients

• alpha_beta_computer_t – Choice of computing of \(\alpha\) and \(\beta\) parameters

• eig_computer_t – Computing of eigenvalues of explicit part (diffusion)

• shampine_trick_caller_t – Computing of Shampine’s trick

Parameters:

• alpha_beta_computer – \(\alpha\) and \(\beta\) computer object

• eig_computer – Eigenvalue computer of explicit part of the problem

• shampine_trick_caller – Shampine’s trick computer

• absolute_tolerance – Absolute tolerance for adaptive time step method (default: ponio::default_config::tol)

• relatove_tolerance – Relative tolerance for adaptive time step method (default: ponio::default_config::tol)

template<std::size_t l = 1, typename value_t = double, typename alpha_beta_computer_t, typename eig_computer_t>
auto ponio::runge_kutta::pirock::pirock_RDA(alpha_beta_computer_t &&alpha_beta_computer, eig_computer_t &&eig_computer)

helper function to build PIROCK algorithm

Note

Without a Shampine’s trick caller, this method is fixed time step.

Template Parameters:

• l – Number of augmented stages (l = 1, 2)

• value_t – Type of coefficients

• alpha_beta_computer_t – Choice of computing of \(\alpha\) and \(\beta\) parameters

• eig_computer_t – Computing of eigenvalues of explicit part (diffusion)

Parameters:

• alpha_beta_computer – \(\alpha\) and \(\beta\) computer object

• eig_computer – Eigenvalue computer of explicit part of the problem

template<std::size_t l = 1, typename value_t = double, typename eig_computer_t>
auto ponio::runge_kutta::pirock::pirock_RDA(eig_computer_t &&eig_computer)

helper function to build PIROCK algorithm

Note

Without a Shampine’s trick caller, this method is fixed time step, and without a \(\alpha\) and \(\beta\) computer, parameters are fixed to \(\beta = 0\).

Template Parameters:

• l – Number of augmented stages (l = 1, 2)

• value_t – Type of coefficients

• eig_computer_t – Computing of eigenvalues of explicit part (diffusion)

Parameters:

eig_computer – Eigenvalue computer of explicit part of the problem

template<std::size_t l = 1, typename value_t = double>
auto ponio::runge_kutta::pirock::pirock_RDA()

helper function to build PIROCK algorithm

Note

Without a Shampine’s trick caller, this method is fixed time step, without a \(\alpha\) and \(\beta\) computer, parameters are fixed to \(\beta = 0\) and without a eigenvalues computer this method use power method to estimate spectral radius of explicit part (diffusion operator).

Template Parameters:

• l – Number of augmented stages (l = 1, 2)

• value_t – Type of coefficients

Following functions are useful for to build a PIROCK method with \(\ell=2\) and \(\alpha = 1\) (with ponio::runge_kutta::pirock::alpha_fixed computer).

template<typename value_t = double, typename eig_computer_t>
auto ponio::runge_kutta::pirock::pirock_RDA_a1(eig_computer_t &&eig_computer)

helper function to build PIROCK algorithm

Note

This function builds a PIROCK algorithm with parameters \(\ell=2\) and \(\alpha = 1\) and specific spectral radius estimator of explicit part (diffusion operator).

Template Parameters:

• value_t – Type of coefficients

• eig_computer_t – Computing of eigenvalues of explicit part (diffusion)

Parameters:

eig_computer – Eigenvalue computer of explicit part of the problem

template<typename value_t = double>
auto ponio::runge_kutta::pirock::pirock_RDA_a1()

helper function to build PIROCK algorithm

Note

This function builds a PIROCK algorithm with parameters \(\ell=2\) and \(\alpha = 1\) and power method to estimate spectral radius.

Template Parameters:

value_t – Type of coefficients

Following functions are useful for to build a PIROCK method with \(\ell=1\) and \(\beta = 0\) (with ponio::runge_kutta::pirock::beta_0 computer).

template<typename value_t = double, typename eig_computer_t>
auto ponio::runge_kutta::pirock::pirock_RDA_b0(eig_computer_t &&eig_computer)

helper function to build PIROCK algorithm

Note

This function builds a PIROCK algorithm with parameters \(\ell=1\) and \(\beta = 0\) and specific spectral radius estimator of explicit part (diffusion operator).

Template Parameters:

• value_t – Type of coefficients

• eig_computer_t – Computing of eigenvalues of explicit part (diffusion)

Parameters:

eig_computer – Eigenvalue computer of explicit part of the problem

template<typename value_t = double>
auto ponio::runge_kutta::pirock::pirock_RDA_b0()

helper function to build PIROCK algorithm

Note

This function builds a PIROCK algorithm with parameters \(\ell=1\) and \(\beta = 0\) and power method to estimate spectral radius.

Template Parameters:

value_t – Type of coefficients

Bibliography

[Abd02]

Assyr Abdulle. Fourth order chebyshev methods with recurrence relation. SIAM Journal on Scientific Computing, 23(6):2041–2054, 2002. doi:10.1137/S1064827500379549.

[AM01]

Assyr Abdulle and Alexei A. Medovikov. Second order chebyshev methods based on orthogonal polynomials. Numerische Mathematik, 2001. doi:10.1007/s002110100292.

[AV13]

Assyr Abdulle and Gilles Vilmart. Pirock: a swiss-knife partitioned implicit–explicit orthogonal runge–kutta chebyshev integrator for stiff diffusion–advection–reaction problems with or without noise. Journal of Computational Physics, 242:869–888, 2013. doi:10.1016/j.jcp.2013.02.009.

[DDM15]

Stéphane Descombes, Max Duarte, and Marc Massot. Operator splitting methods with error estimator and adaptive time-stepping. Application to the simulation of combustion phenomena. In Roland Glowinski, Stanley Osher, and Wotao Yin, editors, Splitting Methods in Communication, Imaging, Science, and Engineering, pages 1–13. Springer International Publishing, 2015. doi:10.1007/978-3-319-41589-5.

[DP80]

J.R. Dormand and P.J. Prince. A family of embedded runge-kutta formulae. Journal of Computational and Applied Mathematics, 6(1):19–26, 1980. doi:10.1016/0771-050X(80)90013-3.

[Law67]

J. Douglas Lawson. Generalized runge-kutta processes for stable systems with large lipschitz constants. SIAM Journal on Numerical Analysis, 4(3):372–380, 1967. doi:10.1137/0704033.

[MBA14]

Chad D. Meyer, Dinshaw S. Balsara, and Tariq D. Aslam. A stabilized runge–kutta–legendre method for explicit super-time-stepping of parabolic and mixed equations. Journal of Computational Physics, 257:594–626, 2014. doi:10.1016/j.jcp.2013.08.021.

[PD81]

P.J. Prince and J.R. Dormand. High order embedded runge-kutta formulae. Journal of Computational and Applied Mathematics, 7(1):67–75, 1981. doi:10.1016/0771-050X(81)90010-3.