Convergence of RK methods¶
We consider the IVP \(u'=\sin[(u+t)^2]\) over \(0\le t \le 4\), with \(u(0)=-1\).
using FundamentalsNumericalComputation
f = (u,p,t) -> sin((t+u)^2)
tspan = (0.0,4.0)
u0 = -1.0
ivp = ODEProblem(f,u0,tspan)
[36mODEProblem[0m with uType [36mFloat64[0m and tType [36mFloat64[0m. In-place: [36mfalse[0m
timespan: (0.0, 4.0)
u0: -1.0
We use a DifferentialEquations solver to construct an accurate approximation to the exact solution.
u_exact = solve(ivp,Tsit5(),reltol=1e-14,abstol=1e-14);
Now we perform a convergence study of our two Runge–Kutta implementations.
n = @. 50*2^(0:5)
err_IE2 = zeros(size(n))
err_RK4 = zeros(size(n))
for (j,n) = enumerate(n)
t,u = FNC.ie2(ivp,n)
err_IE2[j] = maximum( @.abs(u_exact(t)-u) )
t,u = FNC.rk4(ivp,n)
err_RK4[j] = maximum( @.abs(u_exact(t)-u) )
end
pretty_table((n=n,e2=err_IE2,e4=err_RK4),["n","error in IE2","error in RK4"],backend=:html)
| n | error in IE2 | error in RK4 |
|---|---|---|
| 50 | 0.003537842522608692 | 2.07232462955953e-5 |
| 100 | 0.0008914149361871626 | 1.24439566162593e-6 |
| 200 | 0.0002224187237531705 | 7.606554319750103e-8 |
| 400 | 5.556591875843786e-5 | 4.70222111248475e-9 |
| 800 | 1.3887596445827377e-5 | 2.9218255592766695e-10 |
| 1600 | 3.4715876811031166e-6 | 1.8209211916087042e-11 |
The amount of computational work at each time step is assumed to be proportional to the number of stages. Let’s compare on an apples-to-apples basis by using the number of \(f\)-evaluations on the horizontal axis.
using Plots
plot([2n 4n],[err_IE2 err_RK4],m=:o,label=["IE2" "RK4"],
xaxis=(:log10,"f-evaluations"),yaxis=(:log10,"inf-norm error"),
title="Convergence of RK methods",leg=:bottomleft)
plot!(2n,0.01*(n/n[1]).^(-2),l=:dash,label="2nd order")
plot!(4n,1e-6*(n/n[1]).^(-4),l=:dash,label="4th order")
The fourth-order variant is more efficient in this problem over a wide range of accuracy.