restart: with (linalg):with(liesymm):with(difforms):with(plots):

topdynamicsexamples.mws

Examples of a topological thermodynamic System

See Vol. 1, "Non-Equilibrium Thermodynamics"

http://www.lulu.com/kiehn

chapter 3.8.

R. M. Kiehn

updated Sept/05/2005 --Feb 6 2008

A Maple program that uses the Jacobian Matrix (dyadic)

of a 3D velocity field of a dynamical system to compute the

similarity invariants, Xm, Yg, Za, and the Brand invariants

setup(x,y,z,t);defform(x=0,y=0,z=0,u=0,f=0,t=0,Vx=0,Vy=0,Vz=0,D1=0,D2=0,D3=0,Ax=0,Ay=0,Az=0,Phi=0,a=const,b=const,c=const,q=const,Lx=0,Ly=0,Lz=0,A=const,B=const,C=const,D=const,Omega=const,r=const);

STF dynamical system

AA:=[C*z*x-Omega*y,(C-R)*z*y+Omega*x,A-B*z^2-S*x^2-T*y^2];

phi:=simplify((AA[1]^2+AA[2]^2+AA[3]^2));VV:=evalm([AA[1],AA[2],AA[3],phi]);

JAC:=simplify(jacobian(VV,[x,y,z,t])):collect(simplify(expand(charpoly(JAC,q))),q);

S1:=factor(trace(JAC)):S2:=factor(trace(innerprod(JAC,JAC))):S3:=factor(trace(innerprod(JAC,JAC,JAC))):

Xm:=S1;Yg:=factor((1/2)*(Xm*S1-S2));Za:=collect(collect(factor((1/3)*(Yg*S1-Xm*S2+S3)),z),Omega);Tk:=factor(det(JAC));

;BBB:=solve(Xm,B);YgB:=factor(subs(B=BBB,Yg));RRR:=solve(YgB,Omega^2);ZaBR:=factor(subs(B=BBB,Omega^2=RRR,Za));CriticalPoint:=collect(ZaBR-1,Omega);

Thermodynamic crtical point is a surface in {x,y,z} (It is not locally stable)

BBB:=solve(Xm,B);ZaB:=factor(subs(B=BBB,Za));RRR:=solve(ZaB,Omega^2);YYg:=factor(subs(B=BBB,Omega^2=RRR,Yg)); Hopfcrit:=collect(YYg,Omega);

TT:=solve(factor(2*S*x*y*R+2*T*y*C*x-2*S*x*y*C),T);ZZ:=subs(T=TT,collect(collect(collect((4*S*x^2*C*z*R+8*C*z*T*y^2*R-6*C^2*z*T*y^2-6*S*x^2*C^2*z-2*T*y^2*z*R^2+8*C^3*z^3-z^3*R^3-12*C^2*z^3*R+6*C*z^3*R^2)/z,z),x),y));

Chiral - Brand invariants

VORTICITY:=evalm(curl(AA,[x,y,z])):CURLx:=simplify(VORTICITY[1]);CURLy:=simplify(VORTICITY[2]);CURLz:=simplify(VORTICITY[3]);

HELICITY:=collect(factor(innerprod(AA,VORTICITY)),Omega);

KKK:=simplify((2*S*x*y*R+2*T*y*C*x-2*S*x*y*C)/(x*y));

if KKK is zero then helicity depends upon Omega.

enstrophy:=collect(factor(innerprod(VORTICITY,VORTICITY)),Omega);

BRAND:=innerprod(jacobian(AA,[x,y,z]),VORTICITY);

stretch:=collect(factor(innerprod(BRAND,VORTICITY)),Omega);ST1:=factor(coeff(stretch,Omega^2));st0:=factor(subs(Omega=0,stretch));

BI:=collect(factor(innerprod(BRAND,BRAND)),Omega);BI2:=factor(coeff(BI,Omega^2));BI1:=factor(coeff(BI,Omega));BI0:=collect(factor(subs(Omega=0,BI)),z);

Theta:=factor(-2*innerprod(grad(phi,[x,y,z]),VORTICITY));

Theta:=collect(factor(-2*innerprod(grad(phi,[x,y,z]),VORTICITY)),Omega);t2:=factor(coeff(Theta,Omega^2));t1:=factor(coeff(Theta,Omega));

t0:=factor(subs(Omega=0,Theta));

(collect(innerprod(AA,grad(HELICITY,[x,y,z])),Omega));SPHERE:=collect(collect(factor(innerprod(AA,[x,y,z])/z),x^2),y^2);

C-R-T=-A;C-S=-A;B=A;S=C+A;T=C+A-R;T=S-R;

Lorenz dynamical system

AA:=[A*(y-x),-x*z+R*x-y,x*y-B*z];

phi:=simplify((AA[1]^2+AA[2]^2+AA[3]^2));VV:=[AA[1],AA[2],AA[3],phi];

JAC:=simplify(jacobian(evalm(VV),[x,y,z,t]));collect(simplify(expand(charpoly(JAC,q))),q);

S1:=factor(trace(JAC)):S2:=factor(trace(innerprod(JAC,JAC))):S3:=factor(trace(innerprod(JAC,JAC,JAC))):

Xm:=S1;Yg:=factor((1/2)*(Xm*S1-S2));Za:=factor((1/3)*(Yg*S1-Xm*S2+S3));Tk:=factor(det(JAC));

BBB:=solve(Xm-3,B);YgB:=factor(subs(B=BBB,Yg));RRR:=solve(YgB-3,R);ZaBR:=factor(subs(B=BBB,R=RRR,Za));CriticalPoint:=collect(ZaBR-1,A);

Thermodynamic crtical point is a surface in {x,y,z} (It is not locally stable)

BBB:=solve(Xm,B);ZaB:=factor(subs(B=BBB,Za));RRR:=solve(ZaB,R);YYg:=factor(subs(B=BBB,R=RRR,Yg)); Hopfcrit:=collect(YYg,A);

INTER:=collect((numer(CriticalPoint)+numer(Hopfcrit)),A);solve(INTER,A);;

Hopf conditions

AAA:=solve(Xm,A);ZaAAA:=subs(A=AAA,Za);BBA:=solve(ZaAAA,B);YYg:=factor(subs(A=AAA,B=BBA[1],Yg));YYg2:=factor(subs(A=AAA,B=BBA[2],Yg));;

When YYg>0, then get oscillations. Can solve for the zero set of YYg which will divide the space into regions of oscillations and no oscillations. The parameter of the surface is R, can ask for discriminant to see if there is an envelope.

YYgN:=numer(YYg);factor(discrim(YYgN,R));

Chiral - Brand invariants

VORTICITY:=evalm(curl(AA,[x,y,z])):CURLx:=simplify(VORTICITY[1]);CURLy:=simplify(VORTICITY[2]);CURLz:=simplify(VORTICITY[3]);

HELICITY:=collect(factor(innerprod(AA,VORTICITY)),B);

enstrophy:=collect(factor(innerprod(VORTICITY,VORTICITY)),R);

BRAND:=innerprod(jacobian(AA,[x,y,z]),VORTICITY);

stretch:=collect(factor(innerprod(BRAND,VORTICITY)),B);ST1:=factor(coeff(stretch,B));

BI:=collect(factor(innerprod(BRAND,BRAND)),B);BI2:=factor(coeff(BI,B^2));BI1:=factor(coeff(BI,B));BI0:=collect(factor(subs(B=0,BI)),R);

Theta:=collect(factor(-2*innerprod(grad(phi,[x,y,z]),VORTICITY)),R);;Theta:=collect(factor(-2*innerprod(grad(phi,[x,y,z]),VORTICITY)),B);;Theta:=collect(factor(-2*innerprod(grad(phi,[x,y,z]),VORTICITY)),A);;Theta:=collect(factor(-2*innerprod(grad(phi,[x,y,z]),VORTICITY)),z);;factor(-8*A^2*x^2+4*A*x^2*z+4*y^2+4*y*x*A+4*A^2*y^2-4*x^2*z^2+4*A^2*x*y+4*B^2*z^2-4*x^2*y^2-4*x*y*A*B+4*B^2*z*A)/4;

HENON dynamical system

AA:=[A*cos(y)+B*sin(z),B*cos(z)+C*sin(x),C*cos(x)+A*sin(y)];

phi:=simplify((AA[1]^2+AA[2]^2+AA[3]^2));VV:=[AA[1],AA[2],AA[3],phi];

JAC:=simplify(jacobian(evalm(VV),[x,y,z,t]),trig);

S1:=factor(trace(JAC)):S2:=factor(trace(innerprod(JAC,JAC))):S3:=factor(trace(innerprod(JAC,JAC,JAC))):

Xm:=S1;Yg:=factor((1/2)*(Xm*S1-S2));Za:=factor((1/3)*(Yg*S1-Xm*S2+S3));Tk:=factor(det(JAC));

Thermodynamic crtical point is a surface in {x,y,z} (It is not locally stable)

YYg:=factor(subs(A=0,Yg));

Chiral - Brand invariants

VORTICITY:=evalm(curl(AA,[x,y,z])):CURLx:=simplify(VORTICITY[1]);CURLy:=simplify(VORTICITY[2]);CURLz:=simplify(VORTICITY[3]);

HELICITY:=collect(simplify(innerprod(AA,VORTICITY),trig),B);

enstrophy:=collect(simplify(innerprod(VORTICITY,VORTICITY),trig),R);

BRAND:=innerprod(jacobian(AA,[x,y,z]),VORTICITY);

stretch:=simplify(innerprod(BRAND,VORTICITY),trig);

BI:=(simplify(innerprod(BRAND,BRAND),trig));

Theta:=factor(-2*innerprod(grad(phi,[x,y,z]),VORTICITY));

Taylor-Green dynamical system

AA:=[A*sin(x)*cos(y)*cos(z),B*cos(x)*sin(y)*cos(z),C*cos(x)*cos(y)*sin(z)];

phi:=simplify((AA[1]^2+AA[2]^2+AA[3]^2));VV:=[AA[1],AA[2],AA[3],phi];

JAC:=simplify(jacobian(evalm(VV),[x,y,z,t]));

S1:=factor(trace(JAC)):S2:=simplify(trace(innerprod(JAC,JAC)),trig):S3:=simplify(trace(innerprod(JAC,JAC,JAC)),trig):

Xm:=S1;Yg:=simplify((1/2)*(Xm*S1-S2),trig);Za:=factor((1/3)*(Yg*S1-Xm*S2+S3));Tk:=factor(det(JAC));

Chiral - Brand invariants

VORTICITY:=evalm(curl(AA,[x,y,z])):CURLx:=simplify(VORTICITY[1]);CURLy:=simplify(VORTICITY[2]);CURLz:=simplify(VORTICITY[3]);

HELICITY:=collect(factor(innerprod(AA,VORTICITY)),B);

enstrophy:=subs(cos(x)=X,cos(y)=Y,cos(z)=Z,(simplify(innerprod(VORTICITY,VORTICITY),trig)));E1:=collect(enstrophy,A);

BRAND:=innerprod(jacobian(AA,[x,y,z]),VORTICITY);

stretch:=collect(factor(innerprod(BRAND,VORTICITY)),B);ST1:=factor(coeff(stretch,B));

BI:=collect(factor(innerprod(BRAND,BRAND)),B);BI2:=factor(coeff(BI,B^2));BI1:=factor(coeff(BI,B));BI0:=collect(factor(subs(B=0,BI)),R);

Theta:=factor(-2*innerprod(grad(phi,[x,y,z]),VORTICITY));

Moffat dynamical system

AA:=[A*x+Omega*z,B*y+G*x^2+F*x*z,C*z-Omega*x-F*x*y];

phi:=simplify((AA[1]^2+AA[2]^2+AA[3]^2));VV:=evalm([AA[1],AA[2],AA[3],phi]);

JAC:=simplify(jacobian(evalm(VV),[x,y,z,t]));collect(simplify(expand(charpoly(JAC,q))),q);

S1:=factor(trace(JAC)):S2:=factor(trace(innerprod(JAC,JAC))):S3:=factor(trace(innerprod(JAC,JAC,JAC))):

Xm:=S1;Yg:=factor((1/2)*(Xm*S1-S2));Za:=collect(factor((1/3)*(Yg*S1-Xm*S2+S3)),Omega);;Tk:=factor(det(JAC));

BBB:=solve(Xm-3,B);YgB:=factor(subs(B=BBB,Yg));FFF:=solve(YgB-3,F);ZaBR:=factor(subs(B=BBB,F=FFF,Za));CriticalPoint:=collect(ZaBR-1,Omega);

BBB:=solve(Xm,B);ZaB:=factor(subs(B=BBB,Za));GGG:=solve(ZaB,G);YYg:=factor(subs(B=BBB,G=GGG,Yg)); Hopfcrit:=collect(YYg,A);

INTER:=collect((numer(CriticalPoint)+numer(Hopfcrit)),A);solve(INTER,A);;

Chiral - Brand invariants

VORTICITY:=evalm(curl(AA,[x,y,z])):CURLx:=simplify(VORTICITY[1]);CURLy:=simplify(VORTICITY[2]);CURLz:=simplify(VORTICITY[3]);

HELICITY:=collect(factor(innerprod(AA,VORTICITY)),Omega);

enstrophy:=collect(factor(innerprod(VORTICITY,VORTICITY)),Omega);

BRAND:=innerprod(jacobian(AA,[x,y,z]),VORTICITY);

stretch:=collect(factor(innerprod(BRAND,VORTICITY)),Omega);ST1:=factor(coeff(stretch,Omega));st0:=collect(factor(subs(Omega=0,stretch)),F);

BI:=collect(factor(innerprod(BRAND,BRAND)),Omega);BI2:=factor(coeff(BI,Omega^2));BI1:=factor(coeff(BI,Omega));BI0:=collect(factor(subs(Omega=0,BI)),F);

Theta:=collect(factor(-2*innerprod(grad(phi,[x,y,z]),VORTICITY)),Omega);t2:=factor(coeff(Theta,Omega^2));t1:=collect(factor(-coeff(Theta,Omega)),F);;

t0:=collect(factor(subs(Omega=0,Theta)),F);;

Brusselator dynamical system

AA:=[A-Omega*(x-G*x^2*y)-D*x,Omega*(x-G*x^2*y),0];

phi:=simplify((AA[1]^2+AA[2]^2+AA[3]^2));VV:=evalm([AA[1],AA[2],AA[3],phi]);

JAC:=simplify(jacobian(VV,[x,y,z,t]));collect(simplify(expand(charpoly(JAC,q))),q);

S1:=factor(trace(JAC)):S2:=factor(trace(innerprod(JAC,JAC))):S3:=factor(trace(innerprod(JAC,JAC,JAC))):

Xm:=collect(S1,Omega);;Yg:=factor((1/2)*(Xm*S1-S2));Za:=factor((1/3)*(Yg*S1-Xm*S2+S3));Tk:=factor(det(JAC));

Loss of Local Stability

DDD:=solve(Xm,D);YgD:=subs(D=DDD,Yg);solve(YgD,Omega);

thermal Critical Point

DDD:=solve(Xm-1,D);YgD:=subs(D=DDD,Yg-1);solve(YgD,Omega);

Hopf conditions

DDD:=solve(Xm,D);YYg:=subs(D=DDD,Yg);EV:=eigenvalues(JAC):EV1:=factor(subs(D=DDD,EV[3]));EV2:=factor(subs(D=DDD,EV[4]));EEVV:=(YYg)^(1/2);

When YYg>0 get oscillations.

This implies that LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYlLUkjbWlHRiQ2I1EhRictRiM2Jy1JKG1mZW5jZWRHRiQ2JC1GIzYoLUkjbW5HRiQ2JFEiMUYnLyUsbWF0aHZhcmlhbnRHUSdub3JtYWxGJy1JI21vR0YkNi1RKCZtaW51cztGJ0Y6LyUmZmVuY2VHUSZmYWxzZUYnLyUqc2VwYXJhdG9yR0ZDLyUpc3RyZXRjaHlHRkMvJSpzeW1tZXRyaWNHRkMvJShsYXJnZW9wR0ZDLyUubW92YWJsZWxpbWl0c0dGQy8lJ2FjY2VudEdGQy8lJ2xzcGFjZUdRLDAuMjIyMjIyMmVtRicvJSdyc3BhY2VHRlItRiM2KS1GNzYkUSIyRidGOi1GPjYtUTEmSW52aXNpYmxlVGltZXM7RidGOkZBRkRGRkZIRkpGTEZOL0ZRUSYwLjBlbUYnL0ZURmhuLUYsNiVRIkdGJy8lJ2l0YWxpY0dRJXRydWVGJy9GO1EnaXRhbGljRidGWi1GLDYlUSJ4RidGXW9GYG9GWi1GLDYlUSJ5RidGXW9GYG8tRj42LVEiK0YnRjpGQUZERkZGSEZKRkxGTkZQRlMtRiM2JkZqbkZaLUYjNiMtSSVtc3VwR0YkNiVGYm9GVy8lMXN1cGVyc2NyaXB0c2hpZnRHUSIwRidGK0YrRjpGWi1GIzYjLUZgcDYlLUYsNiVRKCZPbWVnYTtGJy9GXm9GQ0Y6RldGYnBGWkZqbkYr < 0

Chiral - Brand invariants

VORTICITY:=evalm(curl(AA,[x,y,z]));CURLx:=simplify(VORTICITY[1]);CURLy:=simplify(VORTICITY[2]);CURLz:=simplify(VORTICITY[3]);

HELICITY:=(factor(innerprod(AA,VORTICITY)));

enstrophy:=(factor(innerprod(VORTICITY,VORTICITY)));

BRAND:=innerprod(jacobian(AA,[x,y,z]),VORTICITY);

stretch:=factor(innerprod(BRAND,VORTICITY));

BI:=(factor(innerprod(BRAND,BRAND)));

Theta:=factor(-2*innerprod(grad(phi,[x,y,z]),VORTICITY));

E:=grad(phi,[x,y,z]);AA;TC:=evalm((crossprod(E,AA)+evalm((phi*VORTICITY))));

TC3:=collect(factor(TC[3]),Omega);;

a1:=factor(coeff(TC3,Omega^3));a2:=factor(coeff(TC3,Omega^2));a3:=factor(coeff(TC3,Omega));

Chiral Separation (KONDEPUDI)

AA:=[(B)*x+2*C*x*y,+D+(B)*y-C*(y^2+x^2)-E*(y^2-x^2),0];

phi:=((AA[1]^2+AA[2]^2+AA[3]^2));VV:=evalm([AA[1],AA[2],AA[3],phi]);

JAC:=(jacobian(VV,[x,y,z,t]));

S1:=factor(trace(JAC)):S2:=factor(trace(innerprod(JAC,JAC))):S3:=factor(trace(innerprod(JAC,JAC,JAC))):

Xm:=collect(S1,Omega);;Yg:=factor((1/2)*(Xm*S1-S2));Za:=factor((1/3)*(Yg*S1-Xm*S2+S3));Tk:=factor(det(JAC));

HOpf conditions

BBB:=solve(Xm,B);

EV:=eigenvalues(JAC);EV3:=subs(B=BBB,EV[3]);EV4:=subs(B=BBB,EV[4]);

YYg:=factor(subs(B=BBB,Yg));CCC:=solve(YYg,C);CCC1:=factor(CCC[1]);CCC2:=factor(CCC[2]);

XX0:=solve(AA[1],x);YY0:=solve(subs(x=XX,AA[2]),y);YY00:=solve(subs(AA[1]),y);XX00:=solve(subs(y=YY00,AA[2]),x);

Fixed Points and Hopf.

BBB:=solve(Xm-1,B);YgB:=factor(subs(B=BBB,Yg));EEE:=solve(YgB-1,E);

There are two thermo critcal points given by the two values of E

Thermodynamic crtical point is a surface in {x,y,z} (It is not locally stable)

BBB:=solve(Xm,B);YYg:=factor(subs(B=BBB,Yg)); Hopfcrit:=simplify(YYg);

EEE[1];

INTER:=factor(subs(E=EEE[2],Hopfcrit));

YYYG1:=factor(subs(x=XX0,y=YY0[1],YYg));YYYG2:=factor(subs(x=XX0,y=YY0[2],YYg));YYYYG1:=factor(subs(x=XX00[1],y=YY00,(YYg)));YYYYG2:=simplify(subs(x=XX00[2],y=YY00,(YYg)));

If B goes to zero, then Yg is always negative implying that it is a minimal surface of the negative gauss curvature type (soap film), not a Hopf bifurcation for the Racemic alpha = 0 fixed point.

BUT if B goes to zero then it is true that Yg is positive depending on CD....Yg=4CD

Chiral - Brand invariants

VORTICITY:=evalm(curl(AA,[x,y,z]));CURLx:=(VORTICITY[1]);CURLy:=(VORTICITY[2]);CURLz:=(VORTICITY[3]);

HELICITY:=(factor(innerprod(AA,VORTICITY)));

enstrophy:=(factor(innerprod(VORTICITY,VORTICITY)));

BRAND:=innerprod(jacobian(AA,[x,y,z]),VORTICITY);

stretch:=factor(innerprod(BRAND,VORTICITY));

BI:=(factor(innerprod(BRAND,BRAND)));

Theta:=factor(-2*innerprod(grad(phi,[x,y,z]),VORTICITY));

EE:=grad(phi,[x,y,z]);AA:=AA;Bphi:=[phi*VORTICITY[1],Phi*VORTICITY[2],phi*VORTICITY[3]];TC:=evalm((crossprod(EE,AA)+Bphi));CP:=crossprod(EE,AA);EA3:=factor(CP[3]);

TC3:=collect(factor(TC[3]),B);;

a2:=factor(coeff(TC3,B^2));a1:=factor(coeff(TC3,B));a0:=collect(factor(subs(B=0,TC3)),C);;

factor(subs(B=BBB,C=CCC1,TC3)):

Dynamo dynamical system

AA:=[-B*x+C*z*y+Omega*y,-B*y+C*z*x-Omega*x,1-x*y+D*z^2];

phi:=simplify((AA[1]^2+AA[2]^2+AA[3]^2));VV:=evalm([AA[1],AA[2],AA[3],phi]);

JAC:=simplify(jacobian(VV,[x,y,z,t]));collect(simplify(expand(charpoly(JAC,q))),q);

S1:=factor(trace(JAC)):S2:=factor(trace(innerprod(JAC,JAC))):S3:=factor(trace(innerprod(JAC,JAC,JAC))):

Xm:=S1;Yg:=factor((1/2)*(Xm*S1-S2));Za:=factor((1/3)*(Yg*S1-Xm*S2+S3));Tk:=factor(det(JAC));

BBB:=solve(Xm-3,B);YgB:=factor(subs(B=BBB,Yg));RRR:=solve(YgB-3,Omega^2);ZaBR:=factor(subs(B=BBB,Omega^2=RRR,Za));CriticalPoint:=collect(ZaBR-1,C);

Thermodynamic crtical point is a surface in {x,y,z} (It is not locally stable)

BBB:=solve(Xm,B);ZaB:=factor(subs(B=BBB,Za));RRR:=solve(ZaB,Omega^2);YYg:=factor(subs(B=BBB,Omega^2=RRR,Yg)); Hopfcrit:=collect(YYg,A);

Chiral - Brand invariants

VORTICITY:=evalm(curl(AA,[x,y,z])):CURLx:=simplify(VORTICITY[1]);CURLy:=simplify(VORTICITY[2]);CURLz:=simplify(VORTICITY[3]);

HELICITY:=collect(factor(innerprod(AA,VORTICITY)),Omega);

enstrophy:=collect(factor(innerprod(VORTICITY,VORTICITY)),Omega);

BRAND:=innerprod(jacobian(AA,[x,y,z]),VORTICITY);

stretch:=collect(factor(innerprod(BRAND,VORTICITY)),Omega);ST1:=factor(coeff(stretch,Omega^2));st0:=factor(subs(Omega=0,stretch));

BI:=collect(factor(innerprod(BRAND,BRAND)),Omega);BI2:=factor(coeff(BI,Omega^2));BI1:=factor(coeff(BI,Omega));BI0:=collect(factor(subs(Omega=0,BI)),z);

Theta:=factor(-2*innerprod(grad(phi,[x,y,z]),VORTICITY));

Theta:=collect(factor(-2*innerprod(grad(phi,[x,y,z]),VORTICITY)),Omega);t2:=factor(coeff(Theta,Omega^2));t1:=factor(coeff(Theta,Omega));

t0:=factor(subs(Omega=0,Theta));

SPINNING TOP-Eulerian Rotator dynamical system

AA:=[A*y*z,B*z*x,C*x*y];

phi:=simplify((AA[1]^2+AA[2]^2+AA[3]^2));VV:=evalm([AA[1],AA[2],AA[3],phi]);

JAC:=simplify(jacobian(VV,[x,y,z,t]));collect(simplify(expand(charpoly(JAC,q))),q);

S1:=factor(trace(JAC)):S2:=factor(trace(innerprod(JAC,JAC))):S3:=factor(trace(innerprod(JAC,JAC,JAC))):

Xm:=S1;Yg:=factor((1/2)*(Xm*S1-S2));Za:=factor((1/3)*(Yg*S1-Xm*S2+S3));Tk:=factor(det(JAC));

Hopfcrit:=factor(subs(z=0,Yg));

Chiral - Brand invariants

VORTICITY:=evalm(curl(AA,[x,y,z])):CURLx:=simplify(VORTICITY[1]);CURLy:=simplify(VORTICITY[2]);CURLz:=simplify(VORTICITY[3]);

HELICITY:=(factor(innerprod(AA,VORTICITY)));

enstrophy:=factor(innerprod(VORTICITY,VORTICITY));

BRAND:=innerprod(jacobian(AA,[x,y,z]),VORTICITY);

stretch:=(factor(innerprod(BRAND,VORTICITY)/(1)));factor(C^2*A+B^2*A+C^2*B+A^2*B+B^2*C+A^2*C);

BI:=simplify(factor(innerprod(BRAND,BRAND)));zz:=(C^2*y^2*A^2*z^2-2*C*y^2*A^2*z^2*B+B^2*z^2*A^2*y^2)*Omega^2+(C^2*x^2*B^2*z^2-2*C*x^2*B^2*z^2*A+A^2*z^2*B^2*x^2)*Omega+B^2*x^2*C^2*y^2-2*B*x^2*C^2*y^2*A+A^2*y^2*C^2*x^2;

factor(coeff(zz,Omega^2));factor(coeff(zz,Omega));

Theta:=factor(-2*innerprod(grad(phi,[x,y,z]),VORTICITY));

Theta:=collect(factor(-2*innerprod(grad(phi,[x,y,z]),VORTICITY)),Omega);t2:=factor(coeff(Theta,Omega^2));t1:=factor(coeff(Theta,Omega));

t0:=factor(subs(Omega=0,Theta));SYMM:=factor(subs(C=0,B=-A,Theta));factor(subs(I1=1,I2=1,I3=1/2,SYMM));

ROSSLER dynamical system

AA:= [W*y-z,W*x-A*y,B+x*z-C*y+T*z];

phi:=simplify((AA[1]^2+AA[2]^2+AA[3]^2));VV:=evalm([AA[1],AA[2],AA[3],phi]);

JAC:=(jacobian(VV,[x,y,z,t]));

S1:=factor(trace(JAC)):S2:=factor(trace(innerprod(JAC,JAC))):S3:=factor(trace(innerprod(JAC,JAC,JAC))):

Xm:=S1;Yg:=factor((1/2)*(Xm*S1-S2));Za:=factor((1/3)*(Yg*S1-Xm*S2+S3));Tk:=factor(det(JAC));

TTT:=solve(Xm-3,T);ZaB:=factor(subs(T=TTT,Za));AAA:=solve(ZaB-3,A);YgBR:=factor(subs(T=TTT,A=AAA,Yg));CriticalPoint:=collect(factor(expand(YgBR-1)),W);

Thermodynamic crtical point is a surface in {x,y,z} (It is not locally stable)

TTT:=solve(Xm,T);ZaB:=factor(subs(T=TTT,Za));CCC:=solve(ZaB,C);YYg:=factor(subs(T=TTT,C=CCC,Yg)); Hopfcrit:=collect(YYg,A);

Chiral - Brand invariants

VORTICITY:=evalm(curl(AA,[x,y,z])):CURLx:=simplify(VORTICITY[1]);CURLy:=simplify(VORTICITY[2]);CURLz:=simplify(VORTICITY[3]);

HELICITY:=collect(factor(innerprod(AA,VORTICITY)),W);

enstrophy:=(factor(innerprod(VORTICITY,VORTICITY)));

BRAND:=innerprod(jacobian(AA,[x,y,z]),VORTICITY);

stretch:=collect(factor(innerprod(BRAND,VORTICITY)),W);

BI:=collect(factor(innerprod(BRAND,BRAND)),W);BI2:=factor(coeff(BI,W^2));BI1:=factor(coeff(BI,W));BI0:=collect(factor(subs(W=0,BI)),z);

Theta:=collect(factor(-2*innerprod(grad(phi,[x,y,z]),VORTICITY)),W);Theta:=collect(factor(-2*innerprod(grad(phi,[x,y,z]),VORTICITY)),A);

Zhabotinsky dynamical system

AA:=[0,0,0];VV:=[0,0,0,0];phi:=0;

AA:=evalm([A*y+B*x-C*x*y-2*D*x*x,-A*y-C*x*y+F*E*z/2,2*B*x-E*z]);

phi:=factor((AA[1]^2+AA[2]^2+AA[3]^2));VV:=evalm([AA[1],AA[2],AA[3],phi]);

JAC:=(jacobian(VV,[x,y,z,t]));

S1:=factor(trace(JAC)):S2:=factor(trace(innerprod(JAC,JAC))):S3:=factor(trace(innerprod(JAC,JAC,JAC))):

Xm:=S1;Yg:=collect(factor((1/2)*(Xm*S1-S2)),E);;Za:=factor((1/3)*(Yg*S1-Xm*S2+S3));Tk:=factor(det(JAC));

Hopf conditions

AAA:=factor(solve(Xm,A));ZaEEE:=factor(subs(A=AAA,Za));FFF:=solve(ZaEEE,F);YgAF:=factor(subs(A=AAA,F=FFF,Yg));

Chiral - Brand invariants

VORTICITY:=(curl(AA,[x,y,z])):CURLx:=simplify(VORTICITY[1]);CURLy:=simplify(VORTICITY[2]);CURLz:=simplify(VORTICITY[3]);

HELICITY:=collect(factor(innerprod(AA,VORTICITY)),F);

enstrophy:=collect(factor(innerprod(VORTICITY,VORTICITY)),F);

BRAND:=innerprod(jacobian(AA,[x,y,z]),VORTICITY);

stretch:=collect(factor(innerprod(BRAND,VORTICITY)),F);ST2:=factor(coeff(stretch,F^2));ST1:=factor(coeff(stretch,F));st0:=(factor(subs(F=0,stretch)));st0e:=collect(st0,E);E0:=factor(coeff(st0e,E));

BI:=collect(factor(innerprod(BRAND,BRAND)),F);BI2:=factor(coeff(BI,F^2));BI1:=factor(coeff(BI,F));BI0:=collect(factor(subs(F=0,BI)),E);

Theta:=collect(factor(-2*innerprod(grad(phi,[x,y,z]),VORTICITY)),F);T1:=factor(coeff(Theta,F));

factor(subs(F=0,Theta));

Pitchfork Hopf OK

AA:= [e*A-G+C*z^2+T*(x^2+y^2)+Omega*y,e*A-G+C*z^2+T*(x^2+y^2)-Omega*x,z*(A+B*z^2+D*(x^2+y^2))];

phi:=factor((AA[1]^2+AA[2]^2+AA[3]^2));VV:=[AA[1],AA[2],AA[3],phi];

JAC:=simplify(jacobian(VV,[x,y,z,t]));collect(simplify(expand(charpoly(JAC,q))),q):

S1:=factor(trace(JAC)):S2:=factor(trace(innerprod(JAC,JAC))):S3:=factor(trace(innerprod(JAC,JAC,JAC))):

Xm:=S1;Yg:=collect(factor((1/2)*(Xm*S1-S2)),Omega);Za:=collect(factor((1/3)*(Yg*S1-Xm*S2+S3)),Omega);Tk:=factor(det(JAC));Yg0:=factor(subs(Omega=0,Yg));Za1:=factor(coeff(Za,Omega));

Hopf conditions

AAA:=factor(solve(Xm,A));ZaAAA:=subs(A=AAA,Za);DDD:=(solve(ZaAAA,D));YYg:=factor(subs(A=AAA,D=DDD,Yg));

Oscillations occur when YYg>0

Fixed points

RRR:=solve(ZaAAA,Omega);HopfCrit:=factor(subs(A=AAA,Omega=RRR[1],Yg));HopfCrit:=factor(subs(A=AAA,Omega=RRR[2],Yg));

Chiral - Brand invariants

VORTICITY:=evalm(curl(AA,[x,y,z])):CURLx:=simplify(VORTICITY[1]);CURLy:=simplify(VORTICITY[2]);CURLz:=simplify(VORTICITY[3]);EE1:=-factor(diff(phi,x));EE2:=-factor(diff(phi,y));EE3:=-factor(diff(phi,z));

HELICITY:=collect(factor(innerprod(AA,VORTICITY)),Omega);H1:=factor(coeff(HELICITY,Omega));H0:=factor(subs(Omega=0,HELICITY));

enstrophy:=collect(factor(innerprod(VORTICITY,VORTICITY)),Omega);E0:=factor(subs(Omega=0,enstrophy));

BRAND:=innerprod(jacobian(AA,[x,y,z]),VORTICITY);

stretch:=collect(factor(innerprod(BRAND,VORTICITY)),Omega);ST2:=factor(coeff(stretch,Omega^2));ST1:=factor(coeff(stretch,Omega));st0:=factor(subs(Omega=0,stretch));

BI:=collect(factor(innerprod(BRAND,BRAND)),Omega);BI2:=factor(coeff(BI,Omega^2));BI1:=factor(coeff(BI,Omega));BI0:=(factor(subs(Omega=0,BI)));

Theta:=collect(factor(-2*innerprod([EE1,EE2,EE3],VORTICITY)),Omega);

Theta:=collect(factor(-2*innerprod(grad(phi,[x,y,z]),VORTICITY)),Omega);t2:=factor(coeff(Theta,Omega^2));t1:=factor(coeff(Theta,Omega));

t0:=factor(subs(Omega=0,Theta));

t11:=(t1/(-8*z));collect(t11,D);t10:=factor(subs(D=0,t11));

beta:=D*x^2+A+B*z^2+D*y^2;betasq:=expand(beta^2);

simplify(t11-expand(-(D*x^2+D*y^2+A+B*z^2)^2));

TRANSCRITICAL HOPF OK

AA:= [x*(A-G+C*z)+Omega*y,y*(A-G+C*z)-Omega*x,A*z+B*z^2+D*(x^2+y^2)];

phi:=simplify((AA[1]^2+AA[2]^2+AA[3]^2));VV:=[AA[1],AA[2],AA[3],phi];

JAC:=simplify(jacobian(VV,[x,y,z,t]));collect(simplify(expand(charpoly(JAC,q))),q);

S1:=factor(trace(JAC)):S2:=factor(trace(innerprod(JAC,JAC))):S3:=factor(trace(innerprod(JAC,JAC,JAC))):

Xm:=S1;Yg:=collect(factor((1/2)*(Xm*S1-S2)),Omega);Za:=collect(factor((1/3)*(Yg*S1-Xm*S2+S3)),Omega);Tk:=factor(det(JAC));Yg0:=factor(subs(Omega=0,Yg));Za1:=factor(coeff(Za,Omega));Za0:=factor(subs(Omega=0,Za));

Hopf conditions

GGG:=factor(solve(Xm,G));ZaGGG:=factor(subs(G=GGG,Za));RRR:=solve(ZaGGG,Omega^2);YYg:=factor(subs(G=GGG,Omega^2=RRR,Yg));

Oscillations occur when YYg>0

RRR:=solve(ZaGGG,Omega^2);HopfCrit:=factor(subs(G=GGG,Omega^2=RRR,Yg));

Chiral - Brand invariants

VORTICITY:=evalm(curl(AA,[x,y,z])):CURLx:=simplify(VORTICITY[1]);CURLy:=simplify(VORTICITY[2]);CURLz:=simplify(VORTICITY[3]);

HELICITY:=collect(factor(innerprod(AA,VORTICITY)),Omega);

enstrophy:=collect(factor(innerprod(VORTICITY,VORTICITY)),Omega);E0:=factor(subs(Omega=0,enstrophy));

BRAND:=innerprod(jacobian(AA,[x,y,z]),VORTICITY);

stretch:=collect(factor(innerprod(BRAND,VORTICITY)),Omega);ST2:=factor(coeff(stretch,Omega^2));st0:=factor(subs(Omega=0,stretch));

BI:=collect(factor(innerprod(BRAND,BRAND)),Omega);BI2:=factor(coeff(BI,Omega^2));BI1:=factor(coeff(BI,Omega));BI0:=(factor(subs(Omega=0,BI)));

Theta:=factor(-2*innerprod(grad(phi,[x,y,z]),VORTICITY));

Theta:=collect(factor(-2*innerprod(grad(phi,[x,y,z]),VORTICITY)),Omega);t2:=factor(coeff(Theta,Omega^2));t1:=factor(coeff(Theta,Omega));

t0:=factor(subs(Omega=1,Theta)/(-8));

SADDLE NODE HOPF OK

AA:= [x*(G+C*z)+Omega*y,y*(G+C*z)-Omega*x,A+B*z^2+D*(x^2+y^2)];

phi:=simplify((AA[1]^2+AA[2]^2+AA[3]^2));VV:=[AA[1],AA[2],AA[3],phi];

JAC:=simplify(jacobian(VV,[x,y,z,t]));collect(simplify(expand(charpoly(JAC,q))),q);

S1:=factor(trace(JAC)):S2:=factor(trace(innerprod(JAC,JAC))):S3:=factor(trace(innerprod(JAC,JAC,JAC))):

Hopf conditions

GGG:=factor(solve(Xm,G));ZaGGG:=subs(G=GGG,Za);DDD:=solve(ZaGGG,D);YYg:=factor(subs(G=GGG,D=DDD,Yg));EV:=eigenvalues(JAC):EV3:=factor(subs(G=GGG,D=DDD,EV[3]));EV2:=simplify(subs(G=GGG,D=DDD,EV[2]));EV4:=factor(subs(G=GGG,D=DDD,EV[4]));factor(subs(G=GGG,D=DDD,EV[1]));factor(EV3*EV4*EV2);

Oscillations occur when YYg>0 Hence if 3(Omega)^2 exceeds (Bz)^2, oscillations occur.

Fixed points

RRR:=solve(ZaGGG,Omega^2);HopfCrit:=factor(subs(G=GGG,Omega^2=RRR,Yg));

Xm;

ZaCB:=(factor(subs(G=0,C=-B,Za)));DDD:=solve(ZaCB,D);(factor(subs(C=-B,G=0,Yg)));

Chiral - Brand invariants

VORTICITY:=evalm(curl(AA,[x,y,z])):CURLx:=simplify(VORTICITY[1]);CURLy:=simplify(VORTICITY[2]);CURLz:=simplify(VORTICITY[3]);

HELICITY:=collect(factor(innerprod(AA,VORTICITY)),Omega);

enstrophy:=collect(factor(innerprod(VORTICITY,VORTICITY)),Omega);E0:=factor(subs(Omega=0,enstrophy));

BRAND:=innerprod(jacobian(AA,[x,y,z]),VORTICITY);

stretch:=collect(factor(innerprod(BRAND,VORTICITY)),Omega);ST2:=factor(coeff(stretch,Omega^2));st0:=factor(subs(Omega=0,stretch));

BI:=collect(factor(innerprod(BRAND,BRAND)),Omega);BI2:=factor(coeff(BI,Omega^2));BI1:=factor(coeff(BI,Omega));BI0:=(factor(subs(Omega=0,BI)));

Theta:=factor(-2*innerprod(grad(phi,[x,y,z]),VORTICITY));

Theta:=collect(factor(-2*innerprod(grad(phi,[x,y,z]),VORTICITY)),Omega);t2:=factor(coeff(Theta,Omega^2));t1:=factor(coeff(Theta,Omega));

t0:=factor(subs(Omega=1,Theta)/8);ThetaHopf:=(factor(subs(G=GGG,D=DDD,Theta)));;FalacoCRIT:=YYg;boundinghopf:=factor(x^2*C^2-2*A*C+2*Omega^2+y^2*C^2+2*B^2*z^2-2*C*z^2*B);

bbb:=(4/3)*solve(YYg,Omega^2)/(z^2);H:=HELICITY;BBB:=solve(H,B);Tension:=subs(B=BBB,bbb);Rotat:=(Tension)^(1/2);

******************

HYSTERESIS HOPF OK

AA:= [x*(-G+C*z)+Omega*y,y*(-G+C*z)-Omega*x,A+B*z+E*z^3+D*(x^2+y^2)];

phi:=simplify((AA[1]^2+AA[2]^2+AA[3]^2));VV:=[AA[1],AA[2],AA[3],phi];

JAC:=simplify(jacobian(VV,[x,y,z,t]));collect(simplify(expand(charpoly(JAC,q))),q);

S1:=factor(trace(JAC)):S2:=factor(trace(innerprod(JAC,JAC))):S3:=factor(trace(innerprod(JAC,JAC,JAC))):

Hopf conditions

GGG:=factor(solve(Xm,G));ZaGGG:=subs(G=GGG,Za);DDD:=factor(solve(ZaGGG,D));YYg:=factor(subs(G=GGG,D=DDD,Yg));

Oscillations occur when YYg>0

The Hopf critical surface is

RRR:=solve(ZaGGG,Omega^2);HopfCrit:=factor(subs(G=GGG,Omega^2=RRR,Yg));

Chiral - Brand invariants

VORTICITY:=evalm(curl(AA,[x,y,z])):CURLx:=simplify(VORTICITY[1]);CURLy:=simplify(VORTICITY[2]);CURLz:=simplify(VORTICITY[3]);

HELICITY:=collect(factor(innerprod(AA,VORTICITY)),Omega);

enstrophy:=collect(factor(innerprod(VORTICITY,VORTICITY)),Omega);E0:=factor(subs(Omega=0,enstrophy));

BRAND:=innerprod(jacobian(AA,[x,y,z]),VORTICITY);

stretch:=collect(factor(innerprod(BRAND,VORTICITY)),Omega);ST2:=factor(coeff(stretch,Omega^2));st0:=factor(subs(Omega=0,stretch));

BI:=collect(factor(innerprod(BRAND,BRAND)),Omega);BI2:=factor(coeff(BI,Omega^2));BI1:=factor(coeff(BI,Omega));BI0:=(factor(subs(Omega=0,BI)));

Theta:=factor(-2*innerprod(grad(phi,[x,y,z]),VORTICITY));

Theta:=collect(factor(-2*innerprod(grad(phi,[x,y,z]),VORTICITY)),Omega);t1:=factor(coeff(Theta,Omega));

t0:=factor(subs(Omega=1,Theta)/8);

ThetaHopf:=factor(subs(G=GGG,D=DDD,Theta));boundingHopf:=factor(9*E^2*z^4-4*C*z^3*E+6*B*E*z^2-4*C*z*B-4*A*C+2*y^2*C^2+4*Omega^2+2*x^2*C^2+B^2);

H:=HELICITY;BBB:=z*solve(H,B);BBB*BBB;ZZZ:=Omega^2-3/4*(BBB*BBB);

YYg;

Generalized Langford

AA:= [x*(G+C*z)+Omega*y,y*(G+C*z)-Omega*x,P*sinh(alpha*z)+A-0*B*z+0*F*z^2+0*E*z^3+D*(x^2+y^2)];

#AA:= [x*(G+C*z)+Omega*y,y*(G+C*z)-Omega*x,P*(sinh(z))+D*(x^2+y^2)];

phi:=simplify((AA[1]^2+AA[2]^2+AA[3]^2));VV:=[AA[1],AA[2],AA[3],phi];

JAC:=simplify(jacobian(VV,[x,y,z,t]));collect(simplify(expand(charpoly(JAC,q))),q);

S1:=factor(trace(JAC)):S2:=factor(trace(innerprod(JAC,JAC))):S3:=factor(trace(innerprod(JAC,JAC,JAC))):

Hopf conditions

GGG:=factor(solve(Xm,G));ZaGGG:=subs(G=GGG,Za);DDD:=factor(solve(ZaGGG,D));YYg:=factor(subs(G=GGG,D=DDD,Yg));

Oscillations occur when YYg>0

The Hopf critical surface is

RRR:=solve(ZaGGG,Omega^2);HopfCrit:=factor(subs(G=GGG,Omega^2=RRR,Yg));Pcrit:=solve(HopfCrit,P);

Chiral - Brand invariants

VORTICITY:=evalm(curl(AA,[x,y,z])):CURLx:=simplify(VORTICITY[1]);CURLy:=simplify(VORTICITY[2]);CURLz:=simplify(VORTICITY[3]);

HELICITY:=collect(factor(innerprod(AA,VORTICITY)),Omega);

enstrophy:=collect(factor(innerprod(VORTICITY,VORTICITY)),Omega);E0:=factor(subs(Omega=0,enstrophy));

BRAND:=innerprod(jacobian(AA,[x,y,z]),VORTICITY);

stretch:=collect(factor(innerprod(BRAND,VORTICITY)),Omega);ST2:=factor(coeff(stretch,Omega^2));st0:=factor(subs(Omega=0,stretch));

BI:=collect(factor(innerprod(BRAND,BRAND)),Omega);BI2:=factor(coeff(BI,Omega^2));BI1:=factor(coeff(BI,Omega));BI0:=(factor(subs(Omega=0,BI)));

Theta:=factor(-2*innerprod(grad(phi,[x,y,z]),VORTICITY));

Theta:=collect(factor(-2*innerprod(grad(phi,[x,y,z]),VORTICITY)),Omega);t1:=factor(coeff(Theta,Omega));

t0:=factor(subs(Omega=1,Theta)/8);

ThetaHopf:=factor(subs(G=GGG,D=DDD,Theta));boundingHopf:=factor(9*E^2*z^4-4*C*z^3*E+6*B*E*z^2-4*C*z*B-4*A*C+2*y^2*C^2+4*Omega^2+2*x^2*C^2+B^2);

H:=HELICITY;BBB:=z*solve(H,P);BBB*BBB;ZZZ:=Omega^2-3/4*(BBB*BBB);

YYg;