[secp256k1.git] / sage / weierstrass_prover.sage

# Prover implementation for Weierstrass curves of the form
# y^2 = x^3 + A * x + B, specifically with a = 0 and b = 7, with group laws
# operating on affine and Jacobian coordinates, including the point at infinity
# represented by a 4th variable in coordinates.

load("group_prover.sage")


class affinepoint:
  def __init__(self, x, y, infinity=0):
    self.x = x
    self.y = y
    self.infinity = infinity
  def __str__(self):
    return "affinepoint(x=%s,y=%s,inf=%s)" % (self.x, self.y, self.infinity)


class jacobianpoint:
  def __init__(self, x, y, z, infinity=0):
    self.X = x
    self.Y = y
    self.Z = z
    self.Infinity = infinity
  def __str__(self):
    return "jacobianpoint(X=%s,Y=%s,Z=%s,inf=%s)" % (self.X, self.Y, self.Z, self.Infinity)


def point_at_infinity():
  return jacobianpoint(1, 1, 1, 1)


def negate(p):
  if p.__class__ == affinepoint:
    return affinepoint(p.x, -p.y)
  if p.__class__ == jacobianpoint:
    return jacobianpoint(p.X, -p.Y, p.Z)
  assert(False)


def on_weierstrass_curve(A, B, p):
  """Return a set of zero-expressions for an affine point to be on the curve"""
  return constraints(zero={p.x^3 + A*p.x + B - p.y^2: 'on_curve'})


def tangential_to_weierstrass_curve(A, B, p12, p3):
  """Return a set of zero-expressions for ((x12,y12),(x3,y3)) to be a line that is tangential to the curve at (x12,y12)"""
  return constraints(zero={
    (p12.y - p3.y) * (p12.y * 2) - (p12.x^2 * 3 + A) * (p12.x - p3.x): 'tangential_to_curve'
  })


def colinear(p1, p2, p3):
  """Return a set of zero-expressions for ((x1,y1),(x2,y2),(x3,y3)) to be collinear"""
  return constraints(zero={
    (p1.y - p2.y) * (p1.x - p3.x) - (p1.y - p3.y) * (p1.x - p2.x): 'colinear_1',
    (p2.y - p3.y) * (p2.x - p1.x) - (p2.y - p1.y) * (p2.x - p3.x): 'colinear_2',
    (p3.y - p1.y) * (p3.x - p2.x) - (p3.y - p2.y) * (p3.x - p1.x): 'colinear_3'
  })


def good_affine_point(p):
  return constraints(nonzero={p.x : 'nonzero_x', p.y : 'nonzero_y'})


def good_jacobian_point(p):
  return constraints(nonzero={p.X : 'nonzero_X', p.Y : 'nonzero_Y', p.Z^6 : 'nonzero_Z'})


def good_point(p):
  return constraints(nonzero={p.Z^6 : 'nonzero_X'})


def finite(p, *affine_fns):
  con = good_point(p) + constraints(zero={p.Infinity : 'finite_point'})
  if p.Z != 0:
    return con + reduce(lambda a, b: a + b, (f(affinepoint(p.X / p.Z^2, p.Y / p.Z^3)) for f in affine_fns), con)
  else:
    return con

def infinite(p):
  return constraints(nonzero={p.Infinity : 'infinite_point'})


def law_jacobian_weierstrass_add(A, B, pa, pb, pA, pB, pC):
  """Check whether the passed set of coordinates is a valid Jacobian add, given assumptions"""
  assumeLaw = (good_affine_point(pa) +
               good_affine_point(pb) +
               good_jacobian_point(pA) +
               good_jacobian_point(pB) +
               on_weierstrass_curve(A, B, pa) +
               on_weierstrass_curve(A, B, pb) +
               finite(pA) +
               finite(pB) +
               constraints(nonzero={pa.x - pb.x : 'different_x'}))
  require = (finite(pC, lambda pc: on_weierstrass_curve(A, B, pc) +
             colinear(pa, pb, negate(pc))))
  return (assumeLaw, require)


def law_jacobian_weierstrass_double(A, B, pa, pb, pA, pB, pC):
  """Check whether the passed set of coordinates is a valid Jacobian doubling, given assumptions"""
  assumeLaw = (good_affine_point(pa) +
               good_affine_point(pb) +
               good_jacobian_point(pA) +
               good_jacobian_point(pB) +
               on_weierstrass_curve(A, B, pa) +
               on_weierstrass_curve(A, B, pb) +
               finite(pA) +
               finite(pB) +
               constraints(zero={pa.x - pb.x : 'equal_x', pa.y - pb.y : 'equal_y'}))
  require = (finite(pC, lambda pc: on_weierstrass_curve(A, B, pc) +
             tangential_to_weierstrass_curve(A, B, pa, negate(pc))))
  return (assumeLaw, require)


def law_jacobian_weierstrass_add_opposites(A, B, pa, pb, pA, pB, pC):
  assumeLaw = (good_affine_point(pa) +
               good_affine_point(pb) +
               good_jacobian_point(pA) +
               good_jacobian_point(pB) +
               on_weierstrass_curve(A, B, pa) +
               on_weierstrass_curve(A, B, pb) +
               finite(pA) +
               finite(pB) +
               constraints(zero={pa.x - pb.x : 'equal_x', pa.y + pb.y : 'opposite_y'}))
  require = infinite(pC)
  return (assumeLaw, require)


def law_jacobian_weierstrass_add_infinite_a(A, B, pa, pb, pA, pB, pC):
  assumeLaw = (good_affine_point(pa) +
               good_affine_point(pb) +
               good_jacobian_point(pA) +
               good_jacobian_point(pB) +
               on_weierstrass_curve(A, B, pb) +
               infinite(pA) +
               finite(pB))
  require = finite(pC, lambda pc: constraints(zero={pc.x - pb.x : 'c.x=b.x', pc.y - pb.y : 'c.y=b.y'}))
  return (assumeLaw, require)


def law_jacobian_weierstrass_add_infinite_b(A, B, pa, pb, pA, pB, pC):
  assumeLaw = (good_affine_point(pa) +
               good_affine_point(pb) +
               good_jacobian_point(pA) +
               good_jacobian_point(pB) +
               on_weierstrass_curve(A, B, pa) +
               infinite(pB) +
               finite(pA))
  require = finite(pC, lambda pc: constraints(zero={pc.x - pa.x : 'c.x=a.x', pc.y - pa.y : 'c.y=a.y'}))
  return (assumeLaw, require)


def law_jacobian_weierstrass_add_infinite_ab(A, B, pa, pb, pA, pB, pC):
  assumeLaw = (good_affine_point(pa) +
               good_affine_point(pb) +
               good_jacobian_point(pA) +
               good_jacobian_point(pB) +
               infinite(pA) +
               infinite(pB))
  require = infinite(pC)
  return (assumeLaw, require)


laws_jacobian_weierstrass = {
  'add': law_jacobian_weierstrass_add,
  'double': law_jacobian_weierstrass_double,
  'add_opposite': law_jacobian_weierstrass_add_opposites,
  'add_infinite_a': law_jacobian_weierstrass_add_infinite_a,
  'add_infinite_b': law_jacobian_weierstrass_add_infinite_b,
  'add_infinite_ab': law_jacobian_weierstrass_add_infinite_ab
}


def check_exhaustive_jacobian_weierstrass(name, A, B, branches, formula, p):
  """Verify an implementation of addition of Jacobian points on a Weierstrass curve, by executing and validating the result for every possible addition in a prime field"""
  F = Integers(p)
  print("Formula %s on Z%i:" % (name, p))
  points = []
  for x in range(0, p):
    for y in range(0, p):
      point = affinepoint(F(x), F(y))
      r, e = concrete_verify(on_weierstrass_curve(A, B, point))
      if r:
        points.append(point)

  for za in range(1, p):
    for zb in range(1, p):
      for pa in points:
        for pb in points:
          for ia in range(2):
            for ib in range(2):
              pA = jacobianpoint(pa.x * F(za)^2, pa.y * F(za)^3, F(za), ia)
              pB = jacobianpoint(pb.x * F(zb)^2, pb.y * F(zb)^3, F(zb), ib)
              for branch in range(0, branches):
                assumeAssert, assumeBranch, pC = formula(branch, pA, pB)
                pC.X = F(pC.X)
                pC.Y = F(pC.Y)
                pC.Z = F(pC.Z)
                pC.Infinity = F(pC.Infinity)
                r, e = concrete_verify(assumeAssert + assumeBranch)
                if r:
                  match = False
                  for key in laws_jacobian_weierstrass:
                    assumeLaw, require = laws_jacobian_weierstrass[key](A, B, pa, pb, pA, pB, pC)
                    r, e = concrete_verify(assumeLaw)
                    if r:
                      if match:
                        print("  multiple branches for (%s,%s,%s,%s) + (%s,%s,%s,%s)" % (pA.X, pA.Y, pA.Z, pA.Infinity, pB.X, pB.Y, pB.Z, pB.Infinity))
                      else:
                        match = True
                      r, e = concrete_verify(require)
                      if not r:
                        print("  failure in branch %i for (%s,%s,%s,%s) + (%s,%s,%s,%s) = (%s,%s,%s,%s): %s" % (branch, pA.X, pA.Y, pA.Z, pA.Infinity, pB.X, pB.Y, pB.Z, pB.Infinity, pC.X, pC.Y, pC.Z, pC.Infinity, e))
  print()


def check_symbolic_function(R, assumeAssert, assumeBranch, f, A, B, pa, pb, pA, pB, pC):
  assumeLaw, require = f(A, B, pa, pb, pA, pB, pC)
  return check_symbolic(R, assumeLaw, assumeAssert, assumeBranch, require)

def check_symbolic_jacobian_weierstrass(name, A, B, branches, formula):
  """Verify an implementation of addition of Jacobian points on a Weierstrass curve symbolically"""
  R.<ax,bx,ay,by,Az,Bz,Ai,Bi> = PolynomialRing(QQ,8,order='invlex')
  lift = lambda x: fastfrac(R,x)
  ax = lift(ax)
  ay = lift(ay)
  Az = lift(Az)
  bx = lift(bx)
  by = lift(by)
  Bz = lift(Bz)
  Ai = lift(Ai)
  Bi = lift(Bi)

  pa = affinepoint(ax, ay, Ai)
  pb = affinepoint(bx, by, Bi)
  pA = jacobianpoint(ax * Az^2, ay * Az^3, Az, Ai)
  pB = jacobianpoint(bx * Bz^2, by * Bz^3, Bz, Bi)

  res = {}

  for key in laws_jacobian_weierstrass:
    res[key] = []

  print("Formula " + name + ":")
  count = 0
  for branch in range(branches):
    assumeFormula, assumeBranch, pC = formula(branch, pA, pB)
    pC.X = lift(pC.X)
    pC.Y = lift(pC.Y)
    pC.Z = lift(pC.Z)
    pC.Infinity = lift(pC.Infinity)

    for key in laws_jacobian_weierstrass:
      res[key].append((check_symbolic_function(R, assumeFormula, assumeBranch, laws_jacobian_weierstrass[key], A, B, pa, pb, pA, pB, pC), branch))

  for key in res:
    print("  %s:" % key)
    val = res[key]
    for x in val:
      if x[0] is not None:
        print("    branch %i: %s" % (x[1], x[0]))

  print()
Commit	Line	Data
03d4611c PW	1	# Prover implementation for Weierstrass curves of the form
	2	# y^2 = x^3 + A * x + B, specifically with a = 0 and b = 7, with group laws
	3	# operating on affine and Jacobian coordinates, including the point at infinity
	4	# represented by a 4th variable in coordinates.
	5
	6	load("group_prover.sage")
	7
	8
	9	class affinepoint:
	10	def __init__(self, x, y, infinity=0):
	11	self.x = x
	12	self.y = y
	13	self.infinity = infinity
	14	def __str__(self):
	15	return "affinepoint(x=%s,y=%s,inf=%s)" % (self.x, self.y, self.infinity)
	16
	17
	18	class jacobianpoint:
	19	def __init__(self, x, y, z, infinity=0):
	20	self.X = x
	21	self.Y = y
	22	self.Z = z
	23	self.Infinity = infinity
	24	def __str__(self):
	25	return "jacobianpoint(X=%s,Y=%s,Z=%s,inf=%s)" % (self.X, self.Y, self.Z, self.Infinity)
	26
	27
	28	def point_at_infinity():
	29	return jacobianpoint(1, 1, 1, 1)
	30
	31
	32	def negate(p):
	33	if p.__class__ == affinepoint:
	34	return affinepoint(p.x, -p.y)
	35	if p.__class__ == jacobianpoint:
	36	return jacobianpoint(p.X, -p.Y, p.Z)
	37	assert(False)
	38
	39
	40	def on_weierstrass_curve(A, B, p):
	41	"""Return a set of zero-expressions for an affine point to be on the curve"""
	42	return constraints(zero={p.x^3 + A*p.x + B - p.y^2: 'on_curve'})
	43
	44
	45	def tangential_to_weierstrass_curve(A, B, p12, p3):
	46	"""Return a set of zero-expressions for ((x12,y12),(x3,y3)) to be a line that is tangential to the curve at (x12,y12)"""
	47	return constraints(zero={
	48	(p12.y - p3.y) * (p12.y * 2) - (p12.x^2 * 3 + A) * (p12.x - p3.x): 'tangential_to_curve'
	49	})
	50
	51
	52	def colinear(p1, p2, p3):
	53	"""Return a set of zero-expressions for ((x1,y1),(x2,y2),(x3,y3)) to be collinear"""
	54	return constraints(zero={
	55	(p1.y - p2.y) * (p1.x - p3.x) - (p1.y - p3.y) * (p1.x - p2.x): 'colinear_1',
	56	(p2.y - p3.y) * (p2.x - p1.x) - (p2.y - p1.y) * (p2.x - p3.x): 'colinear_2',
	57	(p3.y - p1.y) * (p3.x - p2.x) - (p3.y - p2.y) * (p3.x - p1.x): 'colinear_3'
	58	})
	59
	60
	61	def good_affine_point(p):
	62	return constraints(nonzero={p.x : 'nonzero_x', p.y : 'nonzero_y'})
	63
	64
65	def good_jacobian_point(p):
66	return constraints(nonzero={p.X : 'nonzero_X', p.Y : 'nonzero_Y', p.Z^6 : 'nonzero_Z'})
67
68
69	def good_point(p):
70	return constraints(nonzero={p.Z^6 : 'nonzero_X'})
71
72
73	def finite(p, *affine_fns):
74	con = good_point(p) + constraints(zero={p.Infinity : 'finite_point'})
75	if p.Z != 0:
76	return con + reduce(lambda a, b: a + b, (f(affinepoint(p.X / p.Z^2, p.Y / p.Z^3)) for f in affine_fns), con)
77	else:
78	return con
79
80	def infinite(p):
81	return constraints(nonzero={p.Infinity : 'infinite_point'})
82
83
84	def law_jacobian_weierstrass_add(A, B, pa, pb, pA, pB, pC):
85	"""Check whether the passed set of coordinates is a valid Jacobian add, given assumptions"""
86	assumeLaw = (good_affine_point(pa) +
87	good_affine_point(pb) +
88	good_jacobian_point(pA) +
89	good_jacobian_point(pB) +
90	on_weierstrass_curve(A, B, pa) +
91	on_weierstrass_curve(A, B, pb) +
92	finite(pA) +
93	finite(pB) +
94	constraints(nonzero={pa.x - pb.x : 'different_x'}))
95	require = (finite(pC, lambda pc: on_weierstrass_curve(A, B, pc) +
96	colinear(pa, pb, negate(pc))))
97	return (assumeLaw, require)
98
99
100	def law_jacobian_weierstrass_double(A, B, pa, pb, pA, pB, pC):
101	"""Check whether the passed set of coordinates is a valid Jacobian doubling, given assumptions"""
102	assumeLaw = (good_affine_point(pa) +
103	good_affine_point(pb) +
104	good_jacobian_point(pA) +
105	good_jacobian_point(pB) +
106	on_weierstrass_curve(A, B, pa) +
107	on_weierstrass_curve(A, B, pb) +
108	finite(pA) +
109	finite(pB) +
110	constraints(zero={pa.x - pb.x : 'equal_x', pa.y - pb.y : 'equal_y'}))
111	require = (finite(pC, lambda pc: on_weierstrass_curve(A, B, pc) +
112	tangential_to_weierstrass_curve(A, B, pa, negate(pc))))
113	return (assumeLaw, require)
114
115
116	def law_jacobian_weierstrass_add_opposites(A, B, pa, pb, pA, pB, pC):
117	assumeLaw = (good_affine_point(pa) +
118	good_affine_point(pb) +
119	good_jacobian_point(pA) +
120	good_jacobian_point(pB) +
121	on_weierstrass_curve(A, B, pa) +
122	on_weierstrass_curve(A, B, pb) +
123	finite(pA) +
124	finite(pB) +
125	constraints(zero={pa.x - pb.x : 'equal_x', pa.y + pb.y : 'opposite_y'}))
126	require = infinite(pC)
127	return (assumeLaw, require)
128
129
130	def law_jacobian_weierstrass_add_infinite_a(A, B, pa, pb, pA, pB, pC):
131	assumeLaw = (good_affine_point(pa) +
132	good_affine_point(pb) +
133	good_jacobian_point(pA) +
134	good_jacobian_point(pB) +
135	on_weierstrass_curve(A, B, pb) +
136	infinite(pA) +
137	finite(pB))
138	require = finite(pC, lambda pc: constraints(zero={pc.x - pb.x : 'c.x=b.x', pc.y - pb.y : 'c.y=b.y'}))
139	return (assumeLaw, require)
140
141
142	def law_jacobian_weierstrass_add_infinite_b(A, B, pa, pb, pA, pB, pC):
143	assumeLaw = (good_affine_point(pa) +
144	good_affine_point(pb) +
145	good_jacobian_point(pA) +
146	good_jacobian_point(pB) +
147	on_weierstrass_curve(A, B, pa) +
148	infinite(pB) +
149	finite(pA))
150	require = finite(pC, lambda pc: constraints(zero={pc.x - pa.x : 'c.x=a.x', pc.y - pa.y : 'c.y=a.y'}))
151	return (assumeLaw, require)
152
153
154	def law_jacobian_weierstrass_add_infinite_ab(A, B, pa, pb, pA, pB, pC):
155	assumeLaw = (good_affine_point(pa) +
156	good_affine_point(pb) +
157	good_jacobian_point(pA) +
158	good_jacobian_point(pB) +
159	infinite(pA) +
160	infinite(pB))
161	require = infinite(pC)
162	return (assumeLaw, require)
163
164
165	laws_jacobian_weierstrass = {
166	'add': law_jacobian_weierstrass_add,
167	'double': law_jacobian_weierstrass_double,
168	'add_opposite': law_jacobian_weierstrass_add_opposites,
169	'add_infinite_a': law_jacobian_weierstrass_add_infinite_a,
170	'add_infinite_b': law_jacobian_weierstrass_add_infinite_b,
171	'add_infinite_ab': law_jacobian_weierstrass_add_infinite_ab
172	}
173
174
175	def check_exhaustive_jacobian_weierstrass(name, A, B, branches, formula, p):
176	"""Verify an implementation of addition of Jacobian points on a Weierstrass curve, by executing and validating the result for every possible addition in a prime field"""
177	F = Integers(p)
13c88efe	178	print("Formula %s on Z%i:" % (name, p))
03d4611c	179	points = []
13c88efe FC	180	for x in range(0, p):
13c88efe FC	181	for y in range(0, p):
03d4611c PW	182	point = affinepoint(F(x), F(y))
	183	r, e = concrete_verify(on_weierstrass_curve(A, B, point))
	184	if r:
	185	points.append(point)
	186
13c88efe FC	187	for za in range(1, p):
13c88efe FC	188	for zb in range(1, p):
03d4611c PW	189	for pa in points:
03d4611c PW	190	for pb in points:
13c88efe FC	191	for ia in range(2):
13c88efe FC	192	for ib in range(2):
03d4611c PW	193	pA = jacobianpoint(pa.x * F(za)^2, pa.y * F(za)^3, F(za), ia)
03d4611c PW	194	pB = jacobianpoint(pb.x * F(zb)^2, pb.y * F(zb)^3, F(zb), ib)
13c88efe	195	for branch in range(0, branches):
03d4611c PW	196	assumeAssert, assumeBranch, pC = formula(branch, pA, pB)
	197	pC.X = F(pC.X)
	198	pC.Y = F(pC.Y)
	199	pC.Z = F(pC.Z)
	200	pC.Infinity = F(pC.Infinity)
	201	r, e = concrete_verify(assumeAssert + assumeBranch)
	202	if r:
	203	match = False
	204	for key in laws_jacobian_weierstrass:
	205	assumeLaw, require = laws_jacobian_weierstrass[key](A, B, pa, pb, pA, pB, pC)
	206	r, e = concrete_verify(assumeLaw)
	207	if r:
	208	if match:
13c88efe	209	print(" multiple branches for (%s,%s,%s,%s) + (%s,%s,%s,%s)" % (pA.X, pA.Y, pA.Z, pA.Infinity, pB.X, pB.Y, pB.Z, pB.Infinity))
03d4611c PW	210	else:
	211	match = True
	212	r, e = concrete_verify(require)
	213	if not r:
13c88efe FC	214	print(" failure in branch %i for (%s,%s,%s,%s) + (%s,%s,%s,%s) = (%s,%s,%s,%s): %s" % (branch, pA.X, pA.Y, pA.Z, pA.Infinity, pB.X, pB.Y, pB.Z, pB.Infinity, pC.X, pC.Y, pC.Z, pC.Infinity, e))
13c88efe FC	215	print()
03d4611c PW	216
	217
	218	def check_symbolic_function(R, assumeAssert, assumeBranch, f, A, B, pa, pb, pA, pB, pC):
	219	assumeLaw, require = f(A, B, pa, pb, pA, pB, pC)
	220	return check_symbolic(R, assumeLaw, assumeAssert, assumeBranch, require)
	221
	222	def check_symbolic_jacobian_weierstrass(name, A, B, branches, formula):
	223	"""Verify an implementation of addition of Jacobian points on a Weierstrass curve symbolically"""
	224	R.<ax,bx,ay,by,Az,Bz,Ai,Bi> = PolynomialRing(QQ,8,order='invlex')
	225	lift = lambda x: fastfrac(R,x)
	226	ax = lift(ax)
	227	ay = lift(ay)
	228	Az = lift(Az)
	229	bx = lift(bx)
	230	by = lift(by)
	231	Bz = lift(Bz)
	232	Ai = lift(Ai)
	233	Bi = lift(Bi)
	234
	235	pa = affinepoint(ax, ay, Ai)
	236	pb = affinepoint(bx, by, Bi)
	237	pA = jacobianpoint(ax * Az^2, ay * Az^3, Az, Ai)
	238	pB = jacobianpoint(bx * Bz^2, by * Bz^3, Bz, Bi)
	239
	240	res = {}
	241
	242	for key in laws_jacobian_weierstrass:
	243	res[key] = []
	244
13c88efe	245	print("Formula " + name + ":")
03d4611c	246	count = 0
13c88efe	247	for branch in range(branches):
03d4611c PW	248	assumeFormula, assumeBranch, pC = formula(branch, pA, pB)
	249	pC.X = lift(pC.X)
	250	pC.Y = lift(pC.Y)
	251	pC.Z = lift(pC.Z)
	252	pC.Infinity = lift(pC.Infinity)
	253
	254	for key in laws_jacobian_weierstrass:
	255	res[key].append((check_symbolic_function(R, assumeFormula, assumeBranch, laws_jacobian_weierstrass[key], A, B, pa, pb, pA, pB, pC), branch))
	256
	257	for key in res:
13c88efe	258	print(" %s:" % key)
03d4611c PW	259	val = res[key]
	260	for x in val:
	261	if x[0] is not None:
13c88efe	262	print(" branch %i: %s" % (x[1], x[0]))
03d4611c	263
13c88efe	264	print()