Runge–Kutta–Fehlberg method
In mathematics, the Runge–Kutta–Fehlberg method (or Fehlberg method) is an algorithm in numerical analysis for the numerical solution of ordinary differential equations. It was developed by the German mathematician Erwin Fehlberg and is based on the large class of Runge–Kutta methods.
The novelty of Fehlberg's method is that it is an embedded method[definition needed] from the Runge–Kutta family, meaning that identical function evaluations are used in conjunction with each other to create methods of varying order and similar error constants. The method presented in Fehlberg's 1969 paper has been dubbed the RKF45 method, and is a method of order O(h4) with an error estimator of order O(h5).[1] By performing one extra calculation, the error in the solution can be estimated and controlled by using the higher-order embedded method that allows for an adaptive stepsize to be determined automatically.
Butcher tableau for Fehlberg's 4(5) method
Any Runge–Kutta method is uniquely identified by its Butcher tableau. The embedded pair proposed by Fehlberg[2]
| 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 coefficients at the bottom of the table gives the fifth-order accurate method, and the second row gives the fourth-order accurate method.
Implementing an RK4(5) Algorithm
The coefficients found by Fehlberg for Formula 1 (derivation with his parameter α2=1/3) are given in the table below, using array indexing of base 1 instead of base 0 to be compatible with most computer languages:
The coefficients in the below table do not work.
| K | A(K) | B(K,L) | C(K) | CH(K) | CT(K) | ||||
|---|---|---|---|---|---|---|---|---|---|
| L=1 | L=2 | L=3 | L=4 | L=5 | |||||
| 1 | 0 | 1/9 | 47/450 | 1/150 | |||||
| 2 | 2/9 | 2/9 | 0 | 0 | 0 | ||||
| 3 | 1/3 | 1/12 | 1/4 | 9/20 | 12/25 | -3/100 | |||
| 4 | 3/4 | 69/128 | -243/128 | 135/64 | 16/45 | 32/225 | 16/75 | ||
| 5 | 1 | -17/12 | 27/4 | -27/5 | 16/15 | 1/12 | 1/30 | 1/20 | |
| 6 | 5/6 | 65/432 | -5/16 | 13/16 | 4/27 | 5/144 | 6/25 | -6/25 | |
Fehlberg[2] outlines a solution to solving a system of n differential equations of the form: [math]\displaystyle{ \frac{dy_i}{dx} = f_i(x,y_1,y_2, \ldots, y_n), i=1,2,\ldots,n }[/math] to iterative solve for [math]\displaystyle{ y_i(x+h), i=1,2,\ldots,n }[/math] where h is an adaptive stepsize to be determined algorithmically:
The solution is the weighted average of six increments, where each increment is the product of the size of the interval, [math]\displaystyle{ h }[/math], and an estimated slope specified by function f on the right-hand side of the differential equation.
[math]\displaystyle{ \begin{align} k_1&=h\cdot f(x+A(1) \cdot h,y) \\ k_2&=h\cdot f(x+A(2)\cdot h,y+B(2,1)\cdot k_1) \\ k_3&=h\cdot f(x+A(3)\cdot h, y+B(3,1)\cdot k_1+B(3,2)\cdot k_2 ) \\ k_4&=h\cdot f(x+A(4)\cdot h, y+B(4,1)\cdot k_1+B(4,2)\cdot k_2+B(4,3)\cdot k_3 ) \\ k_5&=h\cdot f(x+A(5)\cdot h, y+B(5,1)\cdot k_1+B(5,2)\cdot k_2+B(5,3)\cdot k_3+B(5,4)\cdot k_4 ) \\ k_6&=h\cdot f(x+A(6)\cdot h, y+B(6,1)\cdot k_1+B(6,2)\cdot k_2+B(6,3)\cdot k_3+B(6,4)\cdot k_4+B(6,5) \cdot k_5) \end{align} }[/math]
Then the weighted average is:
[math]\displaystyle{ y(x+h)=y(x) + CH(1) \cdot k_1 + CH(2) \cdot k_2 + CH(3) \cdot k_3 + CH(4) \cdot k_4 + CH(5) \cdot k_5 + CH(6) \cdot k_6 }[/math]
The estimate of the truncation error is: [math]\displaystyle{ \mathrm{TE} = \left|\mathrm{CT}(1) \cdot k_1 + \mathrm{CT}(2) \cdot k_2 + \mathrm{CT}(3) \cdot k_3 + \mathrm{CT}(4) \cdot k_4 + \mathrm{CT}(5) \cdot k_5 + \mathrm{CT}(6) \cdot k_6\right| }[/math]
At the completion of the step, a new stepsize is calculated:[3]
[math]\displaystyle{ h_{\text{new}} = 0.9 \cdot h \cdot \left ( \frac {\varepsilon} {TE} \right )^{1/5} }[/math]
If [math]\displaystyle{ \mathrm{TE} \gt \varepsilon }[/math], then replace [math]\displaystyle{ h }[/math] with [math]\displaystyle{ h_{\text{new}} }[/math] and repeat the step. If [math]\displaystyle{ TE\leqslant\varepsilon }[/math], then the step is completed. Replace [math]\displaystyle{ h }[/math] with [math]\displaystyle{ h_{\text{new}} }[/math] for the next step.
The coefficients found by Fehlberg for Formula 2 (derivation with his parameter α2 = 3/8) are given in the table below, using array indexing of base 1 instead of base 0 to be compatible with most computer languages:
| K | A(K) | B(K,L) | C(K) | CH(K) | CT(K) | ||||
|---|---|---|---|---|---|---|---|---|---|
| L=1 | L=2 | L=3 | L=4 | L=5 | |||||
| 1 | 0 | 25/216 | 16/135 | -1/360 | |||||
| 2 | 1/4 | 1/4 | 0 | 0 | 0 | ||||
| 3 | 3/8 | 3/32 | 9/32 | 1408/2565 | 6656/12825 | 128/4275 | |||
| 4 | 12/13 | 1932/2197 | -7200/2197 | 7296/2197 | 2197/4104 | 28561/56430 | 2187/75240 | ||
| 5 | 1 | 439/216 | -8 | 3680/513 | -845/4104 | -1/5 | -9/50 | -1/50 | |
| 6 | 1/2 | -8/27 | 2 | -3544/2565 | 1859/4104 | -11/40 | 2/55 | -2/55 | |
In another table in Fehlberg,[2] coefficients for an RKF4(5) derived by D. Sarafyan are given:
| K | A(K) | B(K,L) | C(K) | CH(K) | CT(K) | ||||
|---|---|---|---|---|---|---|---|---|---|
| L=1 | L=2 | L=3 | L=4 | L=5 | |||||
| 1 | 0 | 0 | 1/6 | 1/24 | 1/8 | ||||
| 2 | 1/2 | 1/2 | 0 | 0 | 0 | ||||
| 3 | 1/2 | 1/4 | 1/4 | 2/3 | 0 | 2/3 | |||
| 4 | 1 | 0 | -1 | 2 | 1/6 | 5/48 | 1/16 | ||
| 5 | 2/3 | 7/27 | 10/27 | 0 | 1/27 | 27/56 | -27/56 | ||
| 6 | 1/5 | 28/625 | -1/5 | 546/625 | 54/625 | -378/625 | 125/336 | -125/336 | |
See also
- List of Runge–Kutta methods
- Numerical methods for ordinary differential equations
- Runge–Kutta methods
Notes
- ↑ According to Hairer et al. (1993, §II.4), the method was originally proposed in Fehlberg (1969); Fehlberg (1970) is an extract of the latter publication.
- ↑ 2.0 2.1 2.2 2.3 2.4 2.5 (Hairer Nørsett) refer to (Fehlberg 1969)
- ↑ Gurevich, Svetlana (2017). "Appendix A Runge-Kutta Methods". pp. 8–11. https://www.uni-muenster.de/imperia/md/content/physik_tp/lectures/ss2017/numerische_Methoden_fuer_komplexe_Systeme_II/rkm-1.pdf.
References
- Fehlberg, Erwin (1968) Classical fifth-, sixth-, seventh-, and eighth-order Runge-Kutta formulas with stepsize control. NASA Technical Report 287. https://ntrs.nasa.gov/api/citations/19680027281/downloads/19680027281.pdf
- Fehlberg, Erwin (1969) Low-order classical Runge-Kutta formulas with stepsize control and their application to some heat transfer problems. Vol. 315. National aeronautics and space administration.
- Fehlberg, Erwin (1969). "Klassische Runge-Kutta-Nystrom-Formeln funfter und siebenter Ordnung mit Schrittweiten-Kontrolle". Computing 4: 93–106. doi:10.1007/BF02234758.
- Fehlberg, Erwin (1970) Some experimental results concerning the error propagation in Runge-Kutta type integration formulas. NASA Technical Report R-352. https://ntrs.nasa.gov/api/citations/19700031412/downloads/19700031412.pdf
- Fehlberg, Erwin (1970). "Klassische Runge-Kutta-Formeln vierter und niedrigerer Ordnung mit Schrittweiten-Kontrolle und ihre Anwendung auf Wärmeleitungsprobleme," Computing (Arch. Elektron. Rechnen), vol. 6, pp. 61–71. doi:10.1007/BF02241732
- Hairer, Ernst; Nørsett, Syvert; Wanner, Gerhard (1993). Solving Ordinary Differential Equations I: Nonstiff Problems (Second ed.). Berlin: Springer-Verlag. ISBN 3-540-56670-8.
- Sarafyan, Diran (1966) Error Estimation for Runge-Kutta Methods Through Pseudo-Iterative Formulas. Technical Report No. 14, Louisiana State University in New Orleans, May 1966.
Further reading
- Fehlberg, E (1958). "Eine Methode zur Fehlerverkleinerung beim Runge-Kutta-Verfahren". Zeitschrift für Angewandte Mathematik und Mechanik 38 (11/12): 421–426. doi:10.1002/zamm.19580381102. Bibcode: 1958ZaMM...38..421F.
- Fehlberg, E (1964). "New high-order Runge-Kutta formulas with step size control for systems of first and second-order differential equations". Zeitschrift für Angewandte Mathematik und Mechanik 44 (S1): T17–T29. doi:10.1002/zamm.19640441310.
- Fehlberg, E (1972). "Klassische Runge-Kutta-Nystrom-Formeln mit Schrittweiten-Kontrolle fur Differentialgleichungen x.. = f(t,x)". Computing 10: 305–315. doi:10.1007/BF02242243.
- Fehlberg, E (1975). "Klassische Runge-Kutta-Nystrom-Formeln mit Schrittweiten-Kontrolle fur Differentialgleichungen x.. = f(t,x,x.)". Computing 14: 371–387. doi:10.1007/BF02253548.
- Simos, T. E. (1993). "A Runge-Kutta Fehlberg method with phase-lag of order infinity for initial-value problems with oscillating solution". Computers & Mathematics with Applications 25 (6): 95–101. doi:10.1016/0898-1221(93)90303-D..
- Handapangoda, C. C.; Premaratne, M.; Yeo, L.; Friend, J. (2008). "Laguerre Runge-Kutta-Fehlberg Method for Simulating Laser Pulse Propagation in Biological Tissue". IEEE Journal of Selected Topics in Quantum Electronics 1 (14): 105–112. doi:10.1109/JSTQE.2007.913971. Bibcode: 2008IJSTQ..14..105H..
- Simos, T. E. (1995). "Modified Runge–Kutta–Fehlberg methods for periodic initial-value problems". Japan Journal of Industrial and Applied Mathematics 12 (1): 109. doi:10.1007/BF03167384..
- Sarafyan, D. (1994). "Approximate Solution of Ordinary Differential Equations and Their Systems Through Discrete and Continuous Embedded Runge-Kutta Formulae and Upgrading Their Order". Computers & Mathematics with Applications 28 (10–12): 353–384. doi:10.1016/0898-1221(94)00201-0.
- Paul, S.; Mondal, S. P.; Bhattacharya, P. (2016). "Numerical solution of Lotka Volterra prey predator model by using Runge–Kutta–Fehlberg method and Laplace Adomian decomposition method". Alexandria Engineering Journal 55 (1): 613–617. doi:10.1016/j.aej.2015.12.026.
 |  |
