addhelp(hppack,"This package contains the following functions: anlt, clnolst, clnolstA, ellhp, ellvtautovz, ellvztohp, esse, ffindd, findd, hplst, hppolisin, hprint, hptopol, iscond, isinzd, maxd, polisinzd, qfbclgp, qfbclgplift, qfbprimeto, qfbvlift, qfbtotau, vectau, zlst") kill(clnolst) clnolst(K)=local(X);vector(K,X,if(isfundamental(-1*X),qfbclassno(-1*X),0)) addhelp(clnolst,"clnolst(K): Creates a vector of length K whose n-th entry is the class number of Q(sqrt(-n)) if -n is a fundamental discrimant, and 0 otherwise") kill(clnolstA) clnolstA(K)=local(X);vector(K,X,if( (!isfundamental(-1*X))&&((X%4==3)||(X%4==0)) ,qfbclassno(-1*X),0)) addhelp(clnolstA,"clnolstA(K): Creates a vector of length K whose n-th entry is the class number of the order of discriminant -n if -n is not fundamental, and 0 otherwise") kill(esse) esse(N,D)=local(i);for(i=0,4*N,if(issquare(D+4*N*i),return(round(sqrt(D+4*N*i))))) addhelp(esse,"esse(N,D): Finds a square root of D modulo 4N") kill(ffindd) ffindd(N,f=0,n=0,sf=0,C=1,K=30000)= local(t,lt,i,D,i2); i2=1; if( ((f==2)||(f==1))&&( (sf<=0)||(n<=0) ), t=factor(N); lt=matsize(t)[1]; n=prod(i=1,lt,t[i,1]); sf=prod(i=1,lt,if(t[i,2]==1 ,t[i,1],1));); if(n<=0, t=factor(N); lt=matsize(t)[1]; n=prod(i=1,lt,t[i,1]);); D=0; if( (f<0)||(f>4),f=0); if( f==0, for(i=3,K, if( isfundamental(-i)&&(gcd(N,i)==1), if(issquare(Mod(-i,4*N)), if((i2>=C), return(-i);, i2++ ))))); if( f==1, for(i=3,K, if( isfundamental(-i)&&( (gcd(N,i)>1)||(sf==1)), if(issquare(Mod(-i,4*N)), if((i2>=C), return(-i);,i2++))))); if( f==2, t=n/sf; for(i=1,floor(K/sf), if( isfundamental(-i*sf)&&(gcd(t,i)==1), if(issquare(Mod(-i*sf,4*N)), if((i2>=C), return(-i*sf);,i2++))))); if( f==3, for(i=3,K, if( isfundamental(-i)&&issquare(Mod(-i,4*N)), if((i2>=C), return(-i);,i2++)))); if( f==4, for(i=3,K, if( (!isfundamental(-i))&&((i%4==3)||(i%4==0)),if(issquarefree(gcd(i,N)) , if(issquare(Mod(-i,4*N)), if((i2>=C), return(-i);,i2++)))))); if(D==0, print("*** error: no D found")); return(D); addhelp(ffindd,"ffindd(N,{flag=0},{n},{sf},{C=1},{K}): Finds a discriminant D of an imaginary quadratic field, so that D is a square modulo 4N. N is a positive integer (the conductor of the elliptic curve), f is a flag, n is the product of the primes dividing N, sf is the product of the primes dividing N only once, K is a bound on D i.e. how far the function will search. If n and sf are not entered then they are computed. If flag=0 then gcd(D,N)=1, if flag=1 then D is such that gcd(D,N)>1 (when possible, i.e sf>1), and if flag=2 then gcd(D,N) is as big as possible. If flag=3 then no additional condition are made on gcd(D,N). If flag=4 then searches for a discriminant of a non-maximal quadratic imaginary order instead. Note that the D given is the C-th D found") kill(findd) findd(N,clno,f=0,n=0,sf=0,C=1)=local(K,t,lt,i,D,i2); i2=1; tmp=length(clno); K=length(clno); D=0; if( (f<0)||(f>9),f=0); if((f>=2)&&(f<=5)&&((sf<=0)||(n<=0)) ,t=factor(N); lt=matsize(t)[1]; n=prod(i=1,lt,t[i,1]); sf=prod(i=1,lt,if(t[i,2]==1 ,t[i,1],1))); if(n==0, t=factor(N); lt=matsize(t)[1]; n=prod(i=1,lt,t[i,1]);); if(f==0,for(i=3,K,if( (clno[i]!=0)&&(clno[i]<=tmp)&&(gcd(N,i)==1), if(issquare(Mod(-i,4*n)), tmp=clno[i]; D=-i; )))); if(f==1, for(i=3,K, if( (clno[i]!=0)&&(gcd(N,i)==1), if(issquare(Mod(-i,4*n)), if((i2>=C), return(-i);, i2++))))); if(f==2, for(i=3,K, if( (clno[i]!=0)&&(clno[i]<=tmp)&&((gcd(N,i)>1)||(sf==1)), if(issquare(Mod(-i,4*N)), tmp=clno[i]; D=-i; )))); if(f==3, for(i=3,K, if( (clno[i]!=0)&&((gcd(N,i)>1)||(sf==1)), if(issquare(Mod(-i,4*N)), if((i2>=C), return(-i); ,i2++))))); if(f==4,t=n/sf; for(i=1,floor(K/sf), if( (clno[sf*i]!=0)&&(clno[sf*i]<=tmp)&&(gcd(t,i)==1), if(issquare(Mod(-i*sf,4*N)), tmp=clno[i]; D=-i*sf;)))); if(f==5,t=n/sf; for(i=1,floor(K/sf), if( (clno[sf*i]!=0)&&(gcd(t,i)==1), if(issquare(Mod(-i*sf,4*N)), if((i2>=C), return(-i*sf) ; ,i2++ ))))); if(f==6, for(i=3,K, if( (clno[i]!=0)&&(clno[i]<=tmp), if(issquare(Mod(-i,4*N)), tmp=clno[i]; D=-i; )))); if(f==7, for(i=3,K, if( (clno[i]!=0), if(issquare(Mod(-i,4*N)), if((i2>=C), return(-i); ,i2++))))); if(f==8, for(i=3,K, if( (clno[i]!=0)&&(clno[i]<=tmp),if(gcd(-i/maxd(-i),N)==1, if(issquare(Mod(-i,4*N)), tmp=clno[i]; D=-i; ))))); if(f==9, for(i=3,K, if( (clno[i]!=0)&&(gcd(-i/maxd(-i),N)==1), if(issquare(Mod(-i,4*N)), if((i2>=C), return(-i); ,i2++))))); if(D==0, print("*** error: no D found in findd")); return(D); addhelp(findd,"findd(N,clno,{flag=0},{n},{sf},{C=1}): Finds a discriminant D of an imaginary quadratic order, so that D is a square modulo 4N, and the conductor is coprime to N. N is a positive integer (the conductor of the elliptic curve), clno is the result of the clnolst function, f is a flag, n is the product of the primes dividing N, sf is the product of the primes dividing N only once. If n and sf are not entered then they are computed. If the flag is 0 or 1 then gcd(D,N)=1, if the flag is 2 or 3 then gcd(D,N)>1 (when possible i.e. sf>1), and if the flag is 4 or 5 then gcd(D,N) is as big as possible. If the flag is 6 or 7 then no additional condition are made on gcd(D,N). If the flag is 8 or 9 then clno has to come from clnolstA rather than clnolst. If the flag is even then the output is chosen amongst the D's bigger than -length(clno) so that the class number is as small as possible (and in case of ties the smallest D is chosen). If the flag is odd then the output is the C-th D found.") kill(iscond) iscond(nlst)=local(tmp,i1,i2,i3);tmp=[]; for(i1=1,nlst[1][1]-1,tmp=concat(tmp,0)); tmp=concat(tmp,1); for(i2=2,length(nlst), if(nlst[i2-1][1]!=nlst[i2][1], for(i3=nlst[i2-1][1]+1,nlst[i2][1]-1, tmp=concat(tmp,0); ); tmp=concat(tmp,1); ) ); return(tmp); addhelp(iscond,"iscond(nlst): Makes a vector whose n-th entry is one if n is the conductor of an elliptic curve described in nlst, and 0 else. nlst should be a vector whose entry are vectors with a conductor as first component (e.g. lines from the cremona table: [conductor, isogeny class, ...]).") \\kill(xor) \\xor(x,y)=(x||y)&&(!(x&&y)) \\addhelp(xor,"xor(x,y): The exclusive OR of x and y.") kill(qfbclgp) qfbclgp(D)=local(k,a,b,c,tmp); k=floor(sqrt(-D/3)); tmp=[]; for(a=1,k,if(issquare(Mod(D,4*a)),forstep(b=-a+2-bitxor(a%2,D%2),a,2,if( ( (b^2-D)%(4*a)==0),c=(b^2-D)/(4*a);if( ((0<=b)&&(a<=c))||((0<-b)&&(a=h) , a=floor(l/2)*(-1)^l; f=j+s*a+N*a^2; fordiv(f,i3 , if( ((!l) || ((i3!=1)&&(i3!=f)&&(i3<=10+N)))&&(gcd(gcd(N*i3,s+2*N*a),f/i3)==1) , Q=Qfb(N*i3,s+2*N*a,f/i3); for(i4=1,h , if( (!clgp2[i4][2])||(component(Q,1)<=component(clgp2[i4][3],1)) , if(qfbred(Q)==clgp2[i4][1] , clgp2[i4][3]=Q; clgp2[i4][2]=1; ))))); l++; ); if((flag<0)||(flag>1),flag=0); if(flag==0, return(vector(h,i5,clgp2[i5][3]))); if(flag==1, return(vector(h,i5,[clgp2[i5][3],clgp2[i5][1]]))); addhelp(qfbvlift,"qfbvlift(s,D,N,qfbv,{flag=0}): Lifts qfbv (a vector of reduced binary quadratic forms) to equivalent forms Qfb(A,B,C) such that N divides A and 2N divides B-s. The result is a vector containing the lifted forms if the flag is 0, or pairs [reduced form, lifted form] if the flag is 1. These lifts do not always give the smallest A possible.") kill(qfbclgplift) qfbclgplift(s,D,N,f=0)= local(clgp,i,h,Q); h=Qfb(N,s,(s^2-D)/(4*N)); clgp=qfbclgp(D); for(i=1,length(clgp),clgp[i]=qfbprimeto(h,clgp[i])); vector(length(clgp),i, if(f==0, qfbcompraw(h,clgp[i]), Q=qfbcompraw(h,clgp[i]);[qfbred(Q),Q] )) addhelp(qfbclgplift,"qfbclgplift(s,D,N,{flag=0}): Lifts the whole class group of discriminant D to equivalent forms Qfb(A,B,C) such that N divides A and 2N divides B-s. The result is a vector containing the lifted forms if the flag is 0, or pairs [reduced form, lifted form] if the flag is 1. These lifts give the smallest A possible, but is not very fast.") kill(qfbprimeto) qfbprimeto(Q1,Q2)= local(s,D,N,A,B,C,A1,B1,C1,A2,B2,C2,tmp,tmp2,i,i2,i3,a,b,c,d,h); N=component(Q1,1); s=component(Q1,2); h=component(Q1,3); A=component(Q2,1); B=component(Q2,2); C=component(Q2,3); if( (s^2-4*N*h)==(B^2-4*A*C) , D=B^2-4*A*C ,error("The two Qfb in qfbprimeto have different discriminant")); A2=10*N^2; C2=1; B2=0; i2=0; i3=1; if( gcd(gcd(N,A),(B+s)/2)!=1 , until(i2==1 , tmp=qfminim([2*A,B;B,2*C],10^i3,10^i3,0); h=matsize(tmp[3])[2]; a=vector(h,i,tmp[3][1,i]); c=vector(h,i,tmp[3][2,i]); tmp2=vector(h,i,bezout(a[i],c[i])); b=vector(h,i,-tmp2[i][2]); d=vector(h,i,tmp2[i][1]); for(i=1,h , A1=a[i]^2*A+a[i]*c[i]*B+c[i]^2*C; B1= B+2*(a[i]*b[i]*A+b[i]*c[i]*B+c[i]*d[i]*C); if( (gcd(gcd(N,A1),(B1+s)/2)==1)&&(A1<=A2)&&(tmp2[i][3]==1) , A2=A1; B2=B1; C2= (B2^2-D)/(4*A2); ); ); if((A2==10*N^2)&&(C2==1)&&(B2==0) , i3++,i2=1);); A=A2; B=B2; C=C2;); return(Qfb(A,B,C)); addhelp(qfbprimeto,"qfbprimeto(Q1,Q2): Finds a binary quadratic form Q2' equivalent to Q2 but such that the Dirichlet composition between Q1 and Q2 is possible. The X^2 coefficient in Q2 is as small as possible. The method used is slow, for more speed use qfbvlift.") kill(qfbtotau) qfbtotau(Q)= local(A,B,C,D);A=component(Q,1); B=component(Q,2);C=component(Q,3);D=B^2-4*A*C;(-B+0.0+sqrt(D))/(2*A+0.0) addhelp(qfbtotau,"For a binary quadratic form Qfb(A,B,C), of discriminant D<0, gives the associated point in the upper half plane: (-B+sqrt(D))/2*A)") kill(ellvtautovz) ellvtautovz(vt,E)=local(prec,mt,i1,q0,n0,an,vz,i2,q,i3);prec=default(realprecision,,1); mt=vecmin(vector(length(vt),i1,imag(vt[i1]))); q0=exp(-2*Pi*mt); n0=anlt(q0,prec); an=ellan(E,n0); vz=vector(length(vt),i2,q=exp(2*Pi*I*vt[i2]);sum(i3=1,n0,(an[i3]+0.0)*(q^i3)/(i3+0.0))); return(vz) addhelp(ellvtautovz, "ellvtautovz(vt,E): Send a vector vt of points on the modular curve X_0(N)=H/Gamma_0(N) to a vector of points on the complex plane (in the lattice associated to the elliptic curve E).") kill(anlt) anlt(q1,pr)=local(q0,x0,x1,fct,fctp); q0=abs(q1); x0=log(10^(-pr)*(1-q0)/2)/log(q0); x1=x0-x0*log(x0)/(-2*log(q0)*x0+1); until( (x0-x1)<0.5,x0=x1;x1=x0-(x0+log(x0)/(-2*log(q0))-log(10^(-pr)*(1-q0)/2)/log(q0)) / (1+1/(-2*x0*log(q0))) ); return(floor(x1)+1) addhelp(anlt,"anlt(q,prec): Estimates the number of coefficients a_n of the L-function associated to some elliptic curve will be needed so that the sum of the a_i*q^i/i will be precise to the number of digits prec.") kill(ellvztohp) ellvztohp(vz,E)= local(i1);vector(length(vz),i1,ellztopoint(E,vz[i1])) addhelp(ellvztohp,"ellvztohp(vz,E): Sends a vector vz of points on the complex to the elliptic curve E by Weierstrass uniformization") kill(hptopol) hptopol(hp)=local(h,px,py);h=length(hp);pr=default(realprecision,,1);for(i=1,h,if(length(hp[i])==1,hp[i]=[10^(2*pr),10^(3*pr)])); px=prod(i1=1,h,if(abs(hp[i1][1])>=10^(pr-5),X^2,(X-hp[i1][1])));py=prod(i2=1,h,if(abs(hp[i2][2])>=10^(pr-5),Y^2,(Y-hp[i2][2])));[px,py] addhelp(hptopol,"hptopol(hp): Gives the polynomials satisfied by the x and y coordinates of the Heegner points hp. If these points happen to be [0] (the point at infinity) then the polynomial is multiplied by X^2 for each point, so that the degree of the polynomial will be twice the number of Heegner points.") kill(isinzd) isinzd(nb,D,prec)= local(re,im,im2,t,i1,i2,i3); if( abs(nb)>10^(default(realprecision,,1)-prec),print("*** warning: lack of precision in isinzd");return(1&&1) ); re=abs(real(nb+0.0)); im=abs(imag(nb+0.0)); im2=im/sqrt(-D); i2=( abs(2*im2-round(2*im2)) < (10^(-prec))); if( (D%4)==1,t=2,t=1); i1=( abs(t*re-round(t*re)) < (10^(-prec))); if( (D%4)==1, i3=( abs(re+im2-round(re+im2)) < (10^(-prec))),i3=1); (i1&&i2&&i3) addhelp(isinzd,"isinzd(z,D,prec): true(1) if the complex number z is in Z[sqrt(D)], for D<0, up to prec digits. false(0) if not.") kill(polisinzd) polisinzd(pol,D,prec,f=0)=local(n,i);if(f!=0,n=vecmin([poldegree(pol),abs(f)]),n=poldegree(pol));prod(i=0,n,isinzd(polcoeff(pol,poldegree(pol)-i),D,prec)) addhelp(polisinzd,"polisinzd(pol,D,prec,{flag=0}): true(1) if the polynomial pol has its coefficient is in Z[sqrt(D)], for D<0, up to prec digits. false(0) if not. If the flag is non-zero then only the first (flag+1)th coefficients are taken into account (those with higher power of X).") kill(hppolisin) hppolisin(hp,D,prec,f=0)=local(hppol,px,py,i1,i2,i3,tmp); hppol=hptopol(hp); px=hppol[1]; py=hppol[2]; i1=polisinzd(px,D,prec,f); i2=polisinzd(py,D,prec,f); tmp=poldegree(px)-length(hp); if(tmp==length(hp),return(2*tmp),i3=2*tmp); ((i1&&i2)+i3) addhelp(hppolisin,"hppolisin(hp,D,prec,{flag=0}): Checks whether the polynomial associated to the Heegner points hp has coefficients in Z[sqrt(D)], for D<0 (up to prec digits). The output is 0 if false, 1 if true and 2*h(D) if the Heegner point are the point at infinity [0]. If the flag is non-zero then only the first (flag+1)th coefficients in the polynomials associated to hp are taken into account (those with higher power of X). ") kill(vectau) vectau(clno,iscnd,f=0,C=1)=local(tmp,tmp2,N,D,i3,i1,i2,ngrp,vt,t,lt,s,sf); tmp=vector(length(iscnd),N, if(iscnd[N]==1, t=factor(N); lt=matsize(t)[1]; n=prod(i1=1,lt,t[i1,1]); sf=prod(i2=1,lt,if(t[i2,2]==1,t[i2,1],1) ); D=findd(N,clno,f,n,sf,C),0)); tmp2=vector(length(iscnd),N,D=tmp[N];if(D==0,0, s=esse(N,D); ngrp=qfbclgplift(s,D,N); vt=vector(length(ngrp),i3,qfbtotau(ngrp[i3])) ));[tmp,tmp2] addhelp(vectau,"vectau(clno,iscond,{flag=0},{C=1}): clno is a result from the clnolst function and iscond is the result of the iscond function. The output of this function is [Dlst,taulst] where Dlst is a vector whose n-th component is the discriminant satisfying the Heegner hypothesis for n (the discriminants listed are found by findd with additional condition using the flag and C), and taulst is a vector whose n-th component is the tau such that both the lattice <1,tau> and <1,n*tau> are elliptic curves with CM by the order of the n-th discriminant stated in Dlist.") kill(zlst) zlst(nlst,vtau)= local(tmp,i,N,vt,E);tmp=vector(length(nlst),i, if(i==1, N=nlst[i][1]; vt=vtau[N]; E=ellinit(nlst[i][4]); ellvtautovz(vt,E), if( (nlst[i-1][1]==nlst[i][1])&&(nlst[i-1][2]==nlst[i][2]), 0 , N=nlst[i][1]; vt=vtau[N]; E=ellinit(nlst[i][4]); ellvtautovz(vt,E) ))); for(i=1,length(tmp), if(tmp[i]==0, tmp[i]=tmp[i-1]) ); tmp addhelp(zlst, "zlst(nlst,taulst): Takes a cremona table (nlst) and taulst (from the vectau function) as input. The output are the z on the complex plane associated to the tau precedently computed and displayed as a vector whose n-th component gives the z related to the n-th curve in nlst") kill(hplst) hplst(nlst, zlist)=local(tmp,i,E,vz); tmp=vector(length(nlst),i,0); for(i=1,length(nlst), E=ellinit(nlst[i][4]); vz=zlist[i]; tmp[i]=ellvztohp(vz,E); ); return(tmp) addhelp(hplst, "hplst(nlst,zlst): Computes the Heegner points from the point on the complex plane. The result takes the form of a vector whose n-th entry is the Heegner point on the curve given by the n-th entry in nlst (a cremona table vector). zlst is the output of the zlst function.") kill(hprint) hprint(nlst,Dlst,hplist,prec,f=0,f2=0)= local(i,D,N); vector(length(hplist),i , N=nlst[i][1]; if( (f!=0) , D=maxd(Dlst[N]), D=Dlst[N]); hppolisin(hplist[i],D,prec,f2)) addhelp(hprint,"hprint(nlst,Dlst,hplst,prec,{flag=0},{flag2=0}): Checks if the Heegner points on the elliptic curve in hplist are "integral" (i.e. the polynomial of the field over which they are defined is in the order associated to them if flag=0, in the maximal order associated to them if flag=1) up to prec digits. Dlst is as in the output of vectau, hplst is the output of the hplst function, and nlst is a cremona table vector. If the flag2 is non-zero then only the first (flag+1)th coefficients in the polynomials associated to the Heegner points are taken into account (those with higher power of X). ") kill(maxd) maxd(D)= local(t,lt,i2,i); if(D%4==0 , t=factor(D/4); , t=factor(D); ); lt=matsize(t)[1]; i2=prod(i=1,lt,t[i,1]^(t[i,2]%2)); D= i2*if(i2%4==1,1,4); return(D); addhelp(maxd,"maxd(D): For D the discriminant of a quadratic order, gives back the discriminant of the maximal order.") kill(ellhp) ellhp(E,D=0,s=0,f=0,C=1,K=0,h=0,lt=1)= local(N,grp,ngrp,vt,i3,vz,hp); N=ellglobalred(E)[1]; if(K==0,K=vecmax([4*N,30000])) ;if(D==0, D=ffindd(N,f,,,C,K)); if(s==0, s=esse(N,D)); if( (s^2-D)%(4*N)!=0, error("s is not a square root of D modulo 4N")); ngrp=qfbclgplift(s,D,N); vt=vector(length(ngrp),i3,qfbtotau(ngrp[i3])); vz=ellvtautovz(vt,E); hp=ellvztohp(lt*vz,E); if(h==1,[hp,[N,D,s]],hp) addhelp(ellhp,"ellhp(E,{D},{s},{flag=0},{C=1},{K},{info=0},{manin=1}): Gives the Heegner points associated to the discriminant D and its square root (modulo 4N) s for the elliptic curve E. If D is omitted it is computed using ffindd with the options flag, C and K. If s is omitted then it is computed using esse. manin is the Manin constant associated to the elliptic curve; the points on the complex plane are multiplied by this before being sent to a point of E (it is believed to be one for strong Weil curves, but may be non trivial for other curves, e.g. 27A3 in Cremona's table). If info is one then the output is [hp,[N,D,s]] instead of just hp where hp is a vector containing the Heegner points, N is the condustor of E, D and s are as above. ") kill(ellhpq) ellhpq(E,f=0,C=1,K=0,lt=1)= local(N,D,s,ngrp,vt,i3,vz,hp); N=ellglobalred(E)[1]; if(K==0,K=vecmax([4*N,30000])) ;D=ffindd(N,f,,,C,K); s=esse(N,D); if( (s^2-D)%(4*N)!=0, error("s is not a square root of D modulo 4N")); ngrp=qfbclgplift(s,D,N); vt=vector(length(ngrp),i3,qfbtotau(ngrp[i3])); vz=ellvtautovz(vt,E); hp=ellvztohp(lt*vz,E); [hp,[N,D,s]] addhelp(ellhpq,"ellhpq(E,{flag=0},{C=1},{K},{manin=1}): A version of ellhp with less parameters. Gives the Heegner points associated to a discriminant D and its square root (modulo 4N) s for the elliptic curve E. D is computed using ffindd with the options flag, C and K. s is computed using esse. The output is [hp,[N,D,s]] where hp is a vector containing the Heegner points, N is the condustor of E, D and s are as mentionned above. manin is a number by which the points on the complex plane are multiplied before being sent to a point of E, is believed to be one for strong Weil curves.") kill(hpin) hpin(vhprint)= local(tmp,tmp2,tmp3,tmp4,tmp5,tmp6,i,i2,h,D); h=length(vhprint[1]); D=length(vhprint); tmp=[]; for(i=1,h ,tmp2=1; tmp3=0; for(i2=1,D,tmp2=tmp2&&(vhprint[i2][i]%2); tmp3=tmp3||(vhprint[i2][i]%2)); if(tmp2!=tmp3,tmp=concat(tmp,i))); tmp2=tmp; tmp=[]; for(i=1,h,tmp3=0; for(i2=1,D ,tmp3=tmp3||(vhprint[i2][i]>=2)); if(tmp3,tmp=concat(tmp,i))); tmp3=tmp; tmp=[]; for(i=1,h ,tmp4=1; tmp5=1; for(i2=1,D,tmp4=tmp4&&vhprint[i2][i]; tmp5=tmp5&&(vhprint[i2][i]%2)); if(tmp4!=tmp5,tmp=concat(tmp,i))); tmp7=tmp; tmp=[]; for(i=1,h,tmp4=1; for(i2=1,D,tmp4=tmp4&&(vhprint[i2][i]%2);); if(!(tmp4),tmp=concat(tmp,i))); tmp4=tmp; tmp=[]; for(i=1,h,tmp5=1; for(i2=1,D,tmp5=tmp5&&(vhprint[i2][i]%2)); if(tmp5,tmp=concat(tmp,i))); tmp5=tmp; tmp=[]; for(i=1,h,tmp6=1; for(i2=1,D, tmp6=tmp6&&(!(vhprint[i2][i]%2))); if(tmp6,tmp=concat(tmp,i))); tmp6=tmp; tmp=vector(h,i,vector(D,i2,vhprint[i2][i])); [tmp,tmp2,tmp7,tmp3,tmp4,tmp6,tmp5] addhelp(hpin,"hpin(vhprint): [hpin, mech, contrad, spec, non-ent, non-ent seul, ent]") kill(beau) beau(nlst,qqc)= tmp=""; for(i=1,length(qqc),tmp=concat(tmp,concat(nlst[qqc[i]][1],concat(nlst[qqc[i]][2],concat(nlst[qqc[i]][3],", ")))));tmp addhelp(beau, "beau(nlst, vec no de EC): string nom des EC") kill(combo) combo(nlst,vecfc,pr,K)= local(iscnd,qqc,clno,tmp,tmp2,tmp3,tmp4,i,i1,i2,i3,N,mt,vt,n0,q0); h=length(vecfc); iscnd=iscond(nlst); clno=clnolst(K); tmp2=vector(h,i,vectau(clno,iscnd,vecfc[i][1],vecfc[i][2])); qqc=vector(h,i,tmp2[i][2]); tmp2=vector(h,i,tmp2[i][1]); tmp=vector(h,i1,vector(length(nlst),i, 0)); for(i=1,length(tmp[1]), if(i==1 , N=nlst[i][1]; vt=vector(h,i1,qqc[i1][N]); mt=vecmin(vector(h,i1,vecmin(vector(length(vt[i1]),i2,imag(vt[i1][i2]))))); q0=exp(-2*Pi*mt); n0=anlt(q0,pr); E=ellinit(nlst[i][4]); an=ellan(E,n0); for(i1=1,h, tmp[i1][i]=vector(length(vt[i1]),i2,q=exp(2*Pi*I*vt[i1][i2]);n0=anlt(q,pr);sum(i3=1,n0,q^i3*an[i3]/i3))); , if( (nlst[i-1][1]==nlst[i][1])&&(nlst[i-1][2]==nlst[i][2]) , for(i1=1,h, tmp[i1][i]=tmp[i1][i-1] ) , N=nlst[i][1]; vt=vector(h,i1,qqc[i1][N]); mt=vecmin(vector(h,i1,vecmin(imag(vt[i1])))); q0=exp(-2*Pi*mt); n0=anlt(q0,pr); E=ellinit(nlst[i][4]); an=ellan(E,n0); for(i1=1,h, tmp[i1][i]=vector(length(vt[i1]),i2,q=exp(2*Pi*I*vt[i1][i2]);n0=anlt(q,pr);sum(i3=1,n0,q^i3*an[i3]/i3))); ))); tmp3=vector(h,i,hplst(nlst, tmp[i])); [tmp3,tmp2,tmp] addhelp(combo,"combo(nlsdt,vecteur [f,C], prec calc, K): [vec hplst, vec dlst, vec zlst]") kill(cop) cop(N,i)= local(i2,i3);if(i<=0,return(1));i3=0;for(i2=2,(i+1)*N,if((gcd(N,i2)==1),i3++;if(i3>=i,return(i2)))) kill(conclst) conclst(lst1,lst2)=local(tmp,tmp2,i1,i2); tmp2=[]; tmp=concat(lst1,lst2); for(i1=1,vecmax(tmp), for(i2=1,length(tmp), if(tmp[i2]==i1,tmp2=concat(tmp2,i1);break()))); tmp2 kill(combo2) combo2(nlst,prof,pr,pr2,K,f=0,C=1)= local(iscnd,qqc,ngrp,t,lt,n,sf,s,clno,tmp,tmp2,tmp3,tmp4,i,i1,i2,i3,i4,N,mt,vt,n0,q0); h=prof; iscnd=iscond(nlst); clno=clnolst(K); if(C<=0,C=1); tmp2=vector(prof,N, vector(length(iscnd),i,0)); for(N=1,length(iscnd), if(iscnd[N]==1, t=factor(N); lt=matsize(t)[1]; n=prod(i1=1,lt,t[i1,1]); sf=prod(i2=1,lt,if(t[i2,2]==1,t[i2,1],1) ); D=findd(N,clno,f,n,sf,C);for(i3=1,prof,tmp2[i3][N]=D*cop(N,i3-1)^2 ) )); qqc=vector(length(iscnd),N,vector(prof,i4,D=tmp2[i4][N];if(D==0,0, s=esse(N,D); ngrp=qfbclgplift(s,D,N); vt=vector(length(ngrp),i3,qfbtotau(ngrp[i3]))) )); tmp3=vector(prof,i,vector(length(iscnd),i2,qqc[i2][i])); qqc=tmp3;tmp=vector(prof,i1,vector(length(nlst),i, 0)); for(i=1,length(tmp[1]) , if(i==1 , N=nlst[i][1]; vt=vector(h,i1,qqc[i1][N]); mt=vecmin(vector(h,i1,vecmin(vector(length(vt[i1]),i2,imag(vt[i1][i2]))))); q0=exp(-2*Pi*mt); n0=anlt(q0,pr); E=ellinit(nlst[i][4]); an=ellan(E,n0); for(i1=1,h, tmp[i1][i]=vector(length(vt[i1]),i2,q=exp(2*Pi*I*vt[i1][i2]);n0=anlt(q,pr);sum(i3=1,n0,q^i3*an[i3]/i3))); , if( (nlst[i-1][1]==nlst[i][1])&&(nlst[i-1][2]==nlst[i][2]) , for(i1=1,h, tmp[i1][i]=tmp[i1][i-1] ) , N=nlst[i][1]; vt=vector(h,i1,qqc[i1][N]); mt=vecmin(vector(h,i1,vecmin(imag(vt[i1])))); q0=exp(-2*Pi*mt); n0=anlt(q0,pr); E=ellinit(nlst[i][4]); an=ellan(E,n0); for(i1=1,h , tmp[i1][i]=vector(length(vt[i1]),i2,q=exp(2*Pi*I*vt[i1][i2]); n0=anlt(q,pr);sum(i3=1,n0,q^i3*an[i3]/i3))); ))); tmp3=vector(h,i,hplst(nlst, tmp[i])); [tmp3,tmp2,tmp] addhelp(combo2,"combo2(nlst, no de cond diff, prec calc, K, flag, C)=[ vec hplst, vec dlst, vec zlst]") kill(combo3) combo3(nlst,vecfc,pr,K)= local(iscnd,qqc,clno,tmp,tmp2,tmp3,tmp4,i,i1,i2,i3,N,mt,vt,n0,q0); h=length(vecfc); iscnd=iscond(nlst); clno=clnolstA(K); tmp2=vector(h,i,vectau(clno,iscnd,vecfc[i][1],vecfc[i][2])); qqc=vector(h,i,tmp2[i][2]); tmp2=vector(h,i,tmp2[i][1]); tmp=vector(h,i1,vector(length(nlst),i, 0)); for(i=1,length(tmp[1]), if(i==1 , N=nlst[i][1]; vt=vector(h,i1,qqc[i1][N]); mt=vecmin(vector(h,i1,vecmin(vector(length(vt[i1]),i2,imag(vt[i1][i2]))))); q0=exp(-2*Pi*mt); n0=anlt(q0,pr); E=ellinit(nlst[i][4]); an=ellan(E,n0); for(i1=1,h, tmp[i1][i]=vector(length(vt[i1]),i2,q=exp(2*Pi*I*vt[i1][i2]);n0=anlt(q,pr);sum(i3=1,n0,q^i3*an[i3]/i3))); , if( (nlst[i-1][1]==nlst[i][1])&&(nlst[i-1][2]==nlst[i][2]) , for(i1=1,h, tmp[i1][i]=tmp[i1][i-1] ) , N=nlst[i][1]; vt=vector(h,i1,qqc[i1][N]); mt=vecmin(vector(h,i1,vecmin(imag(vt[i1])))); q0=exp(-2*Pi*mt); n0=anlt(q0,pr); E=ellinit(nlst[i][4]); an=ellan(E,n0); for(i1=1,h, tmp[i1][i]=vector(length(vt[i1]),i2,q=exp(2*Pi*I*vt[i1][i2]);n0=anlt(q,pr);sum(i3=1,n0,q^i3*an[i3]/i3))); ))); tmp3=vector(h,i,hplst(nlst, tmp[i])); [tmp3,tmp2,tmp] addhelp(combo3,"combo3(nlsdt,vecteur [f,C], prec calc, K): [vec hplst, vec dlst, vec zlst]") kill(ellsum) ellsum(E,vecfc)=local(tmp); tmp=[0]; for(i=1,length(veccf),tmp=elladd(E,tmp,vecfc[i])); tmp addhelp(ellsum, "ellsum(E, vecpt): Computes the sum of the points on E in the vector vecpt") \\kill(ztozd) kill(hppk) hppk(nlst,N,pr)=local(tmp,tmp2,tmp3,tmp4,tmp5,h,D,i,i1,i2,i3); tmp5=vector(length(N[3][1]),i1 , E=ellinit(nlst[i1][4]); vector(length(N[3]),i3, tmp= sum(i2=1,length(N[3][i3][i1]),N[3][i3][i1][i2]); tmp2= ellztopoint( E , tmp); if( (length(tmp2)==1)||(abs(tmp2[1])>10^default(realprecision,,1)), tmp2=[0];i=2, tmp3= isinzd(tmp2[1],N[2][i3][nlst[i1][1]],pr); tmp4= isinzd(tmp2[2],N[2][i3][nlst[i1][1]],pr); i=tmp3&&tmp4); h=ordzapp(E,tmp,pr); D=(sum(i4=1,length(tmp2),abs(imag(tmp2[i4])))<10^-pr); [D,h,i,tmp2,tmp] )); tmp3=vector(5,i1,vector(length(tmp5),i2,vector(length(tmp5[1]),i3,tmp5[i2][i3][i1]))) addhelp(hppk,"hppk(nlst, res combo, pr)= [vec P_K reel, vec ordre de P_K, vec P_K est ent, vec P_K, vec z_K]") kill(ordzapp) ordzapp(E,tmp,pr)= local(i1,i2,tmp2,tmp3,A1,A2,B1,B2,C1,C2); i1=E[15]; i2=E[16]; A1=bestappr(imag(tmp)/imag(i2),10^floor(pr/2)); if((round(A1)-A1)!=0, B1=component(A1,1); C1=component(A1,2);, B1=A1;C1=1;); tmp2=C1*tmp-B1*i2; if( abs(imag(tmp2))>10^-pr, return(0)); A2=bestappr(real(tmp2)/real(i1),10^floor(pr/2)); if((round(A2)-A2)!=0, B2=component(A2,1); C2=component(A2,2); , B2=A2;C2=1); tmp3=C2*tmp2-B2*i1; if( abs(tmp3)>10^-pr, return(0),return(C1*C2)) addhelp(ordzapp,"ordzapp(E,z,bound): Checks if the point z is a torsion element, on the torus associated to E, of order less than 10^[bound/2]. Returns the order of z if z is torsion, and 0 otherwise.") kill(isinlst) isinlst(i,vecfc)=local(tmp); tmp=0; for(i2=1,length(vecfc),tmp=(tmp||(i==vecfc[i2])); if(tmp==1,return(tmp)) );tmp kill(difflst) difflst(tmp,tmp2)=local(); tmp3=[];for(i=1,vecmax(tmp),if(isinlst(i,tmp)&&(!isinlst(i,tmp2)),tmp3=concat(tmp3,i)));tmp3