Inverse Problem 3D (Ch 13, § 13.2.2 page 215)

Clear["Global`*"]
Off[General::spell1]

<<VectorFieldPlots`;

Let' s consider the given flow curvature manifold of an unknown two-dimensional polynomial dynamical system.

Flow Curvature Manifold 1

(Subscript[ϕ, 1])(x,y, z)=σ (-((r^2) (x^3)+(x ((y^2)+(z^2))-y z (-1+β)) ((x^2)+β)-r x (2 (x^2) z-z (-1+β) β+x y (1+β))-r (x-y) z β σ+(x-y) ((y^2)+(z^2) β) σ) (-y+x (r+x y-z (1+β))+(-y+(2 x-y) (r-z)) σ+(x-y) (σ^2))+(x ((r^2) x+x (y^2)-r (x+y+2 x z)+x z (1+z+β)+y (1-(x^2)+z-z β))+(x-y) y (-1+r-z) σ) ((x^3) y+z (β^3)+(y^2) (1+β) σ+(y^2) (σ^2)-(x^2) z (1+2 β+σ)-x y (1+β+(β^2)-2 z σ+β σ+(σ^2))+r x (-2 y σ+x (1+β+σ)))+(-(x^3) z-(x^2) y (1+β)+r ((x^3)-x z β)+y (-(y^2) σ+z β (1-β+σ))+x (z β (z+β-σ)+(y^2) (β+σ))) ((r^2) x σ-x (y^2) σ-x z (1-(x^2)+β+(β^2)+σ-z σ+(σ^2))+y (-1+z σ (2+σ)+(x^2) (2+β+σ))+r (-(x^3)-y σ (2+σ)+x (1+σ-2 z σ+(σ^2)))))

By using the symmetry of this flow curvature manifold: ϕ (-x, -y, z) = ϕ (x, y, z)

This implies that: (Subscript[b, 020])=0 ;  (Subscript[b, 110]) =0 ; (Subscript[b, 200]) =0 ; (Subscript[c, 010]) =0 ; (Subscript[c, 100]) =0. Thus, the dynamical system may be written as:

System

Parameters

f[x_,y_,z_]:=(Subscript[a, 100])x+ (Subscript[a, 010])y

g[x_,y_,z_]:=(Subscript[b, 100])x+ (Subscript[b, 010])y+(Subscript[b, 101])x z

h[x_,y_,z_]:= (Subscript[c, 110])x y +(Subscript[c, 001])z

Vector field

V={f[x,y,z],g[x,y,z],h[x,y,z]}

Fixed Points

fp=Solve[{f[(Subscript[x, 1]),(Subscript[y, 1]),(Subscript[z, 1])]==0,g[(Subscript[x, 1]),(Subscript[y, 1]),(Subscript[z, 1])]==0,h[(Subscript[x, 1]),(Subscript[y, 1]),(Subscript[z, 1])]==0},{(Subscript[x, 1]),(Subscript[y, 1]),(Subscript[z, 1])}]

fp1=Sort[fp[[1]]];fp2=Sort[fp[[2]]];fp3=Sort[fp[[3]]];

{(Subscript[x, 1]),(Subscript[y, 1]),(Subscript[z, 1])}={Extract[fp1,{1}][[2]],Extract[fp1,{2}][[2]],Extract[fp1,{3}][[2]]};
fixedpoint1=Show[Graphics3D[{PointSize[0.025],Hue[0.33],Point[{(Subscript[x, 1]),(Subscript[y, 1]),(Subscript[z, 1])}]},PlotRange->{{(Subscript[x, min]),(Subscript[x, max])},{(Subscript[y, min]),(Subscript[y, max])},{(Subscript[z, min]),(Subscript[z, max])}}]];

{(Subscript[x, 2]),(Subscript[y, 2]),(Subscript[z, 2])}={Extract[fp2,{1}][[2]],Extract[fp2,{2}][[2]],Extract[fp2,{3}][[2]]};
fixedpoint2=Show[Graphics3D[{PointSize[0.025],Hue[0.33],Point[{(Subscript[x, 2]),(Subscript[y, 2]),(Subscript[z, 2])}]},PlotRange->{{(Subscript[x, min]),(Subscript[x, max])},{(Subscript[y, min]),(Subscript[y, max])},{(Subscript[z, min]),(Subscript[z, max])}}]];

{(Subscript[x, 3]),(Subscript[y, 3]),(Subscript[z, 3])}={Extract[fp3,{1}][[2]],Extract[fp3,{2}][[2]],Extract[fp3,{3}][[2]]};
fixedpoint3=Show[Graphics3D[{PointSize[0.025],Hue[0.33],Point[{(Subscript[x, 3]),(Subscript[y, 3]),(Subscript[z, 3])}]},PlotRange->{{(Subscript[x, min]),(Subscript[x, max])},{(Subscript[y, min]),(Subscript[y, max])},{(Subscript[z, min]),(Subscript[z, max])}}]];

Jacobian

a11=Simplify[D[f[x,y,z],x]]

a12=Simplify[D[f[x,y,z],y]]

a13=Simplify[D[f[x,y,z],z]]

a21=Simplify[D[g[x,y,z],x]]

a22=Simplify[D[g[x,y,z],y]]

a23=Simplify[D[g[x,y,z],z]]

a31=Simplify[D[h[x,y,z],x]]

a32=Simplify[D[h[x,y,z],y]]

a33=Simplify[D[h[x,y,z],z]]

J={{a11,a12,a13},{a21,a22,a23},{a31,a32,a33}}

MatrixForm[J]

Accélération

γ=J . V

Sur-Accélération

Γ=FullSimplify[(J . γ )/.{Dt[x]->f[x,y,z],Dt[y]->g[x,y,z],Dt[z]->h[x,y,z],Dt[(Subscript[a, 200])]->0,Dt[(Subscript[a, 110])]->0,Dt[(Subscript[a, 101])]->0,Dt[(Subscript[a, 011])]->0, Dt[(Subscript[a, 020])]->0,Dt[ (Subscript[a, 002])]->0,Dt[(Subscript[b, 200])]->0,Dt[(Subscript[b, 110])]->0,Dt[(Subscript[b, 101])]->0,Dt[(Subscript[b, 011])]->0, Dt[(Subscript[b, 020])]->0,Dt[ (Subscript[b, 002])]->0,Dt[(Subscript[c, 200])]->0,Dt[(Subscript[c, 110])]->0,Dt[(Subscript[c, 101])]->0,Dt[(Subscript[c, 011])]->0, Dt[(Subscript[c, 020])]->0,Dt[ (Subscript[c, 002])]->0,Dt[ (Subscript[c, 001])]->0,Dt[ (Subscript[a, 100])]->0,Dt[ (Subscript[b, 100])]->0,Dt[ (Subscript[a, 010])]->0,Dt[ (Subscript[b, 010])]->0}]

Fixed Points Stability - Eigenvalues

Eigenvalues[J]/.{x-> (Subscript[x, 1]),y->(Subscript[y, 1]),z->(Subscript[z, 1])}

Flow curvature Manifold 2

(Subscript[ϕ, 2])=Factor[Γ  . Cross[V,γ]]

(Subscript[ϕ, 0xy])=Collect[(Subscript[ϕ, 2]),{x,y,z}]/.{x->0}

(Subscript[ϕ, x0z])=Collect[(Subscript[ϕ, 2]),{x,y,z}]/.{y->0}

(Subscript[ϕ, 00z])=Collect[(Subscript[ϕ, 2]),{x,y,z}]/.{z->0}

(Subscript[ϕ, 1])=-σ ((r^2) (x^6)-r (x^5) y+(x^6) (y^2)-2 r (x^6) z+(x^5) y z+(x^6) (z^2)+(r^2) (x^4) β-(r^3) (x^4) β-r (x^3) y β+(r^2) (x^3) y β-r (x^5) y β+2 (x^4) (y^2) β-r (x^4) (y^2) β-3 r (x^4) z β+2 (r^2) (x^4) z β+2 (x^3) y z β-r (x^3) y z β-(x^5) y z β+2 (x^4) (z^2) β-r (x^4) (z^2) β-r (x^3) y (β^2)+(r^2) (x^3) y (β^2)+(x^2) (y^2) (β^2)-r (x^2) (y^2) (β^2)-r (x^2) z (β^2)+(r^2) (x^2) z (β^2)+r (x^4) z (β^2)+x y z (β^2)-r x y z (β^2)-2 (x^3) y z (β^2)+r (x^3) y z (β^2)+(x^2) (z^2) (β^2)-r (x^2) (z^2) (β^2)+r (x^2) z (β^3)-(r^2) (x^2) z (β^3)-x y z (β^3)+r x y z (β^3)+r (x^3) y σ-2 (r^2) (x^3) y σ+(r^3) (x^3) y σ-(x^2) (y^2) σ+2 r (x^2) (y^2) σ-(r^2) (x^2) (y^2) σ+2 (x^4) (y^2) σ-3 (x^3) (y^3) σ+r (x^3) (y^3) σ-(x^3) y z σ+4 r (x^3) y z σ-3 (r^2) (x^3) y z σ-2 (x^2) (y^2) z σ+2 r (x^2) (y^2) z σ-(x^3) (y^3) z σ-2 (x^3) y (z^2) σ+3 r (x^3) y (z^2) σ-(x^2) (y^2) (z^2) σ-(x^3) y (z^3) σ+r (x^3) y β σ-(r^2) (x^3) y β σ-x (y^3) β σ+r x (y^3) β σ-(x^3) (y^3) β σ+(x^2) (y^4) β σ+r (x^2) z β σ-2 (r^2) (x^2) z β σ+(r^3) (x^2) z β σ-2 r (x^4) z β σ-x y z β σ+2 r x y z β σ-(r^2) x y z β σ+3 r (x^3) y z β σ-(x^2) (y^2) z β σ-r (x^2) (y^2) z β σ+x (y^3) z β σ-(x^2) (z^2) β σ+4 r (x^2) (z^2) β σ-3 (r^2) (x^2) (z^2) β σ+2 (x^4) (z^2) β σ-2 x y (z^2) β σ+2 r x y (z^2) β σ-2 (x^3) y (z^2) β σ+(x^2) (y^2) (z^2) β σ-2 (x^2) (z^3) β σ+3 r (x^2) (z^3) β σ-x y (z^3) β σ-(x^2) (z^4) β σ-(x^2) (y^2) (β^2) σ+r (x^2) (y^2) (β^2) σ+x (y^3) (β^2) σ-r x (y^3) (β^2) σ-r (x^2) z (β^2) σ+(r^2) (x^2) z (β^2) σ+(x^3) y z (β^2) σ+(y^2) z (β^2) σ-r (y^2) z (β^2) σ-(x^2) (y^2) z (β^2) σ+r (x^2) (z^2) (β^2) σ-2 r x y (z^2) (β^2) σ+(y^2) (z^2) (β^2) σ-2 (x^2) (z^3) (β^2) σ+2 x y (z^3) (β^2) σ+x y z (β^3) σ-r x y z (β^3) σ-(y^2) z (β^3) σ+r (y^2) z (β^3) σ-(x^2) (z^2) (β^3) σ+2 x y (z^2) (β^3) σ-(y^2) (z^2) (β^3) σ+(x^2) (y^2) (σ^2)-2 r (x^2) (y^2) (σ^2)+(r^2) (x^2) (y^2) (σ^2)-x (y^3) (σ^2)+2 r x (y^3) (σ^2)-(r^2) x (y^3) (σ^2)+(x^2) (y^4) (σ^2)-x (y^5) (σ^2)+2 (x^2) (y^2) z (σ^2)-2 r (x^2) (y^2) z (σ^2)-2 x (y^3) z (σ^2)+2 r x (y^3) z (σ^2)+(x^2) (y^2) (z^2) (σ^2)-x (y^3) (z^2) (σ^2)+x (y^3) β (σ^2)-r x (y^3) β (σ^2)-(y^4) β (σ^2)+r (y^4) β (σ^2)+x y z β (σ^2)-2 r x y z β (σ^2)+(r^2) x y z β (σ^2)-(y^2) z β (σ^2)+2 r (y^2) z β (σ^2)-(r^2) (y^2) z β (σ^2)-2 (x^2) (y^2) z β (σ^2)+3 x (y^3) z β (σ^2)-(y^4) z β (σ^2)+2 x y (z^2) β (σ^2)-2 r x y (z^2) β (σ^2)-2 (y^2) (z^2) β (σ^2)+2 r (y^2) (z^2) β (σ^2)+x y (z^3) β (σ^2)-(y^2) (z^3) β (σ^2)-x y z (β^2) (σ^2)+r x y z (β^2) (σ^2)+(y^2) z (β^2) (σ^2)-r (y^2) z (β^2) (σ^2)+(x^2) (z^2) (β^2) (σ^2)-2 x y (z^2) (β^2) (σ^2)+(y^2) (z^2) (β^2) (σ^2))

Lets' use the procedure which consists in identifying (Subscript[ϕ, 1])(0,y, z) with (Subscript[ϕ, 2])(0,y,z) and (Subscript[ϕ, 1])(x,0, z) with (Subscript[ϕ, 2])(x,0,z) ; (Subscript[ϕ, 1])(x,y, 0) with (Subscript[ϕ, 2])(x,y,0) and (Subscript[ϕ, 1])(x,0, 0) with (Subscript[ϕ, 2])(x,0,0)  and (Subscript[ϕ, 1])(0,y, 0) with (Subscript[ϕ, 2])(0,y,0)

Inverse Problem Procedure

(Subscript[ϕ, 0xy])=Collect[(Subscript[ϕ, 1])-(Subscript[ϕ, 2]),{x,y,z}]/.{x->0}

(Subscript[ϕ, x0z])=Collect[(Subscript[ϕ, 1])-(Subscript[ϕ, 2]),{x,y,z}]/.{y->0}

(Subscript[ϕ, 00z])=Collect[(Subscript[ϕ, 1])-(Subscript[ϕ, 2]),{x,y,z}]/.{z->0}

(Subscript[ϕ, x00])=Collect[(Subscript[ϕ, 1])-(Subscript[ϕ, 2]),{x,y,z}]/.{y->0,z->0}

(Subscript[ϕ, 0x0])=Collect[(Subscript[ϕ, 1])-(Subscript[ϕ, 2]),{x,y,z}]/.{x->0,z->0}

Solve[{ (β (σ^3)+(Subscript[a, 10]^3) (Subscript[b, 101]^2) (Subscript[c, 1]))==0,(β (σ^3)-(Subscript[a, 10]^3) (Subscript[b, 101]) (Subscript[c, 1]) (Subscript[c, 110]))==0,(β (σ^2)-(Subscript[a, 10]^2) (Subscript[b, 101]^3) (Subscript[c, 1]))==0,(-σ+(Subscript[a, 100]) (Subscript[b, 101]^3) (Subscript[c, 110]^2))==0,(2 r σ+2 (Subscript[a, 100]) (Subscript[b, 100]) (Subscript[b, 101]^2) (Subscript[c, 110]^2))==0,(-(r^2) σ+(Subscript[a, 100]) (Subscript[b, 100]^2) (Subscript[b, 101]) (Subscript[c, 110]^2))==0,(-(r^2) σ+(Subscript[a, 100]) (Subscript[b, 100]^2) (Subscript[b, 101]) (Subscript[c, 110]^2))==0, (-σ-(Subscript[a, 100]) (Subscript[b, 101]^2) (Subscript[c, 110]^3))==0,(β (σ^3)-r β (σ^3)+(Subscript[a, 10]^2) (Subscript[a, 100]) (Subscript[b, 10]) (Subscript[c, 1]) (Subscript[c, 110])-(Subscript[a, 10]^3) (Subscript[b, 100]) (Subscript[c, 1]) (Subscript[c, 110]))==0,(-(r^2) β σ+(r^3) β σ-(Subscript[a, 100]) (Subscript[b, 10]) (Subscript[b, 100]^2) (Subscript[c, 1]) (Subscript[c, 110])+(Subscript[a, 10]) (Subscript[b, 100]^3) (Subscript[c, 1]) (Subscript[c, 110]))==0, (-(r^2) σ+(Subscript[a, 100]) (Subscript[b, 100]^2) (Subscript[b, 101]) (Subscript[c, 110]^2))==0},{(Subscript[a, 100]),(Subscript[a, 010]),(Subscript[b, 100]),(Subscript[b, 010]),(Subscript[b, 101]),(Subscript[c, 110]),(Subscript[c, 001]) }]

FullSimplify[{(Subscript[b, 10])->-(Subscript[c, 110]^7),(Subscript[a, 100])->-(σ/(Subscript[c, 110]^5)),(Subscript[c, 1])->-(β/(Subscript[c, 110]^5)),(Subscript[b, 100])->r (Subscript[c, 110]),(Subscript[a, 10])->σ (Subscript[c, 110]),(Subscript[b, 101])->-(Subscript[c, 110])}/.{(Subscript[c, 110])->1}]

Factor[Expand[V/.{(Subscript[b, 10])->-1,(Subscript[a, 100])->-σ,(Subscript[c, 1])->-β,(Subscript[b, 100])->r,(Subscript[a, 10])->σ,(Subscript[b, 101])->-1,(Subscript[c, 110])->1}]]

Thus the Lorenz vector field is found again.