Understanding Runge-Kutta

(*RKStepValues[Ffn_, p0_, h_]*)

(*A functional alternative to consecutive assignments. *)

Clear@RKStepValues;
RKStepValues[Ffn_,p0_,h_]:=
Module[{p},Fold[Replace[{##},{vs_,sform_}:>
Append[vs,With[{nextpt=p0+(sform/.F[p[i_]]:>vs[[i+1,2]])},
{nextpt,Ffn[nextpt]}]]]&,
{{p0,Ffn[p0]}},
h{1/2 F[p[0]],1/2 F[p[1]], F[p[2]],1/6(F[p[0]]+2 F[p[1]]+2F[p[2]]+F[p[3]])}
]];

(*RKShowSteps[ps_, h_,Ffn_]*)

Clear@RKShowSteps;
RKShowSteps[vs:{{p0_,v0_},{p1_,v1_},{p2_,v2_},{p3_,v3_},{p4_,v4_}},h_]:={{Circle[p0,0.1],{Green,Arrow[{p0,p0+((h v0)/2)}]}},Point[p1],{Green,Arrow[{p1,p1+((h v1)/2)}]},{Pink,Arrow[{p0,p0+((h v1)/2)}]},Point[p2],{Green,Arrow[{p2,p2+((h v2)/2)}]},{Pink,Arrow[{p0,p0+h v2}]},{Green,Arrow[{p3,p3+((h v3)/2)}]},{Orange,Cases[(1/6) {1,2,2,1}Most[vs[[All,2]]],v_:>Arrow[{p0,p0+v}]]},{Orange,Arrow[{p0,p4}]},Circle[p4,0.1`]};

(*PlotNSolve[{{y_ -> (f : InterpolatingFunction[{{t1_, t2}}, ___])}}]*)

Clear[PlotNSolve];PlotNSolve[{{fx:InterpolatingFunction[{{t1_,t2_}},___],fy_InterpolatingFunction}}]:=
ParametricPlot[{fx[t],fy[t]},{t,t1,t2}];
PlotNSolve[_]:={{}};

(*NDSolveODEinRange[Fform_, {x_Symbol, xmin_, xmax_}, {y_Symbol, ymin_, ymax_}, p0_, {t_Symbol, tmin_, tmax_}]*)

Clear@NDSolveODEinRange;
NDSolveODEinRange[Fform_,{x_Symbol,xmin_,xmax_},{y_Symbol,ymin_,ymax_},p0_,{t_Symbol,tmin_,tmax_}]:={x,y}/.NDSolve[Thread[{x'[t],y'[t]}==(Fform/.Thread[{x,y}->{x[t],y[t]}])]~Join~Thread[{x[0],y[0]}==p0],
{x,y},{t,tmin,tmax},
Method->{"EventLocator", "Event":>(xmin ≤x[t] <=xmax && ymin <= y[t] <= ymax),
"EventLocationMethod"->"StepEnd",
"EventAction":> "StopIntegration"
}];

Manipulate[
With[{fform=0.1 x-y,maxh=4},
Replace[{{-1,3},{-2,2}},
prange:{{a_,b_},{c_,d_}}:>
With[{Ffn=(Module[{w},(Function@@{w,{1,fform/.Thread[{x,y}->{w[1],w[2]}]}})/.w[i_]:>w[[i]]]),xrange={a,b},yrange={c,d}},
Graphics[
{{Opacity[0.1],StreamPlot[{1,fform},{x,xrange[[1]],xrange[[2]]},{y,yrange[[1]],yrange[[2]]}][[1]]},
res={{PointSize[0.03],
Replace[{pt,RKStepValues[Ffn,pt,h]},
{{x0_,y0_},vs_}:>
Take[RKShowSteps[vs,h]~Join~
{PlotNSolve[NDSolveODEinRange[{1,fform},{x,xrange[[1]],xrange[[2]]},{y,yrange[[1]],yrange[[2]]},pt,{t,0,h}]
][[1]]},
s]
]
}}}
,Axes->True,AxesLabel->{x,y},PlotRange->prange,
ImageSize->{400,400}
]]]],
{{pt,{0.2,1}},{xrange[[1]],yrange[[1]]},{xrange[[2]],yrange[[2]]},Locator},
{{h,2,"h"},0.01,3,Appearance->"Labeled"},
{{s,1,"steps"},1,12,1,Appearance->"Labeled"},
{{xrange,{-1,3}},ControlType->None},
{{yrange,{-2,2}},ControlType->None},
SynchronousUpdating->False,SaveDefinitions->True,
AutorunSequencing->{2,3}
]

CAPTION

This Demonstration shows the steps involved in computing the Runge–Kutta method of integrating a differential equation and how the approximations behave. in the graphic, drag the locator (from which the calculations start), change the step length h, and move through the steps in the calculation.

DETAILS

Moving the initial point and varying the step size shows how, by sampling from points that contain the expected trajectory, the Runge–Kutta method improves on the Euler and related  methods. Larger step sizes show how the method can err, while decreasing the step size shows the rapid convergence to the actual trajectory.

THIS NOTEBOOK IS THE SOURCE CODE FROM

"Understanding Runge-Kutta" from the Wolfram Demonstrations Project
 http://demonstrations.wolfram.com/UnderstandingRungeKutta/

Contributed by: Gerrard Liddell

A full-function Wolfram Mathematica system (Version 6 or higher) is required to edit this notebook.
GET WOLFRAM MATHEMATICA »

© Wolfram Demonstrations Project & Contributors   |   Terms of Use   |   Make a new version of this Demonstration »