Numerical Integration of the Logistic Equation Using Runge-Kutta Methods

Manipulate[
bifurcation[r0_,r1_,n_] := Quiet@Flatten[Table[limits[[2]][r,n], {r,r0,r1,(r1-r0)/n}], 1];
rmax=limits[[1]];
ListPlot[ bifurcation[rini,Min[endr,rmax],pts], PlotStyle->AbsolutePointSize[0.1],
  AxesOrigin->{rini,ymin}, PlotRange->{{rini,Min[endr,rmax]},{ymin,ymax}}, AxesLabel->{Style["r",12,Italic],Style["(Subscript[y, n])",12]},PlotLabel->Style[Row[{"(Subscript[r, max]) = ", Text[rmax]}] ,18],ImageSize->.65{800,600},ImagePadding->20,Ticks->{Range[rini,If[endr<rmax,(endr),rmax]], Range[0,1,If[ymax-ymin<.3,(ymax-ymin)/10,0.1]]}],
Grid[{{Control[{{limits,limits1,"Method "},{limits1->"Euler",limits2->"RK21",limits3->"RK22",limits4->"RK3"},ControlType->PopupMenu,ControlPlacement->Left}]},{
Control[{{rini,1.5,"r initial"},1.5,2.495,.05,ImageSize->Small,Appearance->"Labeled"}],
Control[{{ymax,1.5,"y max"},.01,1.5,.05,ImageSize->Small,Appearance->"Labeled"}]},
{Control[{{endr,2.5,"r end"},2.5,3.65,.05,ImageSize->Small,Appearance->"Labeled"}],Control[{{ymin,0,"y min"},-0.5,1.5,.05,ImageSize->Small,Appearance->"Labeled"}]}},Alignment->Right],
{{pts,500,"number of points"},500,1000,100,Appearance->"Labeled"}, TrackedSymbols->True,SynchronousUpdating->False,
Initialization:>((* Differential Equations *)
limits1 := {3.0,Quiet@Compile[{{r,_Real},{pts1,_Integer}},Map[{r,#}&, Union[Drop[ NestList[#+r*#*(1-#)&,0.8,pts1],250]]],RuntimeAttributes->{Listable}]};limits2:= {3.65,Quiet@Compile[{{r,_Real},{pts1,_Integer}},Map[{r,#}&, Union[Drop[ NestList[(1/4) #(4-r (-2+r (-1+#)) (-1+#) (-2+r #))&,0.8,pts1],250]]],RuntimeAttributes->{Listable}]};limits3 := {3.4,Quiet@Compile[{{r,_Real},{pts1,_Integer}},Map[{r,#}&, Union[Drop[NestList[(1/2) # (2-r (-1+#) (2+r-2 r #+(r^2) (-1+#) #))&,0.8,pts1],250]]],RuntimeAttributes->{Listable}]};limits4:= {3.55,Quiet@Compile[{{r,_Real},{pts1,_Integer}},Map[{r,#}&, Union[Drop[ NestList[(1/24) #(24+r (-1+#) (-24+r (-12+24 #+r (-4-(-1+#) # (24+r (10-20 #+r (4+(-1+#) #(20+r (4+(-8+r (-1+#)) #)))))))))&,0.8,pts1],250]]],RuntimeAttributes->{Listable}]};
bifurcation[r0_,r1_,n_] := Flatten[Table[limits[[2]][r,n], {r,r0,r1,(r1-r0)/n}], 1]
)
]

The continuous logistic model is described by the differential equation y'=β y (1-y), where β>0 is known as the Malthusian parameter. This Demonstration allows you to integrate this DE numerically with several Runge–Kutta methods with orders from 1 up to 3, with constant step size h and initial value (Subscript[y, 0])=0.8. This procedure illustrates chaos and bifurcation. In particular, by using the Euler method, the map is very similar to the classic logistic map, but other RK methods give different and more interesting maps.

When a Runge–Kutta method is applied to the logistic model, the discrete scheme reads (Subscript[y, n+1])=(Subscript[y, n])+Φ((Subscript[y, n]),r), where r=h β. Then it is possible to explore different areas of the logistic map to analyze its nature by increasing the value of r. For each Runge–Kutta method, there exists a limit value (Subscript[r, max]) where the logistic map is defined in [0, (Subscript[r, max])].

This formulation is modified from that given in [1], p. 941 ff.

Reference

[1] H. Ruskeepää, Mathematica Navigator, 3rd ed., San Diego: Academic Press, 2009.

Contributed by: Luis Rández

Based on programs by: Robert M. Lurie