Application of the Poincaré-Bendixson theorem

Ex. 2: Glycolysis, linear analysis (no parameter values assigned)

Clear[diffX,diffY,x,y,a,b,fp]
diffX=-x+a y+x^2 y;
diffY=b-a y- x^2 y;
fp=Solve[{diffX==0,diffY==0},{x,y}];
A={diffX,diffY};
Print[Style["System matrix",18,Blue]]
MatrixForm[A]
Print[Style["Jacobian matrix",18,Blue]]
D[A,{{x,y}}]//MatrixForm (*Jacobian*)
Jfp=D[A,{{x,y}}]/.fp ; (*lin. systems' matrices about f.p.s*)
Print[Style["Jacobian matric about f.p.",18,Blue]]
Map[MatrixForm][Jfp]
Print[Style["Determinant",18,Blue]]
Simplify[Map[Det][Jfp]]
Print[Style["Trace",18,Blue]]
Map[Tr][Jfp]
Print[Style["Eigenvalues",18,Blue]]
Simplify[Map[Eigenvalues][Jfp]]
Print[Style["F.p. coordinates",18,Blue]]
coordinates={x,y}/.fp

Ex. 2: Glycolysis, linear analysis (parameter values assigned)

Clear[diffX,diffY,x,y,a,b]
domx=3;domy=3;

a=0.07;
b=0.53;

diffX=-x+a y+x^2 y;
diffY=b-a y- x^2 y;

nullclineX=Quiet[Solve[diffX==0,y]];
nullclineY=Quiet[Solve[diffY==0,y]];
ncx=Plot[Evaluate[y/.nullclineX],{x,0,domx},PlotStyle->{Gray, Thick}];
ncy=Plot[Evaluate[y/.nullclineY],{x,0,domx},PlotStyle->Gray];
fixedpoints=Quiet[Solve[{diffX==0,diffY==0},{x,y}]];
A={diffX,diffY};
Print[Style["System matrix",18,Blue]]
MatrixForm[A]
Print[Style["Jacobian matrix",18,Blue]]
D[A,{{x,y}}]//MatrixForm (*Jacobian*)
Jfp=D[A,{{x,y}}]/.fixedpoints ; (*lin. systems' matrices about f.p.s*)
Print[Style["Jacobian matric about f.p.s",18,Blue]]
Map[MatrixForm][Jfp]
Print[Style["Determinants",18,Blue]]
Map[Det][Jfp]
Print[Style["Trace",18,Blue]]
Map[Tr][Jfp]
Print[Style["Eigenvalues",18,Blue]]
Map[Eigenvalues][Jfp]
Print[Style["F.p.s coordinates",18,Blue]]
coordinates={x,y}/.fixedpoints
fixedp=Graphics[{EdgeForm[{Thick,Red}],FaceForm[White],Disk[coordinates[[1]][[All]],{0.04,0.04}]}];
sp=StreamPlot[{diffX,diffY},{x,0,domx},{y,0,domy},FrameLabel->{"ADF concentration" x[t], "F6P concentration" y[t]},ImageSize->400,LabelStyle->Large];
Show[sp,ncx,ncy,fixedp]

Ex. 2: Integrated solution of the glycolysis model (interactive)

Clear[diffX,diffY,x,y]
imagesize=350;domx=2;domy=3;maxT=77;

a=0.07;
b=0.53;(*selected so that limit-cycle exists*)

diffX=-x+a y+x^2 y;
diffY=b-a y- x^2 y;

difX=diffX/.{x->x[t],y->y[t]};
difY=diffY/.{x->x[t],y->y[t]};
nullclineX=Solve[diffX==0,y];
nullclineY=Solve[diffY==0,y];
ncx=Plot[Evaluate[y/.nullclineX],{x,0,domx},PlotStyle->{Gray, Thick}];
ncy=Plot[Evaluate[y/.nullclineY],{x,0,domx},PlotStyle->Gray];
splot=StreamPlot[{diffX,diffY},{x,0,domx},{y,0,domy},FrameLabel->{"ADF" x[t], "F6P" y[t]},ImageSize->imagesize];
fixedpoints=Quiet[Solve[{diffX==0,diffY==0},{x,y}]];
coordinates={x,y}/.fixedpoints;
fixedp=Graphics[{EdgeForm[{Thick,Red}],FaceForm[White],Disk[coordinates[[1]][[All]],{0.04,0.06}]}];
Manipulate[
integration=NDSolve[{
x'[t]==difX,
y'[t]==difY,
x[0]==point[[1]],y[0]==point[[2]]},
{x,y},{t,0,T}];
trajectory=ParametricPlot[
Evaluate[{x[t],y[t]}/.integration],
{t,0,T},PlotStyle->Black];
solx=Plot[
Evaluate[x[t]/.integration],
{t,0,T},PlotRange-> All,AxesLabel-> {t,x[t]},PlotLabel->"ADP", ImageSize->imagesize,LabelStyle-> 16,PlotStyle->Black];
soly=Plot[
Evaluate[y[t]/.integration],
{t,0,T},PlotRange-> All,AxesLabel-> {t,y[t]},PlotLabel->"F6P",ImageSize->imagesize,LabelStyle-> 16,PlotStyle->Black];
Grid[{{Show[splot,ncx,ncy,trajectory,fixedp,LabelStyle-> 16],soly},{solx}},Frame->All],
{{T,0.6 maxT,"Integration time, t"},1,maxT,Appearance->"Labeled"},{{point,{1,1}},Locator}, SaveDefinitions->True,TrackedSymbols:>{T,point}]

Ex. 2: Glycolysis, two initial conditions meet at the limit-cycle

Clear[a,b,x,y,tmax]
l={x[t],y[t]};
tmax=70;a=0.07;b=0.53;

diffX=-x+a y+x^2 y;
diffY=b-a y- x^2 y;

difX=diffX/.{x->x[t],y->y[t]};
difY=diffY/.{x->x[t],y->y[t]};
sol1=NDSolve[{
x'[t]==difX,
y'[t]==difY,
x[0]==1.05,y[0]==2.3},
{x,y},{t,0.0,tmax}];
sol2=NDSolve[{
x'[t]==difX,
y'[t]==difY,
x[0]==0.534,y[0]==1.52},
{x,y},{t,0.0,tmax}];

PPsol1=ParametricPlot[Evaluate[{x[t],y[t]}/.sol1],
{t,0.01,tmax},PlotRange->All,PlotPoints->1000,FrameLabel->TraditionalForm/@l,PlotStyle->{Red},PlotRange-> All,LabelStyle->20,ImageSize->400,Frame->True];
PPsol2=ParametricPlot[Evaluate[{x[t],y[t]}/.sol2],
{t,0.01,tmax},PlotRange->All,PlotPoints->1000,FrameLabel->TraditionalForm/@l,PlotStyle->{Blue},PlotRange-> All,LabelStyle->20,ImageSize->400,Frame->True];
fixedpoints=Quiet[Solve[{diffX==0,diffY==0},{x,y}]];
coordinates={x,y}/.fixedpoints;
fixedp=Graphics[{EdgeForm[{Thick,Black}],FaceForm[White],Disk[coordinates[[1]][[All]],{0.025,0.035}]}];

sol1X=Plot[Evaluate[x[t]/.sol1],{t,0,tmax},PlotRange-> All,FrameLabel-> {t,x[t]},PlotLabel->Style["ADP",20], ImageSize->350,LabelStyle->20,PlotStyle->Red,Frame->True];
sol1Y=Plot[Evaluate[y[t]/.sol1],{t,0,tmax},PlotRange-> All,FrameLabel-> {t,y[t]},PlotLabel->Style["F6P",20],ImageSize->350,LabelStyle->20,PlotStyle->Red,Frame->True];
sol2X=Plot[Evaluate[x[t]/.sol2],{t,0,tmax},PlotRange-> All,FrameLabel-> {t,x[t]}, ImageSize->350,LabelStyle->20,PlotStyle->Blue,Frame->True];
sol2Y=Plot[Evaluate[y[t]/.sol2],{t,0,tmax},PlotRange-> All,FrameLabel-> {t,y[t]},ImageSize->350,LabelStyle->20,PlotStyle->Blue,Frame->True];

Grid[{{Show[PPsol1,PPsol2,fixedp,AspectRatio-> 0.8],Show[sol1Y,sol2Y,AspectRatio-> 0.8]},{Show[sol1X,sol2X,AspectRatio-> 0.8]}},Frame-> All]

Ex. 2: Finding the dividing line between repellers and attractors

Clear[a,b]
Solve[-1-a-(b^2)+((2 (b^2))/(a+(b^2)))==0,b]

Ex. 2: Parameter space corresponding to closed orbits

p=Plot[{Sqrt[(1/2) (1-2 a +Sqrt[1-8 a])],Sqrt[(1/2) (1-2 a -Sqrt[1-8 a])]},{a,0,0.17},
PlotStyle->Blue,AxesLabel->{"a","b"},Filling->{1->{2}},LabelStyle->20,ImageSize->Large];
co=Graphics[Text[Style["Closed orbits",Large],{0.05,0.6}]];
sf=Graphics[Text[Style["Stable f.p.",Large],{0.13,0.2}]];
Show[p,co,sf]

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