<< GInt3.m << thetas.m (**************************************************************************** ToTorus[CurveData,{Px,Py},pres] maps the point {Px,Py} to the analytic Jacobian, where Px and Py are the x and y coordinates of a point on y^2 = pol. pres is the required precision. CurveData should equal InitCurve[pol,pres], with pres at least as big as the precision used here. The analytic Jacobian is understood to be C^2/(1,tau). Options are those for gGInt. ****************************************************************************) ToTorus[CurveData_,{Px_,Py_},pres_,opts___] := Module[{px,py}, If[CurveData[[2]] === Infinity, ToAnnJg[{Px,Py/Sqrt[CurveData[[4]]]},CurveData[[5]],CurveData[[1,1]],CurveData[[1,2]],pres,opts], px = 1/(Px-CurveData[[2]]); py = Py*px^3/Sqrt[CurveData[[3]]*CurveData[[4]]]; ToAnnJg[{px,py},CurveData[[5]],CurveData[[1,1]],CurveData[[1,2]],pres,opts]]] (**************************************************************************** FromTorus[CurveData,z,pres] maps the point z on the analytic Jacobian to algebraic Jacobian. The analytic Jacobian is understood as C^2/(1,tau) and the algebraic Jacobian is represented by the symmetric product of the curve with itself. The function therefore returns two points on the curve y^2 = pol in the form {{P1x,P1y},{P2x,P2y}}. pres is the required precision. CurveData should equal InitCurve[pol,pres], with pres at least as big as the precision used here. Options are those for gGInt and theta. ****************************************************************************) FromTorus[CurveData_,z_,pres_,opts___] := Module[{dum}, dum = ToCF[z,tauFromPerM[CurveData[[1,2]]]]; If[Abs[dum[[1,1]]]+Abs[dum[[2,1]]] < 10^(-pres*8/10),Print["Looks like z is on the lattice!"];{{Null,Null},{Null,Null}}, dum = xyFromz[dum,CurveData[[1,1]],CurveData[[5]],CurveData[[1,2]],pres,opts]; If[dum === {{Null,Null},{Null,Null}},dum, {ToCurve1[CurveData,dum[[1]]],ToCurve1[CurveData,dum[[2]]]}]]] ToCurve1[CurveData_,{px_,py_}] := If[CurveData[[2]] === Infinity, {px,py*Sqrt[CurveData[[4]]]}, {1/px + CurveData[[2]],py*Sqrt[CurveData[[3]]*CurveData[[4]]]/px^3}] (* This table gives the two division point corresponding to each of the five roots *) HlfData = {{0, 0, 1, 0},{1, 0, 1, 0},{1, 0, 0, 1},{1, 1, 0, 1},{1, 1, 0, 0}}; (* This maps the point {Px,Py} to the analytic Jacobian *) ToAnnJg[{Px_,Py_},InfInts_,rts_,PerM_,pres_:16,opts___] := Module[{sgn,CP,Px2=Px,w2i,g,dum,ops = opts}, g = Length[PerM]; If[Abs[Py] < 10^(-pres*9/10), dum = Map[Abs[#-Px]&,rts]; Transpose[{HlfData[[Position[dum,Min[dum]][[1,1]]]]}.Join[{{1,0},{0,1}},Transpose[tauFromPerM[PerM]]]]/2, If[Abs[Py] < 10^(-15) && (IntegrationMethod /. {opts} /. Options[gGIntP]) == Pari, Print["You are very close to a Weierstrass point."]; Print["We will have to do this slowly (no Pari)."]; ops = IntegrationMethod -> Gauss]; sgn = N[Product[Sqrt[Px-rts[[i]]],{i,1,Length[rts]}]/Py]; If[Abs[sgn-Sign[Re[sgn]]] > 10^(-10),Print["The point is not on the curve!"];Abort[]]; sgn = Sign[Re[sgn]]; CP = gGIntP[rts,Path2Inf[Px,rts,InfInts,pres],pres,ops]+InfInts[[2]]; CP = Transpose[{CP}]; w2i = Inverse[Table[Take[PerM[[i]],-g],{i,1,g}]]; ToCF[sgn*w2i.CP,w2i.Table[Take[PerM[[i]],g],{i,1,g}]]]]; Path2Inf[Px_,rts_,InfInts_,pres_] := Module[{dum,tol,thei,they}, dum = Sort[Map[Im,rts]]; tol = Min[Select[Table[dum[[i]]-dum[[i-1]],{i,2,Length[dum]}],# > 10^(-pres*8/10) &],1]/5; dum = Sort[Map[Im,Select[rts,Re[#] > Re[Px]-tol &]]]; Do[If[Abs[dum[[i-1]]-dum[[i]]] < 10^(-pres*8/10),dum = Drop[dum,{i}]],{i,Length[dum],2,-1}]; If[Min[Map[Abs,dum - Im[Px]]] < tol, thei = 0; Do[ If[thei == 0 && Im[Px] < dum[[i]],thei = i], {i,1,Length[dum]}]; If[thei == 0, they = Im[Px]+tol*10, If[thei == 1, they = Im[Px]-tol*10, they = (dum[[thei]]+dum[[thei-1]])/2]]; {Px,Re[Px]+I*they,InfInts[[1]]+I*they,InfInts[[1]]}, If[Abs[Im[Px]] < 10^(-Precision[Px]+10), {Px,InfInts[[1]]}, {Px,InfInts[[1]]+Im[Px]*I,InfInts[[1]]}]]] tauFromPerM[PerM_] := Module[{g}, g = Length[PerM]; Inverse[Table[Take[PerM[[i]],-g],{i,1,g}]].Table[Take[PerM[[i]],g],{i,1,g}]] (* This takes a point on the Jacobian into the usual Fundamental domain *) ToCF[CP_,tau_,tol_:10^(-14)] := Module[{v1,v2}, v1 = MyFloor[Inverse[Im[tau]].Im[CP],tol]; v2 = MyFloor[Re[CP-tau.v1],tol]; CP - v2 - tau.v1]; Attributes[MyFloor] = {Listable}; MyFloor[x_,tol_] := Floor[x+tol]; (* This takes a point in the analytic Jacobian back to the algebraic Jacobian. It is assumed that nz has been ToCF'ed and that nz is not on the lattice. *) xyFromz[nz_,rts_List,InfInts_,PerM_,pres_,opts___] := Module[{nx1,nx2,ny1,ny2,zlst,tau,Stop}, tau = tauFromPerM[PerM]; {nx1,nx2} = xFromz[nz,rts,tau,pres,opts]; If[(Precision[nx1] < 5) || (Precision[nx2] < 5),{{Null,Null},{Null,Null}}, Print["Precision of nx1 and nx2: "]; Print[Precision[nx1]," ",Precision[nx2]]; ny1 = yFromx[nx1,rts,pres][[1]]; ny2 = yFromx[nx2,rts,pres][[1]]; zlst = allzs[{ToAnnJg[{nx1,ny1},InfInts,rts,PerM,16,opts],ToAnnJg[{nx2,ny2},InfInts,rts,PerM,16,opts]},tau]; zlst = Map[Abs[#[[1,1]]-nz[[1,1]]]^2+Abs[#[[2,1]]-nz[[2,1]]]^2 &,zlst]; If[Abs[Min[zlst]] > 10^(-30),Print["*************WARNING***************"];Print["Something is almost certainly wrong, can't find correct y-coordinates"]]; {{{nx1,ny1},{nx2,ny2}},{{nx1,ny1},{nx2,-ny2}}, {{nx1,-ny1},{nx2,ny2}},{{nx1,-ny1},{nx2,-ny2}}}[[Position[zlst,Min[zlst]][[1,1]]]]]] (* Given a point z in the analytic Jacobian corresponding to tau, this finds the two x-coordinates of the point in the algebraic Jacobian corresponding to z *) xFromz[z_,rts_List,tau_,pres_:16,opts___] := xFromU[rts[[1]],Uatak[1,{2,3},z,rts,tau,pres,opts], rts[[2]],Uatak[2,{1,3},z,rts,tau,pres,opts],pres] xFromU[a1_,Uata1_,a2_,Uata2_,pres_] := Module[{disc,b}, If[Abs[Uata1] > 10^(pres*1.8) || Abs[Uata2] > 10^(pres*1.8), Print["You seem to be on the theta divisor. Try to compute twice your point"]; {SetPrecision[1.0,3],1}, b = Uata2 - Uata1 + a1^2 - a2^2; disc=Sqrt[b^2-4*(a2-a1)*(a2*Uata1-a1*Uata2+a1*a2*(a2-a1))]; If[Abs[disc] < 10^(-Round[0.9*Precision[{a1,Uata1,a2,Uata2}]]),disc = 0]; {(-b+disc)/(2*(a2-a1)),(-b-disc)/(2*(a2-a1))}]] Uatak[k_,V_,z_,rts_List,tau_,pres_:16,opts___] := Module[{UV = SymDiv[Append[V,k],UbranchPts]}, (-1)^(4*dp.eta[k][[2]])* (-1)^(4*eta[UV][[2]].eta[k][[1]])* N[(rts[[k]]-rts[[V[[1]]]])*(rts[[k]]-rts[[V[[2]]]]),pres]* theta[eta[UV]+eta[k],{{0},{0}},tau,pres,opts]^2* theta[del+eta[k],z,tau,pres,opts]^2/ (theta[eta[UV],{{0},{0}},tau,pres,opts]^2*theta[del,z,tau,pres,opts]^2)] eta[1] = {{1/2,0},{0,0}}; eta[2] = {{1/2,0},{1/2,0}}; eta[3] = {{0,1/2},{1/2,0}}; eta[4] = {{0,1/2},{1/2,1/2}}; eta[5] = {{0,0},{1/2,1/2}}; eta[list_] := If[list == {}, {{0,0},{0,0}}, Sum[eta[list[[i]]],{i,1,Length[list]}]]; UbranchPts = {1,3,5}; BranchPts = {1,2,3,4,5}; SymDiv[list1_,list2_] := Complement[Union[list1,list2],Intersection[list1,list2]]; del = {{1/2,1/2},{1,1/2}}; dp = {1/2,1/2}; dpp = {1,1/2}; yFromx[x_,rts_List,pres_:16] := Module[{dum}, {dum = Product[N[Sqrt[x-rts[[i]]],pres],{i,1,Length[rts]}],-dum}]; allzs[{z1_,z2_},tau_] := {ToCF[z1+z2,tau],ToCF[z1-z2,tau],ToCF[-z1+z2,tau],ToCF[-z1-z2,tau]}