poisson2.mws  030224

Determine the critical points of the first two moment closure equations, and their stabilities, under Poisson assumption!

>

restart;

>	with(VectorCalculus,Jacobian); with(LinearAlgebra,Eigenvalues);

The rhss of the differential equations for the first two cumulants are eq1 and eq2prel. I reserve the name eq2 for the rhs of the second moment closure equation.

>	eq1:=(a1*N-kappa1)*kappa1-kappa2;

>	eq2prel:=(f1*a1*N-f2*kappa1)*kappa1+(2*a1*N-f2)*kappa2-4*kappa1*kappa2-2*kappa3;

The Poisson assumption means that kappa3=kappa1:

>	eq2:=subs(kappa3=kappa1,eq2prel);

Solve the first equation for kappa2:

>	kap2:=solve(eq1,kappa2);

Substitute this value of kappa2 into the second equation.

>	eq2a:=subs(kappa2=kap2,eq2);

Denote the critical points by (x[1],y[1]), (x[2],y[2]), (x[3],y[3]).

>	x:=solve(eq2a,kappa1);

These are the kappa1-coordinates of the three critical points. The corresponding kappa2-coordinates are:

>	y[1]:=subs(kappa1=x[1],kap2);

>	y[2]:=expand(subs(kappa1=x[2],kap2));

>	y[3]:=expand(subs(kappa1=x[3],kap2));

Asymptotic approximations of the critical points:

>	assume(a1>0);

>	x2asymp:=map(simplify,convert(asympt(x[2],N,2),polynom));

>	x3asymp:=map(simplify,convert(asympt(x[3],N,2),polynom));

>	y2asymp:=map(simplify,convert(asympt(y[2],N,2),polynom));

>	y3asymp:=map(simplify,convert(asympt(y[3],N,2),polynom));

>

a1:='a1':

Study the stability of the three critical points:

>	F:=[eq1,eq2];

>	J:=Jacobian(F,[kappa1,kappa2]);

Determine the Jacobian in each of the three critical points:

>	J1:=subs(kappa1=x[1],kappa2=y[1],J);

>	J2:=subs(kappa1=x[2],kappa2=y[2],J);

>	J3:=subs(kappa1=x[3],kappa2=y[3],J);

The eigenvalues of these 3 matrices are:

>	eg1:=Eigenvalues(J1);

>	eg2:=Eigenvalues(J2);

>	eg3:=Eigenvalues(J3);

Asymptotic approximations of the eigenvalues:

>	assume(a1>0);

>	simplify(convert(asympt(eg1[1],N,1),polynom));

>	simplify(convert(asympt(eg1[2],N,1),polynom));

So both of the eigenvalues of the critical point (0,0) are positive. Hence (0,0) is an unstable node.

>	simplify(convert(asympt(eg2[1],N,1),polynom));

>	simplify(convert(asympt(eg2[2],N,1),polynom));

Both eigenvalues are negative. So this critical point is a stable node.

>	simplify(convert(asympt(eg3[1],N,1),polynom));

>	simplify(convert(asympt(eg3[2],N,1),polynom));

>

The third critical point has one positive and one negative eigenvalue. So it is a saddle point.