Q2
***
allocatemem(2^33)
\p 3000
\r Hilbert.gp

ch="(...)/Theta-Q2/Split/";

MM=read("actionG24");
M=read(Str(ch,"../classeSLDG24-",2));

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

l=7;
beta=3+sqrt(2);
betac=3-sqrt(2);

Den(x,y)=(x^8*y^4 + x^8 + 4*x^6*y^6 - 4*x^6*y^2 + 6*x^4*y^8 - 10*x^4*y^4 + 4*x^2*y^10 - 4*x^2*y^6 + y^12 + y^8 + 16*y^4)^2*(x^12 + 4*x^10*y^2 + 6*x^8*y^4 + x^8 + 4*x^6*y^6 - 4*x^6*y^2 + x^4*y^8 - 10*x^4*y^4 + 16*x^4 - 4*x^2*y^6 + y^8)^2;


MCl=read(Str(ch,"PolQ2Theta-",l,"/classeSLDG24-D2G242G",l));
Tab=vector(8*(l+1),i,2);Tab[8*(l+1)-1]=1;

L1=vector(8*(l+1),i,read(Str(ch,"PolQ2Theta-",l,"/NumPhi/numphi",i)));
L2=vector(8*(l+1),i,read(Str(ch,"PolQ2Theta-",l,"/NumPsi/numpsi",i)));

[b2,b3]=[3+7*I,-5+2*I];
Phi=X^(8*(l+1))+sum(i=1,8*(l+1),substvec(L1[i],[x,y],[b2,b3])*X^(i-1)/Den(b2,b3)^Tab[i]);
Phip=Phi';
Psi=sum(i=1,8*(l+1),substvec(L2[i],[x,y],[b2,b3])/Den(b2,b3)^Tab[i]*X^(i-1))/Phip;

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

l=23;
beta=5+sqrt(2);
betac=5-sqrt(2);

Den(x,y)=( (x^36 + 12*x^34*y^2 + x^32*y^8 + 66*x^32*y^4 + 42*x^32 + 16*x^30*y^10 + 212*x^30*y^6 + 324*x^30*y^2 + 120*x^28*y^12 + 433*x^28*y^8 + 1010*x^28*y^4 + 1253*x^28 + 560*x^26*y^14 + 664*x^26*y^10 + 1108*x^26*y^6 + 11948*x^26*y^2 + 1820*x^24*y^16 + 1214*x^24*y^12 - 1483*x^24*y^8 + 51112*x^24*y^4 - 5567*x^24 + 4368*x^22*y^18 + 2608*x^22*y^14 + 3792*x^22*y^10 + 107524*x^22*y^6 - 50188*x^22*y^2 + 8008*x^20*y^20 + 2761*x^20*y^16 + 83024*x^20*y^12 + 43639*x^20*y^8 - 188718*x^20*y^4 + 20048*x^20 + 11440*x^18*y^22 - 5764*x^18*y^18 + 380496*x^18*y^14 - 440128*x^18*y^10 - 293724*x^18*y^6 + 110976*x^18*y^2 + 12870*x^16*y^24 - 29700*x^16*y^20 + 989611*x^16*y^16 - 1667858*x^16*y^12 + 197263*x^16*y^8 + 253600*x^16*y^4 + 10752*x^16 + 11440*x^14*y^26 - 62028*x^14*y^22 + 1716548*x^14*y^18 - 3687448*x^14*y^14 + 2122600*x^14*y^10 - 40064*x^14*y^6 - 3072*x^14*y^2 + 8008*x^12*y^28 - 81817*x^12*y^24 + 2109130*x^12*y^20 - 5992909*x^12*y^16 + 5730940*x^12*y^12 - 1127072*x^12*y^8 - 17920*x^12*y^4 + 4096*x^12 + 4368*x^10*y^30 - 74624*x^10*y^26 + 1885108*x^10*y^22 - 7531892*x^10*y^18 + 9530216*x^10*y^14 - 3418752*x^10*y^10 - 89088*x^10*y^6 + 1820*x^8*y^32 - 48122*x^8*y^28 + 1234931*x^8*y^24 - 7225716*x^8*y^20 + 11377967*x^8*y^16 - 6341632*x^8*y^12 + 254208*x^8*y^8 + 560*x^6*y^34 - 21688*x^6*y^30 + 586408*x^6*y^26 - 5031916*x^6*y^22 + 10667172*x^6*y^18 - 8571776*x^6*y^14 + 2048000*x^6*y^10 - 32768*x^6*y^6 + 120*x^4*y^36 - 6530*x^4*y^32 + 193876*x^4*y^28 - 2353135*x^4*y^24 + 7648338*x^4*y^20 - 8747568*x^4*y^16 + 3866624*x^4*y^12 + 139264*x^4*y^8 + 16*x^2*y^38 - 1184*x^2*y^34 + 40360*x^2*y^30 - 651768*x^2*y^26 + 3503092*x^2*y^22 - 5304064*x^2*y^18 + 2134016*x^2*y^14 - 32768*x^2*y^10 + y^40 - 98*y^36 + 4003*y^32 - 80066*y^28 + 718945*y^24 - 1281056*y^20 + 1024768*y^16 - 401408*y^12 + 65536*y^8)*(x^40 + 16*x^38*y^2 + 120*x^36*y^4 - 98*x^36 + 560*x^34*y^6 - 1184*x^34*y^2 + 1820*x^32*y^8 - 6530*x^32*y^4 + 4003*x^32 + 4368*x^30*y^10 - 21688*x^30*y^6 + 40360*x^30*y^2 + 8008*x^28*y^12 - 48122*x^28*y^8 + 193876*x^28*y^4 - 80066*x^28 + 11440*x^26*y^14 - 74624*x^26*y^10 + 586408*x^26*y^6 - 651768*x^26*y^2 + 12870*x^24*y^16 - 81817*x^24*y^12 + 1234931*x^24*y^8 - 2353135*x^24*y^4 + 718945*x^24 + 11440*x^22*y^18 - 62028*x^22*y^14 + 1885108*x^22*y^10 - 5031916*x^22*y^6 + 3503092*x^22*y^2 + 8008*x^20*y^20 - 29700*x^20*y^16 + 2109130*x^20*y^12 - 7225716*x^20*y^8 + 7648338*x^20*y^4 - 1281056*x^20 + 4368*x^18*y^22 - 5764*x^18*y^18 + 1716548*x^18*y^14 - 7531892*x^18*y^10 + 10667172*x^18*y^6 - 5304064*x^18*y^2 + 1820*x^16*y^24 + 2761*x^16*y^20 + 989611*x^16*y^16 - 5992909*x^16*y^12 + 11377967*x^16*y^8 - 8747568*x^16*y^4 + 1024768*x^16 + 560*x^14*y^26 + 2608*x^14*y^22 + 380496*x^14*y^18 - 3687448*x^14*y^14 + 9530216*x^14*y^10 - 8571776*x^14*y^6 + 2134016*x^14*y^2 + 120*x^12*y^28 + 1214*x^12*y^24 + 83024*x^12*y^20 - 1667858*x^12*y^16 + 5730940*x^12*y^12 - 6341632*x^12*y^8 + 3866624*x^12*y^4 - 401408*x^12 + 16*x^10*y^30 + 664*x^10*y^26 + 3792*x^10*y^22 - 440128*x^10*y^18 + 2122600*x^10*y^14 - 3418752*x^10*y^10 + 2048000*x^10*y^6 - 32768*x^10*y^2 + x^8*y^32 + 433*x^8*y^28 - 1483*x^8*y^24 + 43639*x^8*y^20 + 197263*x^8*y^16 - 1127072*x^8*y^12 + 254208*x^8*y^8 + 139264*x^8*y^4 + 65536*x^8 + 212*x^6*y^30 + 1108*x^6*y^26 + 107524*x^6*y^22 - 293724*x^6*y^18 - 40064*x^6*y^14 - 89088*x^6*y^10 - 32768*x^6*y^6 + 66*x^4*y^32 + 1010*x^4*y^28 + 51112*x^4*y^24 - 188718*x^4*y^20 + 253600*x^4*y^16 - 17920*x^4*y^12 + 12*x^2*y^34 + 324*x^2*y^30 + 11948*x^2*y^26 - 50188*x^2*y^22 + 110976*x^2*y^18 - 3072*x^2*y^14 + y^36 + 42*y^32+ 1253*y^28 - 5567*y^24 + 20048*y^20 + 10752*y^16 + 4096*y^12) )^2;


MCl=read(Str(ch,"PolQ2Theta-",l,"/classeSLDG24-D2G242G",l));
L1=vector(8*(l+1),i,read(Str(ch,"PolQ2Theta-",l,"/NumPhi/numphi",i)));
L2=vector(8*(l+1),i,read(Str(ch,"PolQ2Theta-",l,"/NumPsi/numpsi",i)));
Tab=vector(8*(l+1),i,2);Tab[8*(l+1)-1]=;1

[b2,b3]=[3+7*I,-5+2*I];
Phi=X^(8*(l+1))+sum(i=1,8*(l+1),substvec(L1[i],[x,y],[b2,b3])*X^(i-1)/Den(b2,b3)^Tab[i]);
Phip=Phi';
Psi=sum(i=1,8*(l+1),substvec(L2[i],[x,y],[b2,b3])/Den(b2,b3)^Tab[i]*X^(i-1))/Phip;

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

l=17;
beta=5+2*sqrt(2);
betac=5-2*sqrt(2);

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

l=41;
beta=7+2*sqrt(2);
betac=7-2*sqrt(2);

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

l=73;
beta=9+2*sqrt(2);
betac=9-2*sqrt(2);

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

l=89;
beta=11+4*sqrt(2);
betac=11-4*sqrt(2);

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

l=97;
beta=13+6*sqrt(2);
betac=13-6*sqrt(2);

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

DD=read(Str(ch,"PolQ2Theta-",l,"/DenQ2Theta",l));
Den(xx,yy)=substvec(DD,[x,y],[xx,yy]);

MCl=read(Str(ch,"PolQ2Theta-",l,"/classeSLDG24-D2G",l));
Tab=vector(2*(l+1),i,2);Tab[2*(l+1)]=1;

L1=vector(2*(l+1),i,read(Str(ch,"PolQ2Theta-",l,"/NumPhi/numphi",i)));
L2=vector(2*(l+1),i,read(Str(ch,"PolQ2Theta-",l,"/NumPsi/numpsi",i)));
[b2,b3]=[3+7*I,-5+2*I];

Phi=X^(2*(l+1))+sum(i=1,2*(l+1)-1,substvec(L1[i],[x,y],[b2,b3])*X^(i-1))/Den(b2,b3)^2+substvec(L1[2*(l+1)],[x,y],[b2,b3])/Den(b2,b3)*X^(2*(l+1)-1);
Phip=Phi';
Psi=sum(i=1,2*(l+1),substvec(L2[i],[x,y],[b2,b3])/Den(b2,b3)^Tab[i]*X^(i-1))/Phip;

Z=taufromthetNivD(b2,b3,M,MM,2)[1];
for(i=1,length(MCl),z=mtauD(MCl[i],Z,2);[b1p,b2p,b3p]=EvalBjp(TauToOmega([z[1]/beta,z[2]/betac],2));print(i"  "round(10^1000*subst(Phi,X,b2p))"  "round(10^1000*(b3p-subst(Psi,X,b2p)))));

for(i=1,length(MCl),z=mtauD(MCl[i],[Z[2],Z[1]],2);[b1p,b2p,b3p]=EvalBjp(TauToOmega([z[1]/beta,z[2]/betac],2));print(i"  "round(10^1000*subst(Phi,X,b2p))"  "round(10^1000*(b3p-subst(Psi,X,b2p)))));