List of Runge–Kutta methods

Runge–Kutta methods are methods for the numerical solution of the ordinary differential equation

${\displaystyle {\frac {dy}{dt}}=f(t,y)}$

which take the form

${\displaystyle y_{n+1}=y_{n}+h\sum _{i=1}^{s}b_{i}k_{i}}$
${\displaystyle k_{1}=f(t_{n},y_{n}),}$
${\displaystyle k_{2}=f(t_{n}+c_{2}h,y_{n}+h(a_{21}k_{1})),}$
${\displaystyle k_{3}=f(t_{n}+c_{3}h,y_{n}+h(a_{31}k_{1}+a_{32}k_{2})),}$
${\displaystyle \vdots }$
${\displaystyle k_{i}=f\left(t_{n}+c_{i}h,y_{n}+h\sum _{j=1}^{i-1}a_{ij}k_{j}\right),}$

Each method listed on this page is defined by its Butcher tableau, which puts the coefficients of the method in a table as follows:

${\displaystyle {\begin{array}{c|cccc}c_{1}&a_{11}&a_{12}&\dots &a_{1s}\\c_{2}&a_{21}&a_{22}&\dots &a_{2s}\\\vdots &\vdots &\vdots &\ddots &\vdots \\c_{s}&a_{s1}&a_{s2}&\dots &a_{ss}\\\hline &b_{1}&b_{2}&\dots &b_{s}\\\end{array}}}$

Explicit methods

The explicit methods are those where the matrix ${\displaystyle [a_{ij}]}$ is lower triangular.

Forward Euler

The Euler method is first order. The lack of stability and accuracy limits its popularity mainly to use as a simple introductory example of a numeric solution method.

${\displaystyle {\begin{array}{c|c}0&0\\\hline &1\\\end{array}}}$

Explicit midpoint method

The (explicit) midpoint method is a second-order method with two stages (see also the implicit midpoint method below):

${\displaystyle {\begin{array}{c|cc}0&0&0\\1/2&1/2&0\\\hline &0&1\\\end{array}}}$

Heun's method

Heun's method is a second-order method with two stages (also known as explicit trapezoid rule):

${\displaystyle {\begin{array}{c|cc}0&0&0\\1&1&0\\\hline &1/2&1/2\\\end{array}}}$

Ralston's method

Ralston's method is a second-order method with two stages and a minimum local error bound:

${\displaystyle {\begin{array}{c|cc}0&0&0\\2/3&2/3&0\\\hline &1/4&3/4\\\end{array}}}$

Generic second-order method

${\displaystyle {\begin{array}{c|ccc}0&0&0\\x&x&0\\\hline &1-{\frac {1}{2x}}&{\frac {1}{2x}}\\\end{array}}}$

Kutta's third-order method

${\displaystyle {\begin{array}{c|ccc}0&0&0&0\\1/2&1/2&0&0\\1&-1&2&0\\\hline &1/6&2/3&1/6\\\end{array}}}$

Classic fourth-order method

The "original" Runge–Kutta method.

${\displaystyle {\begin{array}{c|cccc}0&0&0&0&0\\1/2&1/2&0&0&0\\1/2&0&1/2&0&0\\1&0&0&1&0\\\hline &1/6&1/3&1/3&1/6\\\end{array}}}$

3/8-rule fourth-order method

This method doesn't have as much notoriety as the "classical" method, but is just as classical because it was proposed in the same paper (Kutta, 1901).

${\displaystyle {\begin{array}{c|cccc}0&0&0&0&0\\1/3&1/3&0&0&0\\2/3&-1/3&1&0&0\\1&1&-1&1&0\\\hline &1/8&3/8&3/8&1/8\\\end{array}}}$

Embedded methods

The embedded methods are designed to produce an estimate of the local truncation error of a single Runge-Kutta step, and as result, allow to control the error with adaptive stepsize. This is done by having two methods in the tableau, one with order p and one with order p-1.

The lower-order step is given by

${\displaystyle y_{n+1}^{*}=y_{n}+h\sum _{i=1}^{s}b_{i}^{*}k_{i},}$

where the ${\displaystyle k_{i}}$ are the same as for the higher order method. Then the error is

${\displaystyle e_{n+1}=y_{n+1}-y_{n+1}^{*}=h\sum _{i=1}^{s}(b_{i}-b_{i}^{*})k_{i},}$

which is ${\displaystyle O(h^{p})}$. The Butcher Tableau for this kind of method is extended to give the values of ${\displaystyle b_{i}^{*}}$

${\displaystyle {\begin{array}{c|cccc}c_{1}&a_{11}&a_{12}&\dots &a_{1s}\\c_{2}&a_{21}&a_{22}&\dots &a_{2s}\\\vdots &\vdots &\vdots &\ddots &\vdots \\c_{s}&a_{s1}&a_{s2}&\dots &a_{ss}\\\hline &b_{1}&b_{2}&\dots &b_{s}\\&b_{1}^{*}&b_{2}^{*}&\dots &b_{s}^{*}\\\end{array}}}$

Heun–Euler

The simplest adaptive Runge–Kutta method involves combining Heun's method, which is order 2, with the Euler method, which is order 1. Its extended Butcher Tableau is:

${\displaystyle {\begin{array}{c|cc}0&\\1&1\\\hline &1/2&1/2\\&1&0\end{array}}}$

The error estimate is used to control the stepsize.

Fehlberg RK1(2)

The Fehlberg method[1] has two methods of orders 1 and 2. Its extended Butcher Tableau is:

 0 1/2 1/2 1 1/256 255/256 1/256 255/256 0 1/512 255/256 1/512

The first row of b coefficients gives the first-order accurate solution, and the second row has order two.

Bogacki–Shampine

The Bogacki–Shampine method has two methods of orders 3 and 2. Its extended Butcher Tableau is:

 0 1/2 1/2 3/4 0 3/4 1 2/9 1/3 4/9 2/9 1/3 4/9 0 7/24 1/4 1/3 1/8

The first row of b coefficients gives the third-order accurate solution, and the second row has order two.

Fehlberg

The Runge–Kutta–Fehlberg method has two methods of orders 5 and 4. Its extended Butcher Tableau is:

 0 1/4 1/4 3/8 3/32 9/32 12/13 1932/2197 −7200/2197 7296/2197 1 439/216 −8 3680/513 −845/4104 1/2 -8/27 2 −3544/2565 1859/4104 −11/40 16/135 0 6656/12825 28561/56430 −9/50 2/55 25/216 0 1408/2565 2197/4104 −1/5 0

The first row of b coefficients gives the fifth-order accurate solution, and the second row has order four.

Cash-Karp

Cash and Karp have modified Fehlberg's original idea. The extended tableau for the Cash–Karp method is

 0 1/5 1/5 3/10 3/40 9/40 3/5 3/10 −9/10 6/5 1 −11/54 5/2 −70/27 35/27 7/8 1631/55296 175/512 575/13824 44275/110592 253/4096 37/378 0 250/621 125/594 0 512/1771 2825/27648 0 18575/48384 13525/55296 277/14336 1/4

The first row of b coefficients gives the fifth-order accurate solution, and the second row has order four.

Dormand–Prince

The extended tableau for the Dormand–Prince method is

 0 1/5 1/5 3/10 3/40 9/40 4/5 44/45 −56/15 32/9 8/9 19372/6561 −25360/2187 64448/6561 −212/729 1 9017/3168 −355/33 46732/5247 49/176 −5103/18656 1 35/384 0 500/1113 125/192 −2187/6784 11/84 35/384 0 500/1113 125/192 −2187/6784 11/84 0 5179/57600 0 7571/16695 393/640 −92097/339200 187/2100 1/40

The first row of b coefficients gives the fifth-order accurate solution and the second row gives the fourth-order accurate solution.

Implicit methods

Backward Euler

The backward Euler method is first order. Unconditionally stable and non-oscillatory for linear diffusion problems.

${\displaystyle {\begin{array}{c|c}1&1\\\hline &1\\\end{array}}}$

Implicit midpoint

The implicit midpoint method is of second order. It is the simplest method in the class of collocation methods known as the Gauss methods. It is a symplectic integrator.

${\displaystyle {\begin{array}{c|c}1/2&1/2\\\hline &1\end{array}}}$

Gauss–Legendre methods

These methods are based on the points of Gauss–Legendre quadrature. The Gauss–Legendre method of order four has Butcher tableau:

${\displaystyle {\begin{array}{c|cc}{\frac {1}{2}}-{\frac {\sqrt {3}}{6}}&{\frac {1}{4}}&{\frac {1}{4}}-{\frac {\sqrt {3}}{6}}\\{\frac {1}{2}}+{\frac {\sqrt {3}}{6}}&{\frac {1}{4}}+{\frac {\sqrt {3}}{6}}&{\frac {1}{4}}\\\hline &{\frac {1}{2}}&{\frac {1}{2}}\\&{\frac {1}{2}}+{\frac {1}{2}}{\sqrt {3}}&{\frac {1}{2}}-{\frac {1}{2}}{\sqrt {3}}\\\end{array}}}$

The Gauss–Legendre method of order six has Butcher tableau:

${\displaystyle {\begin{array}{c|ccc}{\frac {1}{2}}-{\frac {\sqrt {15}}{10}}&{\frac {5}{36}}&{\frac {2}{9}}-{\frac {\sqrt {15}}{15}}&{\frac {5}{36}}-{\frac {\sqrt {15}}{30}}\\{\frac {1}{2}}&{\frac {5}{36}}+{\frac {\sqrt {15}}{24}}&{\frac {2}{9}}&{\frac {5}{36}}-{\frac {\sqrt {15}}{24}}\\{\frac {1}{2}}+{\frac {\sqrt {15}}{10}}&{\frac {5}{36}}+{\frac {\sqrt {15}}{30}}&{\frac {2}{9}}+{\frac {\sqrt {15}}{15}}&{\frac {5}{36}}\\\hline &{\frac {5}{18}}&{\frac {4}{9}}&{\frac {5}{18}}\\&-{\frac {5}{6}}&{\frac {8}{3}}&-{\frac {5}{6}}\end{array}}}$

Lobatto methods

There are three main families of Lobatto methods, called IIIA, IIIB and IIIC (in classial mathematical literature, the symbols I and II are reserved for two types of Radau methods). These are named after Rehuel Lobatto. All are implicit methods, have order 2s − 2 and they all have c1 = 0 and cs = 1. Unlike any explicit method, it's possible for these methods to have the order greater than the number of stages. Lobatto lived before the classic fourth-order method was popularized by Runge and Kutta.

Lobatto IIIA methods

The Lobatto IIIA methods are collocation methods. The second-order method is known as the trapezoidal rule:

${\displaystyle {\begin{array}{c|cc}0&0&0\\1&1/2&1/2\\\hline &1/2&1/2\\&1&0\\\end{array}}}$

The fourth-order method is given by

${\displaystyle {\begin{array}{c|ccc}0&0&0&0\\1/2&5/24&1/3&-1/24\\1&1/6&2/3&1/6\\\hline &1/6&2/3&1/6\\&-{\frac {1}{2}}&2&-{\frac {1}{2}}\\\end{array}}}$

This methods are A-stable, but not L-stable and B-stable.

Lobatto IIIB methods

The Lobatto IIIB methods are not collocation methods, but they can be viewed as discontinuous collocation methods (Hairer, Lubich & Wanner 2006, §II.1.4). The second-order method is given by

${\displaystyle {\begin{array}{c|cc}0&1/2&0\\1&1/2&0\\\hline &1/2&1/2\\&1&0\\\end{array}}}$

The fourth-order method is given by

${\displaystyle {\begin{array}{c|ccc}0&1/6&-1/6&0\\1/2&1/6&1/3&0\\1&1/6&5/6&0\\\hline &1/6&2/3&1/6\\&-{\frac {1}{2}}&2&-{\frac {1}{2}}\\\end{array}}}$

Lobatto IIIB methods are A-stable, but not L-stable and B-stable.

Lobatto IIIC methods

The Lobatto IIIC methods also are discontinuous collocation methods. The second-order method is given by

${\displaystyle {\begin{array}{c|cc}0&1/2&-1/2\\1&1/2&1/2\\\hline &1/2&1/2\\&1&0\\\end{array}}}$

The fourth-order method is given by

${\displaystyle {\begin{array}{c|ccc}0&1/6&-1/3&1/6\\1/2&1/6&5/12&-1/12\\1&1/6&2/3&1/6\\\hline &1/6&2/3&1/6\\&-{\frac {1}{2}}&2&-{\frac {1}{2}}\\\end{array}}}$

They are L-stable. They are also algebraically stable and thus B-stable, that makes them suitable for stiff problems.

Lobatto IIIC* methods

The Lobatto IIIC* methods are also known as Lobatto III methods (Butcher, 2008), Butcher’s Lobatto methods (Hairer et al, 1993), and Lobatto IIIC methods (Sun, 2000) in the literature.[2] The second-order method is given by

${\displaystyle {\begin{array}{c|cc}0&0&0\\1&1&0\\\hline &1/2&1/2\\\end{array}}}$

The fourth-order method is given by

${\displaystyle {\begin{array}{c|ccc}0&0&0&0\\1/2&1/4&1/4&0\\1&0&1&0\\\hline &1/6&2/3&1/6\\\end{array}}}$

These methods are not A-stable, B-stable or L-stable. The Lobatto IIIC* method for ${\displaystyle s=2}$ is sometimes called the explicit trapezoidal rule.

Generalized Lobatto methods

One can consider a very general family of methods with three real parameters ${\displaystyle (\alpha _{A},\alpha _{B},\alpha _{C})}$ by considering Lobatto coefficients of the form

${\displaystyle a_{i,j}(\alpha _{A},\alpha _{B},\alpha _{C})=\alpha _{A}a_{i,j}^{A}+\alpha _{B}a_{i,j}^{B}+\alpha _{C}a_{i,j}^{C}+\alpha _{C*}a_{i,j}^{C*}}$,

where

${\displaystyle \alpha _{C*}=1-\alpha _{A}-\alpha _{B}-\alpha _{C}}$.

For example, Lobatto IIID family introduced in (Nørsett and Wanner, 1981), also called Lobatto IIINW, are given by

${\displaystyle {\begin{array}{c|cc}0&1/2&1/2\\1&-1/2&1/2\\\hline &1/2&1/2\\\end{array}}}$

and

${\displaystyle {\begin{array}{c|ccc}0&1/6&0&-1/6\\1/2&1/12&5/12&0\\1&1/2&1/3&1/6\\\hline &1/6&2/3&1/6\\\end{array}}}$

These methods correspond to ${\displaystyle \alpha _{A}=2}$, ${\displaystyle \alpha _{B}=2}$, ${\displaystyle \alpha _{C}=-1}$, and ${\displaystyle \alpha _{C*}=-2}$. The methods are L-stable. They are algebraically stable and thus B-stable.

Radau methods are fully implicit methods (matrix A of such methods can have any structure). Radau methods attain order 2s − 1 for s stages. Radau methods are A-stable, but expensive to implement. Also they can suffer from order reduction. The first order Radau method is similar to backward Euler method.

The third-order method is given by

${\displaystyle {\begin{array}{c|cc}0&1/4&-1/4\\2/3&1/4&5/12\\\hline &1/4&3/4\\\end{array}}}$

The fifth-order method is given by

${\displaystyle {\begin{array}{c|ccc}0&{\frac {1}{9}}&{\frac {-1-{\sqrt {6}}}{18}}&{\frac {-1+{\sqrt {6}}}{18}}\\{\frac {3}{5}}-{\frac {\sqrt {6}}{10}}&{\frac {1}{9}}&{\frac {11}{45}}+{\frac {7{\sqrt {6}}}{360}}&{\frac {11}{45}}-{\frac {43{\sqrt {6}}}{360}}\\{\frac {3}{5}}+{\frac {\sqrt {6}}{10}}&{\frac {1}{9}}&{\frac {11}{45}}+{\frac {43{\sqrt {6}}}{360}}&{\frac {11}{45}}-{\frac {7{\sqrt {6}}}{360}}\\\hline &{\frac {1}{9}}&{\frac {4}{9}}+{\frac {\sqrt {6}}{36}}&{\frac {4}{9}}-{\frac {\sqrt {6}}{36}}\\\end{array}}}$

The ci of this method are zeros of

${\displaystyle P_{s}(2x-1)-P_{s-1}(2x-1)=0,}$

where ${\displaystyle P_{s}}$ is the Legendre polynomial of degree s. The third-order method is given by

${\displaystyle {\begin{array}{c|cc}1/3&5/12&-1/12\\1&3/4&1/4\\\hline &3/4&1/4\\\end{array}}}$

The fifth-order method is given by

${\displaystyle {\begin{array}{c|ccc}{\frac {2}{5}}-{\frac {\sqrt {6}}{10}}&{\frac {11}{45}}-{\frac {7{\sqrt {6}}}{360}}&{\frac {37}{225}}-{\frac {169{\sqrt {6}}}{1800}}&-{\frac {2}{225}}+{\frac {\sqrt {6}}{75}}\\{\frac {2}{5}}+{\frac {\sqrt {6}}{10}}&{\frac {37}{225}}+{\frac {169{\sqrt {6}}}{1800}}&{\frac {11}{45}}+{\frac {7{\sqrt {6}}}{360}}&-{\frac {2}{225}}-{\frac {\sqrt {6}}{75}}\\1&{\frac {4}{9}}-{\frac {\sqrt {6}}{36}}&{\frac {4}{9}}+{\frac {\sqrt {6}}{36}}&{\frac {1}{9}}\\\hline &{\frac {4}{9}}-{\frac {\sqrt {6}}{36}}&{\frac {4}{9}}+{\frac {\sqrt {6}}{36}}&{\frac {1}{9}}\\\end{array}}}$

References

• Hairer, Ernst; Nørsett, Syvert Paul; Wanner, Gerhard (1993), Solving ordinary differential equations I: Nonstiff problems, Berlin, New York: Springer-Verlag, ISBN 978-3-540-56670-0.
• Hairer, Ernst; Wanner, Gerhard (1996), Solving ordinary differential equations II: Stiff and differential-algebraic problems, Berlin, New York: Springer-Verlag, ISBN 978-3-540-60452-5.
• Hairer, Ernst; Lubich, Christian; Wanner, Gerhard (2006), Geometric Numerical Integration: Structure-Preserving Algorithms for Ordinary Differential Equations (2nd ed.), Berlin, New York: Springer-Verlag, ISBN 978-3-540-30663-4.