Construction of the two invariants (Section 6)

We have 4 equations, since the series Y_i cancel both the kernel Ker and the right-hand side R

>    e1:=subs(y=y1,Ker);e2:=subs(y=y2,Ker);

e1 := -x*y1*w*t*((nu-1)*(y1-1)+q*y1)*M(1,y1)-x*y1*t*(x*nu-1)*M(x,1)-x^2*y1*w*t*(nu-1)/(x-1)-x*y1^2*t/(y1-1)+1

e2 := -x*y2*w*t*((nu-1)*(y2-1)+q*y2)*M(1,y2)-x*y2*t*(x*nu-1)*M(x,1)-x^2*y2*w*t*(nu-1)/(x-1)-x*y2^2*t/(y2-1)+1

>    e3:=subs(y=y1,R);e4:=subs(y=y2,R);

e3 := 1-x*y1*w*t*(nu-1)*M(1,y1)/(x-1)-x*y1*t*M(x,1)/(y1-1)

e4 := 1-x*y2*w*t*(nu-1)*M(1,y2)/(x-1)-x*y2*t*M(x,1)/(y2-1)

We first solve e3 and e4 for x and M(x,1). (these are the two simplest equations).

Set  xM(x,1)=S(x).

>    e3b:=subs(M(x,1)=Sx/x,e3);
e4b:=subs(M(x,1)=Sx/x,e4);

e3b := 1-x*y1*w*t*(nu-1)*M(1,y1)/(x-1)-y1*t*Sx/(y1-1)

e4b := 1-x*y2*w*t*(nu-1)*M(1,y2)/(x-1)-y2*t*Sx/(y2-1)

This is Eq. (34)

>    Sxsol1:=factor(solve(subs(isolate(e3b,x),e4b),Sx));

Sxsol1 := (y2-1)*(y1-1)*(-y1*M(1,y1)+y2*M(1,y2))/y1/(y2*M(1,y2)-M(1,y2)-y1*M(1,y1)+M(1,y1))/y2/t

>    Xsol1:=factor(solve(subs(isolate(e3b,Sx),e4b),x));

Xsol1 := (y1-y2)/(y1+y1*w*t*nu*y2^2*M(1,y2)-y1*w*t*nu*y2*M(1,y2)-y1*w*t*y2^2*M(1,y2)+y1*w*t*y2*M(1,y2)-y2-y1^2*w*t*M(1,y1)*nu*y2+y1*w*t*M(1,y1)*nu*y2+y1^2*w*t*M(1,y1)*y2-y1*w*t*M(1,y1)*y2)
Xsol1 := (y1-y2)/(y1+y1*w*t*nu*y2^2*M(1,y2)-y1*w*t*nu*y2*M(1,y2)-y1*w*t*y2^2*M(1,y2)+y1*w*t*y2*M(1,y2)-y2-y1^2*w*t*M(1,y1)*nu*y2+y1*w*t*M(1,y1)*nu*y2+y1^2*w*t*M(1,y1)*y2-y1*w*t*M(1,y1)*y2)

The reciprocal of x looks nicer

>    Xinvsol1:=collect(1/Xsol1,M,factor);

Xinvsol1 := y1*w*t*y2*(nu-1)*(y2-1)/(y1-y2)*M(1,y2)+1-y1*w*t*y2*(nu-1)*(y1-1)/(y1-y2)*M(1,y1)

>   

Let us now work with e1 and e2

We  eliminate M(x,1), and this gives another linear equation in x:

>    eqx1:=collect(op(3,factor(numer(subs(isolate(e1,M(x,1)),e2)))),x,factor);

eqx1 := y1*y2*t*(y1-y2-w*nu*M(1,y2)+w*M(1,y1)*nu-y1*w*y2^2*M(1,y2)+y1^2*w*M(1,y1)*y2+2*y1*w*y2*M(1,y2)-2*y1*w*M(1,y1)*y2+2*w*nu*y2*M(1,y2)-w*M(1,y1)*nu*y2-y2^2*w*M(1,y2)*q-w*nu*y2^2*M(1,y2)+y2*w*M(1,y2...
eqx1 := y1*y2*t*(y1-y2-w*nu*M(1,y2)+w*M(1,y1)*nu-y1*w*y2^2*M(1,y2)+y1^2*w*M(1,y1)*y2+2*y1*w*y2*M(1,y2)-2*y1*w*M(1,y1)*y2+2*w*nu*y2*M(1,y2)-w*M(1,y1)*nu*y2-y2^2*w*M(1,y2)*q-w*nu*y2^2*M(1,y2)+y2*w*M(1,y2...
eqx1 := y1*y2*t*(y1-y2-w*nu*M(1,y2)+w*M(1,y1)*nu-y1*w*y2^2*M(1,y2)+y1^2*w*M(1,y1)*y2+2*y1*w*y2*M(1,y2)-2*y1*w*M(1,y1)*y2+2*w*nu*y2*M(1,y2)-w*M(1,y1)*nu*y2-y2^2*w*M(1,y2)*q-w*nu*y2^2*M(1,y2)+y2*w*M(1,y2...
eqx1 := y1*y2*t*(y1-y2-w*nu*M(1,y2)+w*M(1,y1)*nu-y1*w*y2^2*M(1,y2)+y1^2*w*M(1,y1)*y2+2*y1*w*y2*M(1,y2)-2*y1*w*M(1,y1)*y2+2*w*nu*y2*M(1,y2)-w*M(1,y1)*nu*y2-y2^2*w*M(1,y2)*q-w*nu*y2^2*M(1,y2)+y2*w*M(1,y2...
eqx1 := y1*y2*t*(y1-y2-w*nu*M(1,y2)+w*M(1,y1)*nu-y1*w*y2^2*M(1,y2)+y1^2*w*M(1,y1)*y2+2*y1*w*y2*M(1,y2)-2*y1*w*M(1,y1)*y2+2*w*nu*y2*M(1,y2)-w*M(1,y1)*nu*y2-y2^2*w*M(1,y2)*q-w*nu*y2^2*M(1,y2)+y2*w*M(1,y2...

>    degree(%,x);

1

Second expression of the recirpocal of x

>    Xinvsol2:=collect(1/solve(eqx1,x),M,factor);

Xinvsol2 := (nu*y2-nu-y2+1+q*y2)*y1*w*t*y2/(y1-y2)*M(1,y2)-(y1*nu-y1-nu+1+q*y1)*y1*w*t*y2/(y1-y2)*M(1,y1)+y1*y2*t/(y2-1)/(y1-1)

>   

>   

Compare these two epressions of 1/x : this will give the first  invariant

>    collect((Xinvsol1-Xinvsol2),M,factor);

-y2^2*w*t*y1*q/(y1-y2)*M(1,y2)+y1^2*w*t*y2*q/(y1-y2)*M(1,y1)-(-y1*y2+y1+y2-1+y1*y2*t)/(y2-1)/(y1-1)

>    eq5:=collect((Xinvsol1-Xinvsol2)/y1/y2*(y1-y2),M,factor);

eq5 := -y2*w*t*q*M(1,y2)+y1*w*t*q*M(1,y1)-(-y1*y2+y1+y2-1+y1*y2*t)*(y1-y2)/(y2-1)/(y1-1)/y1/y2

>    convert(op(3,eq5),parfrac,y1);

-(-y2+1+t*y2)/(y2-1)/y2+t/(y1-1)-1/y1

>    convert(%,parfrac,y2);

t/(y1-1)-1/y1-t/(y2-1)+1/y2

>    eq6:=subsop(3=%,eq5);

eq6 := -y2*w*t*q*M(1,y2)+y1*w*t*q*M(1,y1)+t/(y1-1)-1/y1-t/(y2-1)+1/y2

>    factor(eq5/eq6);

1

We can read off  our first  invariant:

>    I10:=y*w*t*q*M(1,y)-1/y+t/(y-1);

I10 := y*w*t*q*M(1,y)-1/y+t/(y-1)

Let us check

>    factor((subs(y=y1,I10)-subs(y=y2,I10))/eq6);

1

We can just as well take

>    I1:=w*t*y*q*M(1,y)-(1-y)/y+t*y/(y-1);

I1 := y*w*t*q*M(1,y)-(1-y)/y+t*y/(y-1)

Indeed..

>    factor(I1-I10);

1+t

>   

>   

Expression of M(1,y) in terms of the invariant I (denoted II). This is Eq. (36).

>    Mysol:=factor(isolate(I1-II,M(1,y)));

Mysol := M(1,y) = -(y^2-2*y+1+y^2*t-II*y^2+II*y)/y^2/(y-1)/w/t/q

We now derive from the equation e1  another expression of S(x):

>    Sxsol2:=factor(solve(subs(M(x,1)=Sx/x,x=1/Xinvsol2, subs(y=y1,Mysol),subs(y=y2,Mysol),e1),Sx));

Sxsol2 := -(2*y1*y2-II*y2^2*y1-II*y2*y1-2*y1^2*y2^2+y1*y2^2*t-3*y1^2*nu*y2^2*II-3*y1^2*q*y2^2*II+2*y1^2*nu*y2^2*t+y1^2*nu*II*y2+q*y1^2*II*y2-y1^2*t*nu*y2+2*y1^2*nu*y2^2+2*y1^2*q*y2^2+3*y1^2*II*y2^2-y1^...
Sxsol2 := -(2*y1*y2-II*y2^2*y1-II*y2*y1-2*y1^2*y2^2+y1*y2^2*t-3*y1^2*nu*y2^2*II-3*y1^2*q*y2^2*II+2*y1^2*nu*y2^2*t+y1^2*nu*II*y2+q*y1^2*II*y2-y1^2*t*nu*y2+2*y1^2*nu*y2^2+2*y1^2*q*y2^2+3*y1^2*II*y2^2-y1^...
Sxsol2 := -(2*y1*y2-II*y2^2*y1-II*y2*y1-2*y1^2*y2^2+y1*y2^2*t-3*y1^2*nu*y2^2*II-3*y1^2*q*y2^2*II+2*y1^2*nu*y2^2*t+y1^2*nu*II*y2+q*y1^2*II*y2-y1^2*t*nu*y2+2*y1^2*nu*y2^2+2*y1^2*q*y2^2+3*y1^2*II*y2^2-y1^...
Sxsol2 := -(2*y1*y2-II*y2^2*y1-II*y2*y1-2*y1^2*y2^2+y1*y2^2*t-3*y1^2*nu*y2^2*II-3*y1^2*q*y2^2*II+2*y1^2*nu*y2^2*t+y1^2*nu*II*y2+q*y1^2*II*y2-y1^2*t*nu*y2+2*y1^2*nu*y2^2+2*y1^2*q*y2^2+3*y1^2*II*y2^2-y1^...

...and compare it with the first expression of S(x) we have obtained:

>    eq7:=map(factor,collect(numer(factor(subs( subs(y=y1,Mysol),subs(y=y2,Mysol),Sxsol1-Sxsol2))),[y1,y2],distributed));

eq7 := -(q-2)*(nu-1)*y1*y2+(nu-1)*y1^2+(nu*II^2-II^2+q*II^2-3*q*II-4*nu*II+4*II+q^2*w*t+q*t*nu-q*t+3*q-q*nu+4*nu-4)*y1^2*y2^2+(nu-1)*y2^2+(2*nu*II-4*nu-2*II+4+q*II-2*q+q*nu)*y1*y2^2+(2*nu*II-4*nu-2*II+...
eq7 := -(q-2)*(nu-1)*y1*y2+(nu-1)*y1^2+(nu*II^2-II^2+q*II^2-3*q*II-4*nu*II+4*II+q^2*w*t+q*t*nu-q*t+3*q-q*nu+4*nu-4)*y1^2*y2^2+(nu-1)*y2^2+(2*nu*II-4*nu-2*II+4+q*II-2*q+q*nu)*y1*y2^2+(2*nu*II-4*nu-2*II+...

Compare with the equation in the paper (just above (37)) :

>    (nu-1)*(y1^2+y2^2-(q-2)*y1*y2)+y1*y2*(y1+y2)*((q+2*nu-2)*(II-2)+q*nu)
+y1^2*y2^2*((q+nu-1)*II^2-(3*q+4*nu-4)*II+q*t*(w*q+nu-1)+2*q-q*(nu-1)+4*nu-4);

(nu-1)*(y1^2+y2^2-(q-2)*y1*y2)+y1*y2*(y1+y2)*((2*nu-2+q)*(II-2)+q*nu)+y1^2*y2^2*((nu-1+q)*II^2-(3*q+4*nu-4)*II+q*t*(w*q+nu-1)+2*q-q*(nu-1)+4*nu-4)
(nu-1)*(y1^2+y2^2-(q-2)*y1*y2)+y1*y2*(y1+y2)*((2*nu-2+q)*(II-2)+q*nu)+y1^2*y2^2*((nu-1+q)*II^2-(3*q+4*nu-4)*II+q*t*(w*q+nu-1)+2*q-q*(nu-1)+4*nu-4)

>    factor(%/eq7);

1

Replace y_i by 1/z_i

>    eq8:=map(factor,collect(numer(factor(subs(y1=1/z1,y2=1/z2,eq7))),[z1,z2],distributed));

eq8 := -(q-2)*(nu-1)*z1*z2+(nu-1)*z2^2+nu*II^2-II^2+q*II^2-3*q*II-4*nu*II+4*II+q^2*w*t+q*t*nu-q*t+3*q-q*nu+4*nu-4+(nu-1)*z1^2+(2*nu*II-4*nu-2*II+4+q*II-2*q+q*nu)*z1+(2*nu*II-4*nu-2*II+4+q*II-2*q+q*nu)*...
eq8 := -(q-2)*(nu-1)*z1*z2+(nu-1)*z2^2+nu*II^2-II^2+q*II^2-3*q*II-4*nu*II+4*II+q^2*w*t+q*t*nu-q*t+3*q-q*nu+4*nu-4+(nu-1)*z1^2+(2*nu*II-4*nu-2*II+4+q*II-2*q+q*nu)*z1+(2*nu*II-4*nu-2*II+4+q*II-2*q+q*nu)*...

>    collect(subs(z1=v1+a,z2=v2+a,eq8/(nu-1)),[v1,v2],distributed);;

v1^2+(2-q)*v1*v2+(-(q-2)*(nu-1)*a+4+2*nu*II-4*nu-2*II+2*(nu-1)*a+q*II-2*q+q*nu)/(nu-1)*v1+v2^2+(-(q-2)*(nu-1)*a+4+2*nu*II-4*nu-2*II+2*(nu-1)*a+q*II-2*q+q*nu)/(nu-1)*v2+(-(q-2)*(nu-1)*a^2-q*t+4*nu+2*(nu...
v1^2+(2-q)*v1*v2+(-(q-2)*(nu-1)*a+4+2*nu*II-4*nu-2*II+2*(nu-1)*a+q*II-2*q+q*nu)/(nu-1)*v1+v2^2+(-(q-2)*(nu-1)*a+4+2*nu*II-4*nu-2*II+2*(nu-1)*a+q*II-2*q+q*nu)/(nu-1)*v2+(-(q-2)*(nu-1)*a^2-q*t+4*nu+2*(nu...
v1^2+(2-q)*v1*v2+(-(q-2)*(nu-1)*a+4+2*nu*II-4*nu-2*II+2*(nu-1)*a+q*II-2*q+q*nu)/(nu-1)*v1+v2^2+(-(q-2)*(nu-1)*a+4+2*nu*II-4*nu-2*II+2*(nu-1)*a+q*II-2*q+q*nu)/(nu-1)*v2+(-(q-2)*(nu-1)*a^2-q*t+4*nu+2*(nu...

>    asol:=solve(coeff(coeff(%,v1,1),v2,0),a);

asol := (2*nu*II-4*nu-2*II+4+q*II-2*q+q*nu)/(nu-1)/(-4+q)

Expression of the constant a given in the paper:

>    1-((q+2*nu-2)*II-q)/(nu-1)/(4-q);

1-((2*nu-2+q)*II-q)/(nu-1)/(4-q)

>    factor(asol-%);

0

>    eq9:=subs(Z=1,map(factor,collect(subs(z1=v1+asol,z2=v2+asol,Z*eq8/(nu-1)),[v1,v2],distributed)));

eq9 := v1^2-(q-2)*v1*v2+v2^2+(-q^2*w*t+q^2*w*t*nu-2*q*t*nu-q*II+q*II^2*nu-4*q*w*t*nu+q+4*w*q*t+q*t+q*t*nu^2-q*II*nu+8*t*nu-4*t-4*t*nu^2-2*nu*II^2+II^2+nu^2*II^2)*q/(nu-1)^2/(-4+q)
eq9 := v1^2-(q-2)*v1*v2+v2^2+(-q^2*w*t+q^2*w*t*nu-2*q*t*nu-q*II+q*II^2*nu-4*q*w*t*nu+q+4*w*q*t+q*t+q*t*nu^2-q*II*nu+8*t*nu-4*t-4*t*nu^2-2*nu*II^2+II^2+nu^2*II^2)*q/(nu-1)^2/(-4+q)

Observe that the coefficient of v1v2 is  -(q-2).

We will have to take the square root of:

>    DD:=collect(normal(-op(4,eq9)/q*(4-q)*(nu-1)^2),[II,t],factor);

DD := (q*nu-2*nu+1+nu^2)*II^2-q*(nu+1)*II+(nu-1)*(-4+q)*(w*q+nu-1)*t+q

(Here, yinv stands for 1/y

>    UU:=(yinv-asol)*sqrt(4-q)*(nu-1)/sqrt(q)/sqrt(C);

UU := (yinv-(2*nu*II-4*nu-2*II+4+q*II-2*q+q*nu)/(nu-1)/(-4+q))*(4-q)^(1/2)*(nu-1)/q^(1/2)/C^(1/2)

The value x where we will take T_m  (x)

>    xsubs:=collect(factor(UU*sqrt(q)/2*sqrt(4-q)),[yinv,II],factor);

xsubs := -1/2*(nu-1)*(-4+q)/C^(1/2)*yinv+1/2*(2*nu-2+q)/C^(1/2)*II+1/2*(-2*q-4*nu+q*nu+4)/C^(1/2)

Expression of x in the paper

>    ((nu-1)*(4-q)*(yinv-1)+(q+2*nu-2)*II-q)/2/sqrt(C);

1/2*((nu-1)*(4-q)*(yinv-1)+(2*nu-2+q)*II-q)/C^(1/2)

>    factor(xsubs-%);

0

>   

>