Runge-Kutta Method

Worksheet Problem 1

Here is a basic implementation, following the worksheet from class:

Clear[y];
f[t_,y_]:=t^2+y
y[0]=2;
Δt=1/8;
t[k_]:=k*Δt
Do[
m=f[t[k],y[k]];
n=f[t[k+1/2],y[k]+m*Δt/2];
q=f[t[k+1/2],y[k]+n*Δt/2];
p=f[t[k+1],y[k]+q*Δt];
y[k+1]=y[k]+(m+2n+2q+p)*Δt/6
,{k,0,8}]

Plot the solution :

ListPlot[Table[{t[k],y[k]},{k,0,8}]]

Here is the approximation of y(1):

N[y[8]]

Here is the exact solution:

Clear[y];
sol=DSolve[{y'[t]==t^2+y[t],y[0]==2},y[t],t]

Evaluate the exact solution at t=1:

exact=y[t]/.sol[[1]]/.t->1

N[exact,20]

Error for Δt=1:

N[Abs[exact-281/48]]

Error for Δt = 1/2:

N[Abs[exact-y[2]]]

%22/%33

A Module Approach

Here, we write a Module containing the Runge-Kutta approximation to (dy/dt)=(t^2)+y, y(0)=2. This will allow us to easily change the number of steps and observe the results.

rkApprox[numsteps_]:=Module[{Δt=1/numsteps,y,t},
f[t_,y_]:=t^2 + y;
y[0]=2;  (* replace this with 2. to use floating-point arithmetic *)
t[k_]:=k*Δt;
Do[
(* compute four slopes *)
m=f[t[k],y[k]];
n=f[ t[k+1/2],y[k]+m*Δt/2];
q=f[t[k+1/2],y[k]+n*Δt/2];
p=f[t[k+1],y[k]+q*Δt];
(* now compute the next value of the approximation *)
y[k+1]=y[k]+(m+2n+2q+p)Δt/6
, {k,0,numsteps-1}];
y[numsteps]
]

Compute the approximation with 1 step:

rkApprox[1]

Convert to a decimal with 10 digits:

N[rkApprox[1],10]

Now, let’s make a table of approximations using the following numbers of steps: 1, 2, 4, 8, ..., (2^10):

TableForm[Table[ {2^k,N[rkApprox[2^k],15]},{k,0,10} ],TableHeadings->{None,{"num steps","approximation"}}]

Make another table with the approximations and errors:

TableForm[Table[ {2^k,N[rkApprox[2^k],15],N[Abs[rkApprox[2^k]-exact]]},{k,0,10} ],TableHeadings->{None,{"num steps","approximation","error"}}]

Compute the ratios of successive errors:

errors = Table[N[Abs[exact-rkApprox[2^k]]],{k,0,10} ]

Table[errors[[k]]/errors[[k+1]],{k,1,Length[errors]-1}]

Observe that many of the ratios are close to 16. In general, when you cut the step size in half, the error is reduced by a factor of about 16.