Contact
CoCalc Logo Icon
StoreFeaturesDocsShareSupport News Sign UpSign In
| Download
Views: 138
For the sake of comparison we first give formulas for Weierstrass curves in projective coordinates; scroll down to see the formulas for Edwards curves
def add_projective(P,Q,A): """addition algorithm for projective points an elliptic curve in short Weierstrass form y^2=x^3+Ax+B -- handles all cases (including points at infinity)""" if P[2] == 0: return Q if Q[2] == 0: return P x1=P[0]; y1=P[1]; z1=P[2]; x2=Q[0]; y2=Q[1]; z2=Q[2] # We use formulas from http://hyperelliptic.org/EFD/g1p/auto-shortw-projective.html#addition-add-1998-cmo-2 x1z2 = x1*z2; x2z1 = x2*z1; y1z2 = y1*z2; y2z1 = y2*z1 if x1z2 != x2z1: # usual case: P and Q are distinct and not opposites z1z2 = z1*z2 u = y2z1-y1z2; uu = u^2 v = x2z1-x1z2; vv = v^2; vvv = v^3 r = vv*x1z2 s = uu*z1z2 - vvv - 2*r x3 = v*s y3 = u*(r-s)-vvv*y1z2 z3 = vvv*z1z2 else: if y1z2 == -y2z1: return (0,1,0) # P = -Q (includes case where P=Q is 2-torsion) # P = Q so we are doubling xx = x1^2; yy = y1^2; zz = z1^2 w = A*zz + 3*xx s = 2*y1*z1; ss = s^2; sss = s*ss r = y1*s; rr = r^2 b = (x1+r)^2-xx-rr h = w^2-2*b x3 = h*s y3 = w*(b-h)-2*rr z3 = sss return (x3,y3,z3) def negate(P): """negates the point P on an elliptic curve in short Weierstrass form""" if P[2] == 0: return (0,1,0) return (P[0],-P[1],P[2]) def dbl_projective(P,A): """doubles a projective non-2-torsion point on an elliptic curve in short Weierstrass form""" x1 = P[0]; y1 = P[1]; z1 = P[2] xx = x1^2; yy = y1^2; zz = z1^2 w = A*zz + 3*xx s = 2*y1*z1; ss = s^2; sss = s*ss r = y1*s; rr = r^2 b = (x1+r)^2-xx-rr h = w^2-2*b x3 = h*s y3 = w*(b-h)-2*rr z3 = sss return (x3,y3,z3) def equal_projective(P,Q): """Test whether two projective points are equal or not""" return P[0]*Q[2] == Q[0]*P[2] and P[1]*Q[2] == Q[1]*P[2] def random_point(A,B): """generates a random affine point on the elliptic curve y^2 = x^2 + Ax + B""" F=A.parent() x0=F.random_element(); t=(x0^3+A*x0+B) while not t.is_square(): x0=F.random_element(); t=(x0^3+A*x0+B) return (x0,sqrt(t),1)
p=random_prime(2^256,lbound=2^255) F=GF(p) A=F(3); B=F.random_element() P0=random_point(A,B) Q0=random_point(A,B) print "starting tests..." timeit('add_projective(P0,Q0,A)',number=50000) timeit('add_projective(P0,P0,A)',number=50000) timeit('dbl_projective(P0,A)',number=50000)
starting tests... 50000 loops, best of 3: 13.6 µs per loop 50000 loops, best of 3: 22.7 µs per loop 50000 loops, best of 3: 20.3 µs per loop
Note that distinct triples of proejctive coordinates may define the same projective point. To compare points for equality we need to use the equal_projective function above (which also works for Edwards curves)
R1=add_projective(P0,Q0,A) R2=add_projective(Q0,P0,A) print "R1 =", R1 print "R2 =", R2 print "R1 == R2:",R1 == R2 print "equal_projective(R1,R2):",equal_projective(R1,R2)
R1 = (88417280862405198524498904230150093073040285267296490962688023695194294232737, 81803801345615446712151118453505899941011682403746951648224039427654521101855, 46633157038940499649062681694319626916323919252431273559123417012920789472386) R2 = (18524196991610979948129892573767159495541479256700257005341500524559971803562, 25137676508400731760477678350411352627570082120249796319805484792099744934444, 60308320815075678823566115109597625652257845271565474408906107206833476563913) R1 == R2: False equal_projective(R1,R2): True
def add_edwards(P,Q,d): """adds two projective points on an elliptic curve in Edwards form (assumes d non-square, so operation is complete)""" # uses unoptimized formula derived in class (could save 1M with optimization) x1=P[0]; y1=P[1]; z1=P[2]; x2=Q[0]; y2=Q[1]; z2=Q[2] x1y2 = x1*y2; x2y1 = x2*y1 r = z1*z2; rr = r^2; s = x1y2+x2y1; t = d*x1y2*x2y1; u = y1*y2 - x1*x2 v = rr+t; w = rr-t x3 = r*s*w y3 = r*u*v z3 = v*w return (x3,y3,z3) def dbl_edwards(P): """doubles a projective point on an elliptic curve in Edwards form""" # uses optimized formula from Bernstein-Lange x1 = P[0]; y1 = P[1]; z1 = P[2] B=(x1+y1)^2; C=x1^2; D=y1^2; E = C+D; J = E-2*z1^2; x3 = (B-E)*J; y3 = E*(C-D); z3 = E*J return (x3,y3,z3) def neg_edwards(P): return (-P[0],P[1],P[2]) def random_edwards_point(d): """generates a random projective point on the Edwards curve x^2 + y^2 = 1+d*x^2*y^2 (assumes d non-square)""" F=d.parent() x0=F.random_element(); t=(1-x0^2)/(1-d*x0^2) while not t.is_square(): x0=F.random_element(); t=(1-x0^2)/(1-d*x0^2) return (x0,sqrt(t),1) def on_edwards_curve (P,d): x0=P[0]; y0=P[1]; z0=P[2] return z0^2*(x0^2+y0^2) == z0^4 + d*x0^2*y0^2 def is_identity_edwards (P): if P[0] == 0 and P[1] == P[2]: return true else: return false
p = 17 F=GF(p) d = F.multiplicative_generator() print "Points on Edwards curve x^2 + y^2 = 1 + %dx^2y^2 over F_%d"%(d,p) for a in F: for b in F: if on_edwards_curve([a,b,1],d): print [a,b]
Points on Edwards curve x^2 + y^2 = 1 + 3x^2y^2 over F_17 [0, 1] [0, 16] [1, 0] [2, 5] [2, 12] [3, 4] [3, 13] [4, 3] [4, 14] [5, 2] [5, 15] [7, 7] [7, 10] [10, 7] [10, 10] [12, 2] [12, 15] [13, 3] [13, 14] [14, 4] [14, 13] [15, 5] [15, 12] [16, 0]
p=random_prime(2^256,lbound=2^255) F=GF(p) d=F(2); while d.is_square(): d+=1; P0=random_edwards_point(d) Q0=random_edwards_point(d) print "starting tests..." timeit('add_edwards(P0,Q0,d)',number=100000) timeit('add_edwards(P0,P0,d)',number=100000) timeit('dbl_edwards(P0)',number=100000)
starting tests... 100000 loops, best of 3: 9.46 µs per loop 100000 loops, best of 3: 9.17 µs per loop 100000 loops, best of 3: 6.4 µs per loop