\\ \\ reciprocal3.gp \\ \\ Last Update: 2004.07.12 \\ \\ find (x,y,z) such that (x+y+z)(1/x+1/y+1/z)=n \\ ec(n)=ellinit([0,n^2-6*n-3,0,16*n,0]); \\ \\ E_n: Y^2Z = X^3+(n^2-6n-3)X^2Z+16nXZ^2 \\ \\ C_n: (x+y+1)(1/x+1/y+1)=n \\ \\ f:(X,Y,Z) ----> (x,y,z) \\ x=(n-1)*X-Y, \\ y=(n-1)*X+Y, \\ z=2(X-4nZ) \\ \\ g:(x,y,z) ----> (X,Y,Z) \\ X=4((x+y)z-(n-1)xy)z, \\ Y=4(n-1)(x-y)((x+y)z+xy), \\ Z=(x+y)z^2 ee(x,y,z,n)=x^3+(n^2-6*n-3)*x^2*z+16*n*x*z^2-y^2*z; cc(x,y,z,n)=(x+y+z)*(1/x+1/y+1/z)-n; \\ phi: C_n(Q)--> E_n(Q) gg(p)= { local(x,y,z,n,X,Y,Z,d); x=p[1]; y=p[2]; z=p[3]; n=(x+y+z)*(1/x+1/y+1/z); \\ print("n=",n); X=4*((x+y)*z-(n-1)*x*y)*z; Y=4*(n-1)*(x-y)*((x+y)*z+x*y); Z=(x+y)*z^2; d=gcd(gcd(X,Y),Z); if(Z<0, d=-d); X=X/d; Y=Y/d; Z=Z/d; [X,Y,Z] } \\ \\ [X:Y:Z] on E_n --> n \\ where n > 0 \\ nn(x,y,z)= { local(a,b2,c,d,d2,n); a=z*x^2; b2=-3*z*x^2 + 8*z^2*x; c=x^3 - 3*z*x^2 - y^2*z; \\ print("[a,b,c]=",[a,2*b2,c]); d=b2^2-a*c; d2=sqrtrational(d); n=nil; if(a==0 && b2!=0, n=-c/(2*b2) ); if(a!=0 && d2!=nil, if(a>0, n=(-b2+d2)/a, n=(-b2-d2)/a ) ); \\ print("n=",n); n } \\ \\ [X:Y:Z] on E_n --> n \\ where n < 0 \\ nneg(x,y,z)= { local(a,b2,c,d,d2,n); a=z*x^2; b2=-3*z*x^2 + 8*z^2*x; c=x^3 - 3*z*x^2 - y^2*z; \\ print("[a,b,c]=",[a,2*b2,c]); d=b2^2-a*c; \\ print("d=",d); d2=sqrtrational(d); \\ print("d2=",d2); n=nil; if(a==0 && b2!=0, n=-c/(2*b2) ); if(a!=0 && d2!=nil, if(a>0, n=(-b2-d2)/a, n=(-b2+d2)/a ) ); \\ print("n=",n); n } \\ \\ phi^(-1): E_n(Q)-->C_n(Q) \\ p[x,y] or p[x,y,z] \\ ff(p,n)= { local(x,y,z,X,Y,Z,d); \\ print("ff:p=",p); if(p==[0],return(nil)); X=p[1]; Y=p[2]; if(length(p)==3, Z=p[3], Z=1 ); if(Z==1, d=lcm(denominator(X),denominator(Y))/gcd(numerator(X),numerator(Y)); Z=Z*d; X=X*d; Y=Y*d ); \\ print("[X,Y,Z]=",[X,Y,Z]); \\ n=nn(X,Y,Z); \\ print("n=",n); x=(n-1)*X-Y; y=(n-1)*X+Y; z=2*(X-4*n*Z); d=gcd(gcd(x,y),z); if(z<0, d=-d); x=x/d; y=y/d; z=z/d; [x,y,z] } \\ \\ [X:Y:Z]-->[X':Y':Z'] \\ \\ X'=min(X,Y,Z); Y'=min(max(X,Y),Z) ; Z'=max(X,Y,Z); \\ |X'| <= |Y'| <= |Z'| \\ so(p)= { local(x,y,z,w); x=p[1]; y=p[2]; z=p[3]; \\ x<=y if(abs(x)>=abs(y), w=x; x=y; y=w ); \\ y<=z if(abs(y)>=abs(z), w=y; y=z; z=w ); \\ x<=y if(abs(x)>=abs(y), w=x; x=y; y=w ); \\ z>0 if(z<=0, z=-z; x=-x; y=-y ); [x,y,z] } \\ ffso(p,n)= { local(q); q=ff(p,n); if(q!=nil, so(q), nil ) } \\ ffsoh(p,n,hlimit)= { local(q); q=ff(p,n); if(q!=nil && hh(q) < hlimit, so(q), nil ) } \\ \\ p1: generator of E_n(Q) where n > 0 \\ f1(p1,m)= { local(en,x,y,z,p,q,r,n,tn); if(length(p1)==3, z=p1[3], z=1 ); x=p1[1]/z; y=p1[2]/z; n=nn(x,y,1); if(n==nil,return); print("n=",n); en=ec(n); tn=[4*n,4*n*(n-1)]; p=[x,y]; for(i=1,m, q=ellpow(en,p,i); \\ print1("[",i,"]p1="); print(ffso(q,n)); r=elladd(en,q,tn); \\ print1("[",i,"]p1+tn="); print(ffso(r,n)) ) } \\ \\ p1,p2: genarator of E_n(Q) where n > 0 \\ f2(p1,p2,m)= { local(en,x1,y1,z1,x2,y2,z2,p,q1,q2,q,r,s,n,n2,tn); if(length(p1)==3, z1=p1[3], z1=1 ); x1=p1[1]/z1; y1=p1[2]/z1; if(length(p2)==3, z2=p2[3], z2=1 ); x2=p2[1]/z2; y2=p2[2]/z2; n=nn(x1,y1,1); n2=nn(x2,y2,1); if(n!=n2,return); print("n=",n); en=ec(n); tn=[4*n,4*n*(n-1)]; q1=[x1,y1]; q2=[x2,y2]; for(i=0,m, p=ellpow(en,q1,i); for(j=0,m, q=ellpow(en,q2,j); r=elladd(en,p,q); if(i!=0 || j!=0, \\ print1("[",i,"]p1+[",j,"]p2="); print(ffso(r,n)); s=elladd(en,r,tn); \\ print1("[",i,"]p1+[",j,"]p2+tn="); print(ffso(s,n)) ); if(j!=0, q=ellpow(en,q,-1); r=elladd(en,p,q); if(i!=0 || j!=0, \\ print1("[",i,"]p1+[",-j,"]p2="); print(ffso(r,n)); s=elladd(en,r,tn); \\ print1("[",i,"]p1+[",-j,"]p2+tn="); print(ffso(s,n)) ) ) ) ) } \\ \\ p1,p2: genarator of E_n(Q) where n > 0 \\ f2l(p1,p2,m,hlimit)= { local(en,x1,y1,z1,x2,y2,z2,p,q1,q2,q,r,s,n,n2,tn); if(length(p1)==3, z1=p1[3], z1=1 ); x1=p1[1]/z1; y1=p1[2]/z1; if(length(p2)==3, z2=p2[3], z2=1 ); x2=p2[1]/z2; y2=p2[2]/z2; n=nn(x1,y1,1); n2=nn(x2,y2,1); if(n!=n2,return); print("n=",n); en=ec(n); tn=[4*n,4*n*(n-1)]; q1=[x1,y1]; q2=[x2,y2]; for(i=0,m, p=ellpow(en,q1,i); for(j=0,m, q=ellpow(en,q2,j); r=elladd(en,p,q); if(i!=0 || j!=0, \\ print1("[",i,"]p1+[",j,"]p2="); print_nonnil(ffsoh(r,n,hlimit)); s=elladd(en,r,tn); \\ print1("[",i,"]p1+[",j,"]p2+tn="); print_nonnil(ffsoh(s,n,hlimit)) ); if(j!=0, q=ellpow(en,q,-1); r=elladd(en,p,q); if(i!=0 || j!=0, \\ print1("[",i,"]p1+[",-j,"]p2="); print_nonnil(ffsoh(r,n,hlimit)); s=elladd(en,r,tn); \\ print1("[",i,"]p1+[",-j,"]p2+tn="); print_nonnil(ffsoh(s,n,hlimit)) ) ) ) ) } \\ \\ p1,p2,p3: genarator of E_n(Q) where n > 0 \\ f3(p1,p2,p3,m)= { local(en,x1,y1,z1,x2,y2,z2,x3,y3,z3,p,q1,q2,q3,q,r,s,u,v,n,n2,n3,tn); if(length(p1)==3, z1=p1[3], z1=1 ); x1=p1[1]/z1; y1=p1[2]/z1; if(length(p2)==3, z2=p2[3], z2=1 ); x2=p2[1]/z2; y2=p2[2]/z2; if(length(p3)==3, z3=p3[3], z3=1 ); x3=p3[1]/z3; y3=p3[2]/z3; n=nn(x1,y1,1); n2=nn(x2,y2,1); if(n!=n2,return); n3=nn(x3,y3,1); if(n!=n3,return); print("n=",n); en=ec(n); tn=[4*n,4*n*(n-1)]; q1=[x1,y1]; q2=[x2,y2]; q3=[x3,y3]; for(i=0,m, p=ellpow(en,q1,i); for(j=0,m, q=ellpow(en,q2,j); r=elladd(en,p,q); for(k=-m,m, s=ellpow(en,q3,k); u=elladd(en,r,s); if(i!=0 || j!=0 || k!=0, \\ print1("[",i,"]p1+[",j,"]p2+[",k,"]p3="); print(ffso(s,n)) ); v=elladd(en,u,tn); \\ print1("[",i,"]p1+[",j,"]p2+[",k,"]p3+tn="); print(ffso(v,n)) ); if(j!=0, q=ellpow(en,q,-1); r=elladd(en,p,q); for(k=-m,m, s=ellpow(en,q3,k); u=elladd(en,r,s); if(i!=0 || j!=0 || k!=0, \\ print1("[",i,"]p1+[",-j,"]p2+[",k,"]p3="); print(ffso(u,n)) ); v=elladd(en,u,tn); \\ print1("[",i,"]p1+[",-j,"]p2+[",k,"]p3+tn="); print(ffso(v,n)) ) ) ) ) } \\ \\ p1,p2,p3: genarator of E_n(Q) where n > 0 \\ f3l(p1,p2,p3,m,hlimit)= { local(en,x1,y1,z1,x2,y2,z2,x3,y3,z3,p,q1,q2,q3,q,r,s,u,v,n,n2,n3,tn); if(length(p1)==3, z1=p1[3], z1=1 ); x1=p1[1]/z1; y1=p1[2]/z1; if(length(p2)==3, z2=p2[3], z2=1 ); x2=p2[1]/z2; y2=p2[2]/z2; if(length(p3)==3, z3=p3[3], z3=1 ); x3=p3[1]/z3; y3=p3[2]/z3; n=nn(x1,y1,1); n2=nn(x2,y2,1); if(n!=n2,return); n3=nn(x3,y3,1); if(n!=n3,return); print("n=",n); en=ec(n); tn=[4*n,4*n*(n-1)]; q1=[x1,y1]; q2=[x2,y2]; q3=[x3,y3]; for(i=0,m, p=ellpow(en,q1,i); for(j=0,m, q=ellpow(en,q2,j); r=elladd(en,p,q); for(k=-m,m, s=ellpow(en,q3,k); u=elladd(en,r,s); if(i!=0 || j!=0 || k!=0, \\ print1("[",i,"]p1+[",j,"]p2+[",k,"]p3="); print_nonnil(ffsoh(s,n,hlimit)) ); v=elladd(en,u,tn); \\ print1("[",i,"]p1+[",j,"]p2+[",k,"]p3+tn="); print_nonnil(ffsoh(v,n,hlimit)) ); if(j!=0, q=ellpow(en,q,-1); r=elladd(en,p,q); for(k=-m,m, s=ellpow(en,q3,k); u=elladd(en,r,s); if(i!=0 || j!=0 || k!=0, \\ print1("[",i,"]p1+[",-j,"]p2+[",k,"]p3="); print_nonnil(ffsoh(u,n,hlimit)) ); v=elladd(en,u,tn); \\ print1("[",i,"]p1+[",-j,"]p2+[",k,"]p3+tn="); print_nonnil(ffsoh(v,n,hlimit)) ) ) ) ) } \\ \\ print_nonnil(p) \\ print p if p != nil \\ print_nonnil(p)= { local(x,y,z); if(p!=nil, x=p[1]; y=p[2]; z=p[3]; print("[", x, " : ", y ," : ", z, "]") ); } \\ \\ p1: generator of E_n(Q) where n < 0 \\ f1n(p1,m)= { local(en,x,y,z,p,q,r,n,tn); if(length(p1)==3, z=p1[3], z=1 ); x=p1[1]/z; y=p1[2]/z; n=nneg(x,y,1); if(n==nil,return); print("n=",n); en=ec(n); tn=[4*n,4*n*(n-1)]; p=[x,y]; for(i=1,m, q=ellpow(en,p,i); \\ print1("[",i,"]p1="); print_nonnil(ffso(q,n)); r=elladd(en,q,tn); \\ print1("[",i,"]p1+tn="); print_nonnil(ffso(r,n)) ) } \\ \\ p1,p2: genrator of E_n(Q) where n < 0 \\ f2n(p1,p2,m)= { local(en,x1,y1,z1,x2,y2,z2,p,q1,q2,q,r,s,n,n2,tn); if(length(p1)==3, z1=p1[3], z1=1 ); x1=p1[1]/z1; y1=p1[2]/z1; if(length(p2)==3, z2=p2[3], z2=1 ); x2=p2[1]/z2; y2=p2[2]/z2; n=nneg(x1,y1,1); n2=nneg(x2,y2,1); if(n!=n2,return); print("n=",n); en=ec(n); tn=[4*n,4*n*(n-1)]; q1=[x1,y1]; q2=[x2,y2]; for(i=0,m, p=ellpow(en,q1,i); for(j=0,m, q=ellpow(en,q2,j); r=elladd(en,p,q); if(i!=0 || j!=0, \\ print1("[",i,"]p1+[",j,"]p2="); print_nonnil(ffso(r,n)); s=elladd(en,r,tn); \\ print1("[",i,"]p1+[",j,"]p2+tn="); print_nonnil(ffso(s,n)) ); if(j!=0, q=ellpow(en,q,-1); r=elladd(en,p,q); if(i!=0 || j!=0, \\ print1("[",i,"]p1+[",-j,"]p2="); print_nonnil(ffso(r,n)); s=elladd(en,r,tn); \\ print1("[",i,"]p1+[",-j,"]p2+tn="); print_nonnil(ffso(s,n)) ) ) ) ) } \\ \\ p1,p2: genrator of E_n(Q) where n < 0 \\ f2nl(p1,p2,m,hlimit)= { local(en,x1,y1,z1,x2,y2,z2,p,q1,q2,q,r,s,n,n2,tn); if(length(p1)==3, z1=p1[3], z1=1 ); x1=p1[1]/z1; y1=p1[2]/z1; if(length(p2)==3, z2=p2[3], z2=1 ); x2=p2[1]/z2; y2=p2[2]/z2; n=nneg(x1,y1,1); n2=nneg(x2,y2,1); if(n!=n2,return); print("n=",n); en=ec(n); tn=[4*n,4*n*(n-1)]; q1=[x1,y1]; q2=[x2,y2]; for(i=0,m, p=ellpow(en,q1,i); for(j=0,m, q=ellpow(en,q2,j); r=elladd(en,p,q); if(i!=0 || j!=0, \\ print1("[",i,"]p1+[",j,"]p2="); print_nonnil(ffsoh(r,n,hlimit)); s=elladd(en,r,tn); \\ print1("[",i,"]p1+[",j,"]p2+tn="); print_nonnil(ffsoh(s,n,hlimit)) ); if(j!=0, q=ellpow(en,q,-1); r=elladd(en,p,q); if(i!=0 || j!=0, \\ print1("[",i,"]p1+[",-j,"]p2="); print_nonnil(ffsoh(r,n,hlimit)); s=elladd(en,r,tn); \\ print1("[",i,"]p1+[",-j,"]p2+tn="); print_nonnil(ffsoh(s,n,hlimit)) ) ) ) ) } \\ \\ p1,p2,p3: genarator of E_n(Q) where n < 0 \\ f3n(p1,p2,p3,m)= { local(en,x1,y1,z1,x2,y2,z2,x3,y3,z3,p,q1,q2,q3,q,r,s,u,v,n,n2,n3,tn); if(length(p1)==3, z1=p1[3], z1=1 ); x1=p1[1]/z1; y1=p1[2]/z1; if(length(p2)==3, z2=p2[3], z2=1 ); x2=p2[1]/z2; y2=p2[2]/z2; if(length(p3)==3, z3=p3[3], z3=1 ); x3=p3[1]/z3; y3=p3[2]/z3; n=nneg(x1,y1,1); n2=nneg(x2,y2,1); if(n!=n2,return); n3=nneg(x3,y3,1); if(n!=n3,return); print("n=",n); en=ec(n); tn=[4*n,4*n*(n-1)]; q1=[x1,y1]; q2=[x2,y2]; q3=[x3,y3]; for(i=0,m, p=ellpow(en,q1,i); for(j=-m,m, q=ellpow(en,q2,j); r=elladd(en,p,q); for(k=-m,m, s=ellpow(en,q3,k); u=elladd(en,r,s); if(i!=0 || j!=0 || k!=0, \\ print1("[",i,"]p1+[",j,"]p2+[",k,"]p3="); print_nonnil(ffso(s,n)) ); v=elladd(en,u,tn); \\ print1("[",i,"]p1+[",j,"]p2+[",k,"]p3+tn="); print_nonnil(ffso(v,n)) ) ) ) } \\ \\ p1,p2,p3: genarator of E_n(Q) where n < 0 \\ f3nl(p1,p2,p3,m,hlimit)= { local(en,x1,y1,z1,x2,y2,z2,x3,y3,z3,p,q1,q2,q3,q,r,s,u,v,n,n2,n3,tn); if(length(p1)==3, z1=p1[3], z1=1 ); x1=p1[1]/z1; y1=p1[2]/z1; if(length(p2)==3, z2=p2[3], z2=1 ); x2=p2[1]/z2; y2=p2[2]/z2; if(length(p3)==3, z3=p3[3], z3=1 ); x3=p3[1]/z3; y3=p3[2]/z3; n=nneg(x1,y1,1); n2=nneg(x2,y2,1); if(n!=n2,return); n3=nneg(x3,y3,1); if(n!=n3,return); print("n=",n); en=ec(n); tn=[4*n,4*n*(n-1)]; q1=[x1,y1]; q2=[x2,y2]; q3=[x3,y3]; for(i=0,m, p=ellpow(en,q1,i); for(j=-m,m, q=ellpow(en,q2,j); r=elladd(en,p,q); for(k=-m,m, s=ellpow(en,q3,k); u=elladd(en,r,s); if(i!=0 || j!=0 || k!=0, \\ print1("[",i,"]p1+[",j,"]p2+[",k,"]p3="); print_nonnil(ffsoh(s,n,hlimit)) ); v=elladd(en,u,tn); \\ print1("[",i,"]p1+[",j,"]p2+[",k,"]p3+tn="); print_nonnil(ffsoh(v,n,hlimit)) ) ) ) } torsion(p)= { local(q,r,n); n=nn(p[1],p[2],1); if(n!=nil, en=ec(n); for(j=0,5, q=ellpow(en,[4*n,4*n*(n-1)],j); r=elladd(en,p,q); print(ff(r)); if(r[2]!=0,print(ff(ellpow(en,r,-1)))) ) ) } torsion2(p)= { local(q,r,n); n=nn(p[1],p[2],1); print(ffso(p)); if(n!=nil, en=ec(n); q=[4*n,4*n*(n-1)]; r=elladd(en,p,q); print(ffso(r)) ) } \\ sqrt of rational number \\ return sqrt(x) if x >= 0 and issuqare(x) \\ nil if otherwise sqrtrational(x) = { local(n,d); n=numerator(x); d=denominator(x); if(issquare(n) && issquare(d), sqrtint(n)/sqrtint(d), nil ) } \\ height of [x:y:z]=log(max{|x|,|y|,|z|}) hh(p)= { local(x,y,z); x=abs(p[1]); y=abs(p[2]); z=abs(p[3]); log(max(max(x,y),z)) } \\ \\ analytic rank r of elliptic curve E \\ ar(e,s)= { local(L1,L2); L1=elllseries(e,s); L2=elllseries(e,s^2-s+1); (L2-L1)/((s-1)*L1) } \\ \\ L'(e,1) where s=1+\epsilon \\ lseries1(e2,s)= { local(L1,L2,e); e=ellchangecurve(e2,ellglobalred(e2)[2]); L1=elllseries(e,1); L2=elllseries(e,s); (L2-L1)/(s-1) }