One dimensional ODEs

This notebook illustrates some Mathematica commands for solving and visualizing one dimensional differential equations of the form x’[t] = f[t,x].   We begin by plotting the slope field in the (t,x)-plane.  Equivalently, we can plot the vectors {1, f[t,x] } with first coordinate 1 and second coordinate f[t,x] using the Mathematica functions VectorPlot or StreamPlot.  The latter seems to look nicer.  Here is the logistic ODE.

f[t_,x_]:= x*(1-x);
VectorPlot[{1,f[t,x]},{t,-1,2},{x,-1,2},FrameLabel->{t,x}]
StreamPlot[{1,f[t,x]},{t,-1,2},{x,-1,2},FrameLabel->{t,x}]

And here is a nonautonomous (but periodic in t) example.  We will give the plot a name so we can use it later.

f[t_,x_]:= 0.5*x*(1-x)+Cos[t]/3;
strplot = StreamPlot[{1,f[t,x]},{t,0,4*Pi},{x,-1,2},StreamPoints->100,FrameLabel->{t,x}]

To find approximate solutions of the differential equation, we can use the functions NDSolve or NDSolveValue.  We need to write down the ODE and the initial condition using the correct syntax.  The output of NDSolve is a substitution rule.  The output of NDSolveValue is the actual solution function x[t].  This cell can be used for other problems by editing the first four lines.

f[t_,x_]:= 0.5*x*(1-x)+Cos[t]/3;
t0=0;
elapsedtime=4*Pi;
x0 = 0;

ode= x'[t]==f[t,x[t]] 
ic = x[t0]==x0 
xsol = NDSolveValue[{ode, ic}, x,{t,t0,t0+elapsedtime}]

solplot= Plot[xsol[t],{t,0,4*Pi},PlotStyle->Red];
Show[strplot,solplot]

Poincaré Maps

Our nonautonomous example is periodic with period 2 Pi.  We will compute the Poincare map which takes an initial state x0 at time t=0 to a final state xPoincare at t=2 Pi.  We define a function which, given x0, solve the ODE and evaluates the solution at time 2 Pi.

xPoincare[x0_] := Module[{xsol, ic = x[0]==x0},
xsol = NDSolveValue[{ode, ic}, x,{t,0,2*Pi}];
xsol[2*Pi]]

Since solutions with x0 < 0 might not exist at t = 2 Pi we will only consider values x0 >= 0.  Here is a plot.

Plot[{x0,xPoincare[x0]},{x0,0,2},PlotRange->{0,2}]

The Poincare map has an attracting fixed point near x0=1.  Here is a way to find it by iterating the Poincare map.

xfixed = Nest[xPoincare,1,20]
xfixed//FullForm

Using this as an initial condition should give a periodic orbit.

ic = x[t0]==xfixed 
xsol = NDSolveValue[{ode, ic}, x,{t,t0,t0+elapsedtime}];
solplot= Plot[xsol[t],{t,0,4*Pi},PlotStyle->Red];
Show[strplot,solplot]

Bifurcation Plots

It is easy to use Mathematica to plot bifurcation diagrams for one-dimensional ODEs depending on a parameter, a.  Just plot the level curve f[x,a] = 0 in the (a,x) plane and add decorations if desired.  Here is an example.

f[x_,a_]:= x^3-x + a
cplot = ContourPlot[f[x,a]==0,{a,-2,2},{x,-2,2},ContourStyle->{Thick, Red}]

To see the flows you could add a stream plot using the vertical vectorfield {0, f[x,a]}.

splot = StreamPlot[{0,f[x,a]},{a,-2,2},{x,-2,2},StreamPoints->30,FrameLabel->{a,x}] 
Show[splot,cplot]