jacobian.py 2.4 KB

12345678910111213141516171819202122232425262728293031323334353637383940414243444546474849505152535455565758596061626364656667686970717273747576777879808182838485868788899091929394959697
  1. from typing import (
  2. Tuple,
  3. )
  4. from eth_keys.constants import (
  5. IDENTITY_POINTS,
  6. SECPK1_A as A,
  7. SECPK1_N as N,
  8. SECPK1_P as P,
  9. )
  10. def inv(a: int, n: int) -> int:
  11. if a == 0:
  12. return 0
  13. lm, hm = 1, 0
  14. low, high = a % n, n
  15. while low > 1:
  16. r = high // low
  17. nm, new = hm - lm * r, high - low * r
  18. lm, low, hm, high = nm, new, lm, low
  19. return lm % n
  20. def to_jacobian(p: Tuple[int, int]) -> Tuple[int, int, int]:
  21. o = (p[0], p[1], 1)
  22. return o
  23. def jacobian_double(p: Tuple[int, int, int]) -> Tuple[int, int, int]:
  24. if not p[1]:
  25. return (0, 0, 0)
  26. ysq = (p[1] ** 2) % P
  27. S = (4 * p[0] * ysq) % P
  28. M = (3 * p[0] ** 2 + A * p[2] ** 4) % P
  29. nx = (M**2 - 2 * S) % P
  30. ny = (M * (S - nx) - 8 * ysq**2) % P
  31. nz = (2 * p[1] * p[2]) % P
  32. return (nx, ny, nz)
  33. def jacobian_add(
  34. p: Tuple[int, int, int], q: Tuple[int, int, int]
  35. ) -> Tuple[int, int, int]:
  36. if not p[1]:
  37. return q
  38. if not q[1]:
  39. return p
  40. U1 = (p[0] * q[2] ** 2) % P
  41. U2 = (q[0] * p[2] ** 2) % P
  42. S1 = (p[1] * q[2] ** 3) % P
  43. S2 = (q[1] * p[2] ** 3) % P
  44. if U1 == U2:
  45. if S1 != S2:
  46. return (0, 0, 1)
  47. return jacobian_double(p)
  48. H = U2 - U1
  49. R = S2 - S1
  50. H2 = (H * H) % P
  51. H3 = (H * H2) % P
  52. U1H2 = (U1 * H2) % P
  53. nx = (R**2 - H3 - 2 * U1H2) % P
  54. ny = (R * (U1H2 - nx) - S1 * H3) % P
  55. nz = (H * p[2] * q[2]) % P
  56. return (nx, ny, nz)
  57. def from_jacobian(p: Tuple[int, int, int]) -> Tuple[int, int]:
  58. z = inv(p[2], P)
  59. return ((p[0] * z**2) % P, (p[1] * z**3) % P)
  60. def jacobian_multiply(a: Tuple[int, int, int], n: int) -> Tuple[int, int, int]:
  61. if a[1] == 0 or n == 0:
  62. return (0, 0, 1)
  63. if n == 1:
  64. return a
  65. if n < 0 or n >= N:
  66. return jacobian_multiply(a, n % N)
  67. if (n % 2) == 0:
  68. return jacobian_double(jacobian_multiply(a, n // 2))
  69. elif (n % 2) == 1:
  70. return jacobian_add(jacobian_double(jacobian_multiply(a, n // 2)), a)
  71. else:
  72. raise Exception("Invariant: Unreachable code path")
  73. def fast_multiply(a: Tuple[int, int], n: int) -> Tuple[int, int]:
  74. return from_jacobian(jacobian_multiply(to_jacobian(a), n))
  75. def fast_add(a: Tuple[int, int], b: Tuple[int, int]) -> Tuple[int, int]:
  76. return from_jacobian(jacobian_add(to_jacobian(a), to_jacobian(b)))
  77. def is_identity(p: Tuple[int, int, int]) -> bool:
  78. return p in IDENTITY_POINTS