\p 100 \\ Global variable D: -D \\ is the discriminant of an imaginary quadratic field K. D = 4; \\ The elliptic curve E is the CM curve of conductor \\ 32, with CM by K. Here D=4. E = ellinit([0,0,0,-1,0]); \\ These are the periods of the elliptic curve Omegap = E[15]; Omegam = E[16]; \\ The function teta(M) returns the theta series of weight 5 \\ attached to the discriminant -4. This is a form on \\ Gamma_1(4), with odd character. \\ It is stored as a vector v of length M^2, \\ where v[n] is the nth \\ fourier coefficient of the theta series. \\ teta(M,\ res, q, a, b, nrm, cf ) = \ res = vector(M^2);\ a = 1;\ while(a<=M,\ b = 0 ;\ while((a^2+b^2)<=M^2,\ nrm = (a^2+b^2);\ cf = (a^4 - 6*a^2*b^2 + b^4);\ res[nrm] = res[nrm] + cf;\ b = b+1);\ a = a+1);\ return(res); f = teta(100); \\ The function evl(f,tau) evaluates the form f at the argument \\ tau in the upper half-plane. \\ It was mainly used as a debugging tool to check that the theta \\ series computed by teta satisfies the required invariance \\ properties. evl(f,tau,\ q, res, M, j) = \ q = exp(2*Pi*I*tau);\ res=0;\ M = length(f);\ for(j=1,M,res=res+f[j]*q^j);\ return(res); \\ prim0(f,tau): Returns the integral from (i.infinity) to tau \\ of f(z) dz, for f a modular form of weight 5. \\ prim0(f,tau,\ res,M,n,q) = \ res=0;\ q = exp(2*Pi*I*tau);\ M = length(f);\ for(n=1,M, res=res+f[n]/n*q^n);\ res=res/(2*Pi*I);\ return(res); \\ prim1(f,tau): Returns the integral from (i.infinity) to tau \\ of z.f(z)dz, for f a modular form of weight 3. \\ prim1(f,tau,\ res,M,n,q) = \ res=0;\ q = exp(2*Pi*I*tau);\ M = length(f);\ for(n=1,M, res=res-f[n]/n^2*q^n);\ res=res/(2*Pi*I)^2;\ res=res+tau*prim0(f,tau);\ return(res); \\ per0(f,g): Returns the period integral of f(z)dz from \\ (i.infty) to g(i.infty), \\ Here g is taken to be a matrix in Gamma_0(D) with a(g) a square \\ modulo D. \\ per0(f,g,\ a,b,c,d,res) = \ a = g[1,1]; b=g[1,2];\ c = g[2,1]; d=g[2,2];\ res = prim0(f,a/c+ I/c);\ res = res - c*prim1(f,-d/c+I/c);\ res = res - d*prim0(f,-d/c+I/c);\ return(4*Pi^2*res); \\ per1(f,g): Returns the period integral of z.f(z)dz from \\ (i.infty) to g(i.infty), \\ Here g is taken to be a matrix in Gamma_0(D) with a(g) a square \\ modulo D. \\ per1(f,g,\ a,b,c,d,res) = \ a = g[1,1]; b=g[1,2];\ c = g[2,1]; d=g[2,2];\ res = prim1(f,a/c+ I/c);\ res = res - a*prim1(f,-d/c+I/c);\ res = res - b*prim0(f,-d/c+I/c);\ return(4*Pi^2*res); \\ The function shp returns the Stark-Heegner point \\ asssociated to the lattice generated by \\ omega1 and omega2. This lattice is assumed to be \\ positively oriented. \\ We expect the resulting invariant to be a global \\ point if [omega1,omega2] is the Neron lattice \\ of an elliptic curve with CM by O_K. shp(omega1,omega2,\ tau, vol,res) = \ tau = omega1/omega2;\ vol = (omega1*conj(omega2)-conj(omega1)*omega2);\ res = vol^(-1)*(conj(omega2)*prim1(f,tau) - conj(omega1)*prim0(f,tau));\ return(4*Pi^2*res); \\ [omega1, omega2] generate the Neron lattice of E. omega1 = Omegap*(D+sqrt(-D))/(2*D);\ omega2 = Omegap;\ \\ Compute the Stark-Heegner point J =3* shp( omega1, omega2 ); J = (J+2*Omegam)/sqrt(-D); print("Stark-Heegner point of conductor 1"); P = ellztopoint(E,J)