Nonlinear dynamics coursework analysis tools

The coursework consists of two parts. The first part requires you to analyse a 2-D system and the second part requires the analysis of a 3-D system. Completed coursework may be handed over in two parts at any time during the semester.

Part 1: Analysis of a 2-D system

Perform linear analysis:

Find a fixed point or points of your system.

Linearise your system about the fixed point or points.

Plot the linearized phase portrait/s (strong suggestion).

Determine the type of the linearized fixed point or points.

Determine if/how changes in the control parameter values influence the dynamics (type of fixed point/s) of your system.

Perform nonlinear analysis of the full homogeneous system using a computer:

Compare the type of the nonlinear fixed point or points with the corresponding linearized fixed point or points.

Plot the nonlinear phase portrait/s. Compare it/them with the corresponding linearized one/s.

Explain any discrepancies between the nonlinear and linearized systems if any occur.

Perform nonlinear analysis of the non-autonomous system (if applicable) using a computer:

Analyse the influence of explicitly time dependant part of the system on the system dynamics.

Study the nonlinear non-autonomous phase trajectories. Compare them against the corresponding homogeneous phase portrait.

Explain the obtained and presented results.

Draw overall conclusions and comment on the presented results.

Input format

The formulas must be entered according to the Wolfram Mathematica syntax rules, read: https://reference.wolfram.com/language/tutorial/InputAndOutputInNotebooks.html

A short primer on the input rules:

Raising to the power: x^2

Multiplication of two values: a*x or a x (separated by a space)

Addition of two values: a+x

Subtraction of two values: a-x

Entering functions, function argument brackets: Sin[x], Tan[a+x]

Selected functions: Exp[x], Log[x] (natural logarithm), Log10[x] (base-10 logarithm), ArcTan[x], etc.

Linear analysis tool

Clear[Dx1,Dy1,x,y,a,b,fpg,fpx,fpy]
Manipulate[
updatelin;
If[fpg==1,
fp=Quiet[Solve[{Dx1==0,Dy1==0},{x,y}]],
fp={{x->fpx,y->fpy}}];
A={Dx1,Dy1};
(*Jacobian matrix*)
JacobianM=Quiet[Simplify[D[A,{{x,y}}]]]//MatrixForm ;
(*Lin. systems' matrices about f.p.s*)
Jfp=Quiet[Simplify[D[A,{{x,y}}]/.fp]];
LinSystems=Quiet[Map[MatrixForm][Jfp]];
(*Determinants*)
Determinants=Quiet[Simplify[Map[Det][Jfp]]];
(*Traces*)
Traces=Quiet[Simplify[Map[Tr][Jfp]]];
(*Eigenvalues*)
Evalues=Quiet[Simplify[Map[Eigenvalues][Jfp]]];
(*F.p coordinates*)
If[fpg==1,
coordinates=Quiet[Simplify[{x,y}/.fp]],
coordinates=Flatten[Quiet[Simplify[{x,y}/.fp]]]];
(*Degenerate nodes*)
deno=Quiet[Simplify[Traces*Traces-4*Determinants]];
Pane[
Text@Grid[{
{Style[ "Linear analysis results:\n",12,Bold]},
{Style[StringForm["Fixed point/s coordinates  ((SuperStar[x]), (SuperStar[y])) = ``",coordinates],12]},
{Style[StringForm["Jacobian matrix  J = ``",JacobianM],12]},
{Style[StringForm["Jacobian matrix evaluated at fixed point/s  J(Subscript[|, ((SuperStar[x]),(SuperStar[y]))]) = ``",LinSystems],12]},
{Style[StringForm["Eigenvalue/s  (Subscript[λ, 1]), (Subscript[λ, 2])  = ``",Evalues],12]},
{Style[StringForm["Determinant/s  Δ = ``",Determinants],12]},
{Style[StringForm["Trace/s  τ = ``",Traces],12]},
{Style[StringForm["(τ^2)- 4Δ = ``",deno],12]}
},Alignment->Left
],
{750,All}],

Style["Nonlinear system (supports symbolic input):",12,Bold],
{{updatelin,0},ControlType->None},
{{Dx1,-x+a y+x^2 y,"(OverDot[x]) = f(x,y) ="},InputField,ImageSize->Large},
{{Dy1,b-a y-x^2 y , "(OverDot[y]) = g(x,y) ="},InputField,ImageSize->Large},
Row[{Button["Evaluate the system",updatelin=1-updatelin]}],
Delimiter,
Style["Coordinates of the fixed point (in case given):",12,Bold],
Row[{ Control@{{fpx,, "(SuperStar[x]) ="},InputField,ImageSize->70},
Control@{{fpy,, "        (SuperStar[y]) ="},InputField,ImageSize->70} }],
Row[{Style["Fill the above input only if the coordinates are given:",11],Control@{{fpg,1,""},{1-> "Fixed point is NOT given,",2-> "Fixed point IS given."}}}],
ControlPlacement->Top,
SaveDefinitions->True,
TrackedSymbols:>{updatelin,fpg}]

2-D phase portrait plotting tool

Clear[Dx,Dy,maxT,maxX,xmin,xmax,ymin,ymax,natDx,natDy,fpx2,fpy2,fpg2]
Manipulate[
update;
difX=Dx/.{x->x[t],y->y[t]};
difY=Dy/.{x->x[t],y->y[t]};
splot=StreamPlot[{Dx,Dy},{x,xmin,xmax},{y,ymin,ymax},FrameLabel->{x[t],y[t]},ImageSize->Medium,LabelStyle->16];
If[fpg2==1,
fixedpoints=Quiet[Solve[{Dx==0,Dy==0},{x,y},Reals]],
fixedpoints={{x->fpx2,y->fpy2}}];
coordinates={x,y}/.fixedpoints;
fixedp=Graphics[{PointSize[0.025],Red,Point[coordinates]}];
integration=Quiet[NDSolve[{
x'[t]==difX+Evaluate[natDx/.{x->x[t],y->y[t]}],
y'[t]==difY+Evaluate[natDy/.{x->x[t],y->y[t]}],
x[0]==point[[1]],y[0]==point[[2]]},{x,y},{t,0,T}]];
trajectory=Quiet[ParametricPlot[Evaluate[{x[t],y[t]}/.integration],{t,0,T},PlotStyle->Black]];
solx=Quiet[Plot[Evaluate[x[t]/.integration],{t,0,T},PlotRange-> All,AxesLabel-> {t,x[t]},PlotLabel->"Displacement x(t)",ImageSize->Medium,LabelStyle->16,PlotStyle->Black]];
soly=Quiet[Plot[Evaluate[y[t]/.integration],{t,0,T},PlotRange-> All,AxesLabel-> {t,y[t]},PlotLabel->"Speed y(t) = (OverDot[x])(t)",ImageSize->Medium,LabelStyle->16,PlotStyle->Black]];
Grid[{{Show[splot,trajectory,fixedp],soly},{solx}},Frame->All],

Style["Governing equations:",12,Bold],
{{update,0},ControlType->None},
Row[{Control@{{Dx,y,"(OverDot[x]) = f(x,y) ="},InputField,ImageSize->401},
Control@{{natDx,0,"     (non-autonomous term)  +"},InputField,ImageSize->100}}],
Row[{Control@{{Dy,-2 x-y, "(OverDot[y]) = g(x,y) ="},InputField,ImageSize->400},
Control@{{natDy,0.1 x Cos[t],"     (non-autonomous term)  +"},InputField,ImageSize->100}}],
Row[{Button["Evaluate the system",update=1-update] }],
Delimiter,
Style["Coordinates of the fixed point (in case given):",12,Bold],
Row[{ Control@{{fpx2,, "(SuperStar[x]) ="},InputField,ImageSize->70},
Control@{{fpy2,, "        (SuperStar[y]) ="},InputField,ImageSize->70} }],
Row[{Style["Fill the above input only if the coordinates are given:",11],Control@{{fpg2,1,""},{1-> "Fixed point is NOT given,",2-> "Fixed point IS given."}}}],
Delimiter,
Style["Phase portrait and time-domain solution ranges:",12, Bold],
Row[{Control@{{maxX,15,"Maximum phase plane range, |(Subscript[x, min])|, |(Subscript[x, max])|, |(Subscript[y, min])|, |(Subscript[y, max])| ="},InputField,ImageSize->Small}}],
Row[{Control@{{xmin,-maxX/15,"(Subscript[x, min])"},-maxX,xmax-0.05,0.1,Appearance->{"UpArrow","Labeled"}},
Control@{{xmax,maxX/15,"(Subscript[x, max])"},xmin+0.05,maxX,0.1,Appearance->{"DownArrow","Labeled"}}}],
Row[{Control@{{ymin,-maxX/15,"(Subscript[y, min])"},-maxX,ymax-0.05,0.1,Appearance->{"UpArrow","Labeled"}},Control@{{ymax,maxX/15,"(Subscript[y, max])"},ymin+0.05,maxX,0.1,Appearance->{"DownArrow","Labeled"}}}],
Row[{Control@{{maxT,20,"Maximum integration time (Subscript[t, max]) ="},InputField,ImageSize->Small}}],
{{T,0.5 maxT,"(Subscript[t, max])"},0.01,maxT,Appearance->"Labeled",ImageSize->Large},{{point,{0.7,0.7}},Locator},
ControlPlacement->Top,
SaveDefinitions->True,TrackedSymbols:>{update, maxT,maxX,T,point,xmin,xmax,ymin,ymax,fpg2}]

Part 2: Analysis of a 3-D system

Compare the known properties of strange attractors against your system.

Draw conclusions based on your analysis results.

3-D phase trajectory plotting tool

Clear[Dx,Dy,Dz,maxTime,x0,y0,z0,x02,y02,z02,nt,nat3Dx,nat3Dy,nat3Dz,T]
Manipulate[
update3d;
div=Quiet[Div[{Dx,Dy,Dz},{x,y,z}]];
difX=Dx/.{x->x[t],y->y[t],z->z[t]};
difY=Dy/.{x->x[t],y->y[t],z->z[t]};
difZ=Dz/.{x->x[t],y->y[t],z->z[t]};
fixedpoints=Quiet[Solve[{Dx==0,Dy==0,Dz==0},{x,y,z},Reals]];
coordinates={x,y,z}/.fixedpoints;
Gfp=Graphics3D[{PointSize[0.013],Red,Point[coordinates]}];
integ[ic__]:=Quiet[NDSolve[{
x'[t]==difX+nat3Dx,
y'[t]==difY+nat3Dy,
z'[t]==difZ+nat3Dz,
x[0]==ic[[1]],y[0]==ic[[2]],z[0]==ic[[3]]},{x,y,z},{t,0,T},MaxSteps->Infinity]];
dtext=Text@Style[StringForm["The divergence of the (autonomous) vector field  ∇·(OverDot[(OverVector[x])]) = ``\n",div],12];

If[nt==1,
icoord={{x0,y0,z0}};
plot=Quiet[ParametricPlot3D[Evaluate[{x[t],y[t],z[t]}/.integ/@icoord],{t,0,T},PlotPoints->700,AxesLabel->{x[t],y[t],z[t]},LabelStyle->16,PlotRange->All,PlotStyle->Thin, ImageSize->415]];
plotx=Quiet[Plot[Evaluate[x[t]/.integ/@icoord],{t,0,T},Frame->True,FrameLabel->{None,x[t]},LabelStyle->14,PlotRange->All,AspectRatio->0.2, ImageSize->415]];
ploty=Quiet[Plot[Evaluate[y[t]/.integ/@icoord],{t,0,T},Frame->True,FrameLabel->{None,y[t]},LabelStyle->14,PlotRange->All,AspectRatio->0.2, ImageSize->415]];
plotz=Quiet[Plot[Evaluate[z[t]/.integ/@icoord],{t,0,T},Frame->True,FrameLabel->{t,z[t]},LabelStyle->14,PlotRange->All,AspectRatio->0.2, ImageSize->415]];
Grid[{{dtext},{Show[plot,Gfp]},{plotx},{ploty},{plotz}}],

icoord={{x0,y0,z0},{x02,y02,z02}};
plot=Quiet[ParametricPlot3D[Evaluate[{x[t],y[t],z[t]}/.integ/@icoord],{t,0,T},PlotPoints->700,AxesLabel->{x[t],y[t],z[t]},LabelStyle->16,PlotRange->All,PlotStyle->Thin, ImageSize->415]];
plotx=Quiet[Plot[Evaluate[x[t]/.integ/@icoord],{t,0,T},Frame->True,FrameLabel->{None,x[t]},LabelStyle->14,PlotRange->All,AspectRatio->0.2, ImageSize->415]];
ploty=Quiet[Plot[Evaluate[y[t]/.integ/@icoord],{t,0,T},Frame->True,FrameLabel->{None,y[t]},LabelStyle->14,PlotRange->All,AspectRatio->0.2, ImageSize->415]];
plotz=Quiet[Plot[Evaluate[z[t]/.integ/@icoord],{t,0,T},Frame->True,FrameLabel->{t,z[t]},LabelStyle->14,PlotRange->All,AspectRatio->0.2, ImageSize->415]];
Grid[{{dtext},{Show[plot,Gfp]},{plotx},{ploty},{plotz}}]],

Style["Governing equations:",12,Bold],
{{update3d,0},ControlType->None},
Row[{Control@{{Dx,10 (y-x),"(OverDot[x]) = f(x,y,z) ="},InputField,ImageSize->167},
Control@{{nat3Dx,0,"   +"},InputField,ImageSize->76}}],
Row[{Control@{{Dy,28 x-y-x z, "(OverDot[y]) = g(x,y,z) ="},InputField,ImageSize->166},
Control@{{nat3Dy,Pi^4 Cos[Pi t^2],"   +"},InputField,ImageSize->76}}],
Row[{Control@{{Dz,x y-8 z/3, "(OverDot[z]) = h(x,y,z) ="},InputField,ImageSize->166},
Control@{{nat3Dz,0,"   +"},InputField,ImageSize->76}}],
Row[{Button["Evaluate the system",update3d=1-update3d] }],

Delimiter,
Style["Initial conditions and number of trajectories:",12,Bold],
Row[{Control@{{nt,1,""},{1-> "One trajectory, (Subscript[(Overscript[x, ->]), 1])(t)",2-> "Two trajectories, (Subscript[(Overscript[x, ->]), 1])(t) and (Subscript[(Overscript[x, ->]), 2])(t)"}}}],
Row[{Control@{{x0,0.1,"(Subscript[x, 0]) = (Subscript[x, 1])(0) ="},InputField,ImageSize->60},Control@{{x02,0.2,"                    (Subscript[x, 2])(0) ="},InputField,ImageSize->60}}],
Row[{Control@{{y0,0.1,"(Subscript[y, 0]) = (Subscript[y, 1])(0) ="},InputField,ImageSize->60},Control@{{y02,0.2,"                    (Subscript[y, 2])(0) ="},InputField,ImageSize->60}}],
Row[{Control@{{z0,0.1,"(Subscript[z, 0]) = (Subscript[z, 1])(0) ="},InputField,ImageSize->60},Control@{{z02,0.2,"                    (Subscript[z, 2])(0) ="},InputField,ImageSize->60}}],
Row[{Button["Evaluate the initial conditions",update3d=1-update3d] }],

Delimiter,
Style["Integration time range:",12, Bold],
Row[{Control@{{maxTime,30,"Maximum integration time (Subscript[t, max]) ="},InputField,ImageSize->Small}}] ,
{{T,0.5 maxTime,"(Subscript[t, max])"},0.01,maxTime,Appearance->"Labeled",ImageSize->Medium},
ControlPlacement->Left,
SaveDefinitions->True,
TrackedSymbols:>{update3d,maxTime,T,nt}]

Notebook version 3.2 (21.05.2023)
CC-BY 2023, Dmitri Kartofelev, PhD
URL: https://www.ioc.ee/~dima/

Creation of this notebook was supported by IT Academy program of Information Technology Foundation for Education (HITSA).