↜ Back to index Introduction to Numerical Analysis 1

# Lecture 4

In today’s lecture we will look at Runge–Kutta methods:

## Runge-Kutta methods

We consider the ODE

\left\{ \begin{aligned} x'(t) &= f(x(t), t), && t > 0,\\ x(0) &= x_0. \end{aligned} \right.

The Euler method that we covered in the previous lecture has a very slow convergence: the global truncation error is proportional to $$h$$, the time-step. We can improve this significantly by using more function evaluations of the right-hand side $$f(x,t)$$ to approximate the value of the solution $$x$$ at time $$t + h$$ from the value at $$t$$.

This is the idea of Runge-Kutta methods. The extra work that these function evaluations require is usually worth the improved convergence rate and smaller error achieved with fewer time steps (larger $$h$$).

### Midpoint method

Let us set $$t_i = ih$$. In the Euler method we use the value of the right-hand side $$f(x,t)$$ at $$(x_{i-1}, t_{i-1})$$ to find the approximate solution $$x_i$$. However, the function $$f(x(t),t)$$ changes with time and therefore this is not very optimal. Instead, we should try to use the value of $$f(x,t)$$ at the middle of the time step $$f(x(t_{i-1} + \frac h2), t_{i-1} + \frac h2)$$. Unfortunately, we do not know the value $$x(t_{i-1} + \frac h2)$$. But we can approximate it by using the Euler method for time step of length $$\frac h2$$. This gives us $$x(t_{i-1} + \frac h2) \sim x_{i-1} + \frac h2 f(x_{i-1}, t_{i-1})$$ and we obtain the formula for the midpoint method

$x_i = x_{i-1} + h f(x_{i-1} + \frac h2 f(x_{i-1}, t_{i-1}), t_{i-1} + \frac h2)$

See wikipedia for more details. This is an explicit method of order 2.

Exercise: Apply the midpoint method to the ODE

\left\{ \begin{aligned} x' &= x(1+t), && t > 0,\\ x(0) &= 1. \end{aligned} \right.

The exact solution is $$x(t) = \exp(t (t + 2) / 2)$$.

The code is implemented in the file midpoint.f90.

Exercise: Run the code in midpoint.f90 and use gnuplot to plot the resulting solution. Compare it to the exact solution $$\exp(x)$$.

Exercise: Use the code in midpoint-error.f90 to find the global truncation error at $$t = 1$$ for various $$h$$. Use gnuplot to estimate the order of the midpoint method. (Recall the logarithmic scale from last lecture.)

### Classical fourth order Runge-Kutta method

If we perform 4 functional evaluations for each time step, we can increase the order of the method to 4. This leads to the classical fourth order Runge-Kutta method:

\begin{aligned} k_1 &= f(x_{i-1}, t),\\ k_2 &= f(x_{i-1} + \frac h2 k_1, t + \frac h2),\\ k_3 &= f(x_{i-1} + \frac h2 k_2, t + \frac h2),\\ k_4 &= f(x_{i-1} + h k_3, t + h),\\ x_i &= x_{i-1} + \frac h6 (k_1 + 2 k_2 + 2 k_3 + k_4). \end{aligned}

to be figured out...