sage/gen_split_lambda_constants.sage

   1 """ Generates the constants used in secp256k1_scalar_split_lambda.
   2
   3 See the comments for secp256k1_scalar_split_lambda in src/scalar_impl.h for detailed explanations.
   4 """
   5
   6 load("secp256k1_params.sage")
   7
   8 def inf_norm(v):
   9     """Returns the infinity norm of a vector."""
  10     return max(map(abs, v))
  11
  12 def gauss_reduction(i1, i2):
  13     v1, v2 = i1.copy(), i2.copy()
  14     while True:
  15         if inf_norm(v2) < inf_norm(v1):
  16             v1, v2 = v2, v1
  17         # This is essentially
  18         #    m = round((v1[0]*v2[0] + v1[1]*v2[1]) / (inf_norm(v1)**2))
  19         # (rounding to the nearest integer) without relying on floating point arithmetic.
  20         m = ((v1[0]*v2[0] + v1[1]*v2[1]) + (inf_norm(v1)**2) // 2) // (inf_norm(v1)**2)
  21         if m == 0:
  22             return v1, v2
  23         v2[0] -= m*v1[0]
  24         v2[1] -= m*v1[1]
  25
  26 def find_split_constants_gauss():
  27     """Find constants for secp256k1_scalar_split_lamdba using gauss reduction."""
  28     (v11, v12), (v21, v22) = gauss_reduction([0, N], [1, int(LAMBDA)])
  29
  30     # We use related vectors in secp256k1_scalar_split_lambda.
  31     A1, B1 = -v21, -v11
  32     A2, B2 = v22, -v21
  33
  34     return A1, B1, A2, B2
  35
  36 def find_split_constants_explicit_tof():
  37     """Find constants for secp256k1_scalar_split_lamdba using the trace of Frobenius.
  38
  39     See Benjamin Smith: "Easy scalar decompositions for efficient scalar multiplication on
  40     elliptic curves and genus 2 Jacobians" (https://eprint.iacr.org/2013/672), Example 2
  41     """
  42     assert P % 3 == 1 # The paper says P % 3 == 2 but that appears to be a mistake, see [10].
  43     assert C.j_invariant() == 0
  44
  45     t = C.trace_of_frobenius()
  46
  47     c = Integer(sqrt((4*P - t**2)/3))
  48     A1 = Integer((t - c)/2 - 1)
  49     B1 = c
  50
  51     A2 = Integer((t + c)/2 - 1)
  52     B2 = Integer(1 - (t - c)/2)
  53
  54     # We use a negated b values in secp256k1_scalar_split_lambda.
  55     B1, B2 = -B1, -B2
  56
  57     return A1, B1, A2, B2
  58
  59 A1, B1, A2, B2 = find_split_constants_explicit_tof()
  60
  61 # For extra fun, use an independent method to recompute the constants.
  62 assert (A1, B1, A2, B2) == find_split_constants_gauss()
  63
  64 # PHI : Z[l] -> Z_n where phi(a + b*l) == a + b*lambda mod n.
  65 def PHI(a,b):
  66     return Z(a + LAMBDA*b)
  67
  68 # Check that (A1, B1) and (A2, B2) are in the kernel of PHI.
  69 assert PHI(A1, B1) == Z(0)
  70 assert PHI(A2, B2) == Z(0)
  71
  72 # Check that the parallelogram generated by (A1, A2) and (B1, B2)
  73 # is a fundamental domain by containing exactly N points.
  74 # Since the LHS is the determinant and N != 0, this also checks that
  75 # (A1, A2) and (B1, B2) are linearly independent. By the previous
  76 # assertions, (A1, A2) and (B1, B2) are a basis of the kernel.
  77 assert A1*B2 - B1*A2 == N
  78
  79 # Check that their components are short enough.
  80 assert (A1 + A2)/2 < sqrt(N)
  81 assert B1 < sqrt(N)
  82 assert B2 < sqrt(N)
  83
  84 G1 = round((2**384)*B2/N)
  85 G2 = round((2**384)*(-B1)/N)
  86
  87 def rnddiv2(v):
  88     if v & 1:
  89         v += 1
  90     return v >> 1
  91
  92 def scalar_lambda_split(k):
  93     """Equivalent to secp256k1_scalar_lambda_split()."""
  94     c1 = rnddiv2((k * G1) >> 383)
  95     c2 = rnddiv2((k * G2) >> 383)
  96     c1 = (c1 * -B1) % N
  97     c2 = (c2 * -B2) % N
  98     r2 = (c1 + c2) % N
  99     r1 = (k + r2 * -LAMBDA) % N
 100     return (r1, r2)
 101
 102 # The result of scalar_lambda_split can depend on the representation of k (mod n).
 103 SPECIAL = (2**383) // G2 + 1
 104 assert scalar_lambda_split(SPECIAL) != scalar_lambda_split(SPECIAL + N)
 105
 106 print('  A1     =', hex(A1))
 107 print(' -B1     =', hex(-B1))
 108 print('  A2     =', hex(A2))
 109 print(' -B2     =', hex(-B2))
 110 print('         =', hex(Z(-B2)))
 111 print(' -LAMBDA =', hex(-LAMBDA))
 112
 113 print('  G1     =', hex(G1))
 114 print('  G2     =', hex(G2))