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)
ODEProblem with uType Float64 and tType Float64. In-place: false
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) )

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.