\\ This file contains the PARI calculations that were performed \\ in section 5 of the paper \\ "Integration on H_p x H and arithmetic applications" \\ \\ \\ The calculation can be divided into three parts: \\ \\ Part I: Initialization of the data concerning X_0(11) \\ Part 2: Initialization of the data on the real quadratic order \\ Part 3: Computing the points \\ Global variables: ggp and ggn. \\ Calculations on the curve will be done modulo ggp^ggn. ggp = 11; ggn = 3; \\ ********************************************************************* \\ ** PART I : The curve X_0(11) \\ ********************************************************************* \\ ggmsymbol = ([C, 0, A, B, B-A, -A, -A, -B, A-B, A, 0, -C]);\ \\ gge = initell([0,-1,1,-10,-20]);\ ggmsymbol = ([[0,0,1], [0,0,0], [1,0,0], [0,1,0], [-1,1,0], [-1,0,0], [-1,0,0], [0,-1,0], [1,-1,0], [1,0,0], [0,0,0], [0,0,-1]]);\ \\ ********************************************************************* \\ ** Modular symbol routines \\ ********************************************************************* \\ \\ path(a,b,c,d): computes the symbol { a/b, c/d } as in Cremona, i.e. \\ the path from (a/b) to (c/d) path(a,b,c,d)= pathinfty(c/d) - pathinfty(a/b) \\ pathinfty(a/b): computes the symbol : { infty, c/d } pathinfty(x,\ c,ll,res,qleft,qmiddle,qright,sgn,j) = \ c = cf(x) ;\ ll=length(c);\ res=0;\ qleft = 1; qmiddle = 0;\ sgn = 1;\ for(j=1,ll,\ qright = c[j]*qmiddle + qleft;\ if(sgn==1,res = res - smb(qmiddle,qright),res = res + smb(qright,qmiddle));\ sgn=-sgn;\ qleft=qmiddle; qmiddle = qright);\ res; smb(c,d)= \ if( (d% ggp)==0, ggmsymbol[ggp+1], ggmsymbol[((c/d)% ggp)+1] ); \\ ********************************************************************** \\ ** Initialization routine: InitForms(D), where D is a (not necessarily \\ ** fundamental) discriminant. This sets up certain global \\ ** variables: \\ ** ggD: The discriminant D of the order. \\ ** ggw: The element such that [1,gw] is a basis for the order. \\ ** ggUnit: The fundamental unit of norm 1 in the order of \\ ** discriminant D. \\ ** ggMUnit: The matrix giving the multiplication by the unit, in the \\ ** basis [1, sqrt{D}] \\ ** ggH: The class number of the order. \\ ** ggClassGroup: The class group of the order, as an array of \\ ** representative forms. \\ ** ggi: an index which divides 12 = 11+1, and is used to make the unit \\ ** congruent to a scalar in O^*. \\ ** ggsd: the square root of d \\ InitForms(d,\ temp,structure,generators,j,u,v)=\ ggD = d;\ ggw = quadgen(d);\ ggi = 1;\ ggunit = unit(d)^ggi;\ u = real(ggunit); v= imag(ggunit);\ if( (ggD % 4) ==1, u = u+ v/2; v=v/2, v=v/2);\ ggMunit = [ u, ggd * v;\ v, u ];\ temp=buchreal(d,1,0.5,0.5,20);\ ggh = temp[1];\ if(norm(ggunit)==1,ggh=ggh*2,);\ structure=temp[2];\ generators=length(structure);\ if(generators==0,gggen=pf(ggd,2),gggen=temp[3][1]);\ gggen=reduce(stripd(gggen));\ ggclassgroup=temp[3];\ if((d%4)==0,ggsd=2*w,ggsd=2*w-1);\ for(j=1,generators,ggclassgroup[j]=stripd(ggclassgroup[j])); \\ stripd: A utility routine to convert a real quadratic form \\ qfr(a,b,c,d) to the usual representation [a,b,c] stripd(form)=\ [compo(form,1),compo(form,2),compo(form,3)]; \\ convertqfr: Performs the converse conversion,i.e,. turns a \\ quadratic form of the type [a,b,c] into qfr(a,b,c,0). convertqfr(form)=\ qfr(form[1],form[2],form[3],0.0); \\ reduces the form f with positive discriminant. reduce(f)=\ stripd(redreal(convertqfr(f)) ); \\ next(f): \\ returns the form adjacent to f in the cycle of reduced forms. next(f)=\ stripd(rhoreal(convertqfr(f)) ); \\ eqf(f1,f2): \\ Returns 1 if f1 is equivalent to f2 in the class group, 0 otherwise. eqf(f1,f2,res,ff)=\ f1=reduce(f1);\ ff=f2;\ until(ff==f1 || ff==f2, ff=next(ff));\ if(ff==f1,res=1,res=0);\ res \\ The folllowing routine implements the Guassian composition of \\ binary quadrartic forms. compose(f1,f2)=\ stripd(redreal(comprealraw( convertqfr(f1),convertqfr(f2) ))); \\ pwr(f1,n): \\ Returns the form f1 to the power n, in the class group. pwr(f1,n)=\ stripd(redreal(powrealraw(convertqfr(f1),n))); \\ utility routine to find the exact order of f in the class group. \\ Needs the form f, and a bound n on the order of f. Tries all divisors of \\ n. exactorder(f,n,id,d,ff,res)=\ res=0;\ fordiv(n,d,if(res==0,ff=pwr(f,d);\ print1(d);print1(" ");res=eqf(gid,ff);print(res);\ if(res==0,,res=d),\ ));\ res \\ Given forms f1 and f2 of exact order n, returns a value k such that \\ f1=k f2. If no such exists, returns k=0. multipl(f1,f2,n,ff,res)=\ ff=f2;\ res=0;\ for(k=1,n,if(res==0,res=eqf(ff,f1)*k;ff=cmp(ff,f2),));\ res \\ List manipulation routines: last extracts the last element from a \\ list v, and rest returns the rest of the list. last(v)=compo(v,length(v)); rest(v,n,res)=\ res=[];\ n=length(v)-1;\ for(j=1,n,res=concat(res,[compo(v,j)]));\ res \\ Generate some quadratic forms. \\ Input is a prime p (first coeff) and j. try(p)= reduce(stripd(pf(ggd,p))) \\ ************************************************************************ \\ ** Routines for computing the geodesic cycles and the local points. \\ ************************************************************************ \\ \\ \\ IPsi(F): takes as input a binary quadratic form F=(a,b,c) of discriminant \\ D, and returns the period attached to the embdding Psi sending \\ sqrt(D) to the matrix | b 2c |, and hence \\ | -2a -b |. \\ w to | (b+1)/2 c | \\ | -a (-b+1)/2 | IPsi(F,\ a,b,c,matw,gammaPsi,zPsi) =\ a=F[1]; b=F[2]; c=F[3];\ d= (b^2-4*a*c);\ if( (d%4)==0, matw=[b,2*c;-2*a,-b],\ matw = [(b+1)/2,c; -a, (-b+1)/2 ]);\ gammaPsi = real(ggUnit)*[1,0;0,1] + imag(ggunit)*matw;\ zPsi = (- b + ggsd )/(2*a);\ dmint(zPsi, gammaPsi[1,1]/gammaPsi[2,1]));\ \\ The function dmint(z,q) implements the double integral introduced in the \\ paper, in an important special case: that where z is defined over the quadratic \\ unramified extension of Q_p, and is of the form a+bw, where a is a rational \\ number in Zp and b is in Zp*. \\ The integral computed is dmint(z-bar,z; infinity,q). dmint(z,q,\ res,p,n,u,alpha,lhs,rhs,Nalpha,M,v) =\ res=1;\ p=ggp; n=ggn;\ for(u=0,p^n-1,\ alpha = u+w;\ lhs=alpha-conj(alpha);\ rhs = (alpha*z)-conj(alpha*z);\ Nalpha = (rhs/lhs) % (p^n);\ M = [1,-Nalpha; 0,p^n];\ res=(res* (alpha/conj(alpha))^pathinfty( (q-Nalpha)/p^n) ) % p^n);\ for(v=0,p^(n-1)-1,\ alpha = 1+p*v*w;\ lhs = (alpha-conj(alpha));\ rhs = (alpha*z)-conj(alpha*z);\ Nalpha = ((-lhs/rhs) % (p^n));\ res=(res*(alpha/conj(alpha))^path(-Nalpha/p^n,(-Nalpha*q-1)/(p^n*q))) % (p^n) );\ res)