Numerical Methods Advanced: Runge-Kutta — Theory
1. Initial value problem
A first-order initial value problem has the form:
\[
\frac{dy}{dx}=f(x,y),
\qquad
y(x_0)=y_0.
\]
The goal is to estimate the value of \(y\) at later points:
\[
x_1=x_0+h,
\qquad
x_2=x_0+2h,
\qquad
\ldots,
\qquad
x_n=x_0+nh.
\]
The number \(h\) is the step size.
2. Why numerical methods are needed
Many differential equations do not have simple exact formulas. Numerical methods approximate
the solution curve point by point.
\[
(x_0,y_0)\to (x_1,y_1)\to (x_2,y_2)\to \cdots \to (x_n,y_n).
\]
More accurate methods try to estimate the shape of the curve inside each step, not only at the starting point.
3. Euler method
Euler’s method uses the slope at the beginning of the step:
\[
k_1=f(x_n,y_n).
\]
Then it updates:
\[
y_{n+1}=y_n+h k_1.
\]
Euler’s method is simple, but the error can grow quickly because it assumes the slope stays constant across the whole step.
4. Heun method
Heun’s method improves Euler’s method by averaging two slopes.
\[
k_1=f(x_n,y_n),
\]
\[
k_2=f(x_n+h,y_n+hk_1).
\]
Then:
\[
y_{n+1}=y_n+\frac{h}{2}(k_1+k_2).
\]
This is also called the improved Euler method.
5. Midpoint method
The midpoint method uses a slope at the middle of the interval.
\[
k_1=f(x_n,y_n),
\]
\[
k_2=f\left(x_n+\frac{h}{2},y_n+\frac{h}{2}k_1\right).
\]
Then:
\[
y_{n+1}=y_n+h k_2.
\]
Midpoint is usually more accurate than Euler because it uses information from inside the step.
6. Fourth-order Runge-Kutta method
The classical fourth-order Runge-Kutta method, usually called RK4, uses four slopes:
\[
k_1=f(x_n,y_n),
\]
\[
k_2=f\left(x_n+\frac{h}{2},y_n+\frac{h k_1}{2}\right),
\]
\[
k_3=f\left(x_n+\frac{h}{2},y_n+\frac{h k_2}{2}\right),
\]
\[
k_4=f(x_n+h,y_n+h k_3).
\]
The update formula is:
\[
y_{n+1}=y_n+\frac{h}{6}(k_1+2k_2+2k_3+k_4).
\]
The middle slopes \(k_2\) and \(k_3\) receive double weight because they usually represent the average behavior
of the solution across the interval better than the endpoint slopes alone.
7. Why RK4 is more accurate
Euler’s method uses one slope:
\[
y_{n+1}=y_n+h f(x_n,y_n).
\]
RK4 uses a weighted average of four slopes:
\[
\frac{k_1+2k_2+2k_3+k_4}{6}.
\]
This makes RK4 much better at following curved solution paths.
For many smooth problems, RK4 gives high accuracy even when Euler visibly drifts.
8. Error order
The global error of a method describes how the final error changes when the step size \(h\) is reduced.
9. Half-step error estimate
A practical way to estimate RK4 error is to compare one full step with two half steps.
\[
y_{\text{full}}=\text{RK4 step using }h,
\]
\[
y_{\text{half}}=\text{two RK4 steps using }\frac{h}{2}.
\]
For RK4, an error estimate is:
\[
E_{\text{RK4}}\approx \frac{|y_{\text{half}}-y_{\text{full}}|}{15}.
\]
This is not the same as exact error, but it gives a useful indication of local numerical accuracy.
10. Worked example: \(y'=x+y,\ y(0)=1\)
Consider:
\[
\frac{dy}{dx}=x+y,
\qquad
y(0)=1.
\]
This equation has exact solution:
\[
y=2e^x-x-1.
\]
Using \(h=0.1\), Euler and RK4 both advance from \(x=0\) to \(x=1\), but RK4 follows the exact curve much more closely.
11. Accuracy comparison graphs
A solution graph compares the numerical curves:
- Euler curve,
- Heun curve,
- midpoint curve,
- RK4 curve,
- exact curve if available.
An error graph compares:
\[
|E_{\text{Euler}}|=|y_{\text{Euler}}-y_{\text{exact}}|,
\qquad
|E_{\text{RK4}}|=|y_{\text{RK4}}-y_{\text{exact}}|.
\]
When no exact curve is available, a RK4 step estimate can still be shown.
13. Common mistakes
- Using \(k_i\) as increments instead of slopes: in this convention, \(k_i\) are slopes and the factor \(h\) is outside the weighted average.
- Forgetting the initial condition: every numerical curve must start at \((x_0,y_0)\).
- Using too large a step size: even RK4 can become inaccurate if \(h\) is too large.
- Comparing methods with different \(h\): use the same step size when comparing Euler and RK4 fairly.
- Assuming numerical means exact: RK4 is accurate, but it is still an approximation unless the exact solution is known.
- Ignoring units: graph axes should show both numbered values and units.
