Runge–Kutta methods

From Wikipedia, the free encyclopedia

In numerical analysis, the Runge–Kutta methods are an important family of implicit and explicit iterative methods for the approximation of solutions of ordinary differential equations. These techniques were developed around 1900 by the German mathematicians C. Runge and M.W. Kutta.

See the article on numerical ordinary differential equations for more background and other methods.

Contents

[edit] The classical fourth-order Runge–Kutta method

One member of the family of Runge–Kutta methods is so commonly used that it is often referred to as "RK4" or simply as "the Runge–Kutta method".

Let an initial value problem be specified as follows.

y' = f(t, y), \quad y(t_0) = y_0

Then, the RK4 method for this problem is given by the following equation:

y_{n+1} = y_n + {h \over 6} (k_1 + 2k_2 + 2k_3 + k_4)

where

k_1 = f \left( t_n, y_n \right)
k_2 = f \left( t_n + {h \over 2}, y_n + {h \over 2} k_1 \right)
k_3 = f \left( t_n + {h \over 2}, y_n + {h \over 2} k_2 \right)
k_4 = f \left( t_n + h, y_n + hk_3 \right)

Thus, the next value (yn+1) is determined by the present value (yn) plus the product of the size of the interval (h) and an estimated slope. The slope is a weighted average of slopes:

  • k1 is the slope at the beginning of the interval;
  • k2 is the slope at the midpoint of the interval, using slope k1 to determine the value of y at the point tn + h/2 using Euler's method;
  • k3 is again the slope at the midpoint, but now using the slope k2 to determine the y-value;
  • k4 is the slope at the end of the interval, with its y-value determined using k3.

In averaging the four slopes, greater weight is given to the slopes at the midpoint:

\mbox{slope} = \frac{k_1 + 2k_2 + 2k_3 + k_4}{6}.

The RK4 method is a fourth-order method, meaning that the error per step is on the order of h5, while the total accumulated error has order h4.

Note that the above formulas are valid for both scalar- and vector-valued functions (i.e., y can be a vector).

[edit] Explicit Runge–Kutta methods

The family of explicit Runge–Kutta methods is a generalization of the RK4 method mentioned above. It is given by

y_{n+1} = y_n + h\sum_{i=1}^s b_i k_i,

where

k_1 = f(t_n, y_n), \,
k_2 = f(t_n+c_2h, y_n+a_{21}hk_1), \,
k_3 = f(t_n+c_3h, y_n+a_{31}hk_1+a_{32}hk_2), \,
\vdots
k_s = f(t_n+c_sh, y_n+a_{s1}hk_1+a_{s2}hk_2+\cdots+a_{s,s-1}hk_{s-1}).

(Note: the above equations have different but equivalent definitions in different texts).

To specify a particular method, one needs to provide the integer s (the number of stages), and the coefficients aij (for 1 ≤ j < is), bi (for i = 1, 2, ..., s) and ci (for i = 2, 3, ..., s). These data are usually arranged in a mnemonic device, known as a Runge–Kutta tableau:

0
c2 a21
c3 a31 a32
\vdots \vdots \ddots
cs as1 as2 \cdots as,s − 1
b1 b2 \cdots bs − 1 bs

The Runge–Kutta method is consistent if

\sum_{j=1}^{i-1} a_{ij} = c_i\ \mathrm{for}\ i=2, \ldots, s.

There are also accompanying requirements if we require the method to have a certain order p, meaning that the truncation error is O(hp+1). These can be derived from the definition of the truncation error itself. For example, a 2-stage method has order 2 if b1 + b2 = 1, b2c2 = 1/2, and b2a21 = 1/2.

[edit] Examples

The RK4 method falls in this framework. Its tableau is:

0
1/2 1/2
1/2 0 1/2
1 0 0 1
1/6 1/3 1/3 1/6

However, the simplest Runge–Kutta method is the (forward) Euler method, given by the formula yn + 1 = yn + hf(tn,yn). This is the only consistent explicit Runge–Kutta method with one stage. The corresponding tableau is:

0
1

An example of a second-order method with two stages is provided by the midpoint method

y_{n+1} = y_n + hf\left(t_n+\frac{h}{2},y_n+\frac{h}{2}f(t_n, y_n)\right).

The corresponding tableau is:

0
1/2 1/2
0 1


[edit] Usage

The following is an example usage of a two-stage explicit Runge–Kutta method, viz.,

0
2/3 2/3
1/4 3/4

to solve the initial-value problem

y' = (\tan{y})+1,\quad y(1)=1,\ t\in [1, 1.1]

with step size h=0.025.

The tableau above yields the equivalent corresponding equations below defining the method:

u_1 = y_n \,
u_2 = y_n + 2/3hf(t_n, u_1) \,
y_{n+1} = y_n + h(1/4f(t_n,u_1)+3/4f(t_n+2/3h,u_2))\,
t0 = 1
y0 = 1
t1 = 1.025
u1 = y0 = 1 f(t0,u1) = 2.557407725 u2 = y0 + 2 / 3hf(t0,u1) = 1.042623462
y1 = y0 + h(1 / 4 * f(t0,u1) + 3 / 4 * f(t0 + 2 / 3h,u2)) = 1.066869388
t2 = 1.05
u1 = y1 = 1.066869388 f(t1,u1) = 2.813524695 u2 = y1 + 2 / 3hf(t1,u1) = 1.113761467
y2 = y1 + h(1 / 4 * f(t1,u1) + 3 / 4 * f(t1 + 2 / 3h,u2)) = 1.141332181
t3 = 1.075
u1 = y2 = 1.141332181 f(t2,u1) = 3.183536647 u2 = y2 + 2 / 3hf(t2,u1) = 1.194391125
y3 = y2 + h(1 / 4 * f(t2,u1) + 3 / 4 * f(t2 + 2 / 3h,u2)) = 1.227417567
t4 = 1.1
u1 = y3 = 1.227417567 f(t3,u1) = 3.796866512 u2 = y3 + 2 / 3hf(t3,u1) = 1.290698676
y4 = y3 + h(1 / 4 * f(t3,u1) + 3 / 4 * f(t3 + 2 / 3h,u2)) = 1.335079087

The numerical solutions correspond to the underlined values. Note that f(ti,u1) has been calculated to avoid recalculation in the yis.

[edit] References

  • George E. Forsythe, Michael A. Malcolm, and Cleve B. Moler. Computer Methods for Mathematical Computations. Englewood Cliffs, NJ: Prentice-Hall, 1977. (See Chapter 6.)
  • Ernst Hairer, Syvert Paul Nørsett, and Gerhard Wanner. Solving ordinary differential equations I: Nonstiff problems, second edition. Berlin: Springer Verlag, 1993. ISBN 3-540-56670-8.
  • William H. Press, Brian P. Flannery, Saul A. Teukolsky, William T. Vetterling. Numerical Recipes in C. Cambridge, UK: Cambridge University Press, 1988. (See Sections 16.1 and 16.2.)
  • Runge-Kutta 4th-order method Notes PPT Matlab Mathematica Maple Mathcad at Holistic Numerical Methods Institute