From: Niels Möller <nisse@lysator.liu.se>
Date: Mon, 8 Nov 2021 16:41:54 +0000 (+0100)
Subject: Implement secp256r1 square root, based on patch by Wim Lewis.
X-Git-Tag: nettle_3.8_release_20220602~80
X-Git-Url: http://git.ipfire.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=bc07754f0eac9854d39a4524444fde887f55a82d;p=thirdparty%2Fnettle.git

Implement secp256r1 square root, based on patch by Wim Lewis.
---

diff --git a/ChangeLog b/ChangeLog
index 30aa255f..57b10eaa 100644
--- a/ChangeLog
+++ b/ChangeLog
@@ -6,6 +6,7 @@
 	Update all curve definitions.
 	* ecc-secp192r1.c (ECC_SECP192R1_SQRT_ITCH): New constant.
 	(ecc_secp192r1_sqrt): New function.
+	* ecc-secp256r1.c (ecc_secp256r1_sqrt): New function.
 
 	* testsuite/ecc-sqrt-test.c (test_sqrt): New function.
 	(test_sqrt_ratio): Renamed function (was test_modulo).
diff --git a/ecc-secp256r1.c b/ecc-secp256r1.c
index 3bdbdc21..e1a14b90 100644
--- a/ecc-secp256r1.c
+++ b/ecc-secp256r1.c
@@ -263,6 +263,70 @@ ecc_secp256r1_inv (const struct ecc_modulo *p,
   ecc_mod_pow_2k_mul (p, rp, rp, 15, a15m1, tp);/* a^{2^{239} - 2^{207} + 2^{175} + 2^{79} - 1} */
   ecc_mod_pow_2k_mul (p, rp, rp, 15, a15m1, tp);/* a^{2^{254} - 2^{222} + 2^{190} + 2^{94} - 1} */
   ecc_mod_pow_2k_mul (p, rp, rp, 2, ap, tp); 	/* a^{2^{256} - 2^{224} + 2^{192} + 2^{96} - 3} */
+
+#undef a5m1
+#undef t0
+#undef a15m1
+#undef a32m1
+#undef tp
+}
+
+/* To guarantee that inputs to ecc_mod_zero_p are in the required range. */
+#if ECC_LIMB_SIZE * GMP_NUMB_BITS != 256
+#error Unsupported limb size
+#endif
+
+#define ECC_SECP256R1_SQRT_ITCH (3*ECC_LIMB_SIZE)
+
+static int
+ecc_secp256r1_sqrt (const struct ecc_modulo *m,
+		    mp_limb_t *rp,
+		    const mp_limb_t *cp,
+		    mp_limb_t *scratch)
+{
+  /* This computes the square root modulo p256 using the identity:
+
+     sqrt(c) = c^(2^254 − 2^222 + 2^190 + 2^94)  (mod P-256)
+
+     which can be seen as a special case of Tonelli-Shanks with e=1.
+
+     It would be nice to share part of the addition chain between inverse and sqrt.
+
+     We need
+
+       p-2 = 2^{256} - 2^{224} + 2^{192} + 2^{96} - 3 (inverse)
+
+     and
+
+       (p+1)/4 = 2^{254} − 2^{222} + 2^{190} + 2^{94} (sqrt)
+
+     which we can both get conveniently from
+
+       (p-3)/4 = 2^{254} − 2^{222} + 2^{190} + 2^{94} - 1
+
+     But addition chain for 2^{94} - 1 appears to cost a few more mul
+     operations than the current, separate, chains. */
+
+#define t0 scratch
+#define tp (scratch + ECC_LIMB_SIZE)
+
+  ecc_mod_sqr        (m, rp, cp, tp);		/* c^2 */
+  ecc_mod_mul        (m, t0, rp, cp, tp);	/* c^3 */
+  ecc_mod_pow_2kp1   (m, rp, t0, 2, tp);	/* c^(2^4 - 1) */
+  ecc_mod_pow_2kp1   (m, t0, rp, 4, tp);	/* c^(2^8 - 1) */
+  ecc_mod_pow_2kp1   (m, rp, t0, 8, tp);	/* c^(2^16 - 1) */
+  ecc_mod_pow_2kp1   (m, t0, rp, 16, tp);	/* c^(2^32 - 1) */
+  ecc_mod_pow_2k_mul (m, rp, t0, 32, cp, tp);	/* c^(2^64 - 2^32 + 1) */
+  ecc_mod_pow_2k_mul (m, t0, rp, 96, cp, tp);	/* c^(2^160 - 2^128 + 2^96 + 1) */
+  ecc_mod_pow_2k     (m, rp, t0, 94,     tp);	/* c^(2^254 - 2^222 + 2^190 + 2^94) */
+
+  ecc_mod_sqr (m, t0, rp, tp);
+  ecc_mod_sub (m, t0, t0, cp);
+
+  return ecc_mod_zero_p (m, t0);
+#undef t0
+#undef tp
+
 }
 
 const struct ecc_curve _nettle_secp_256r1 =
@@ -273,7 +337,7 @@ const struct ecc_curve _nettle_secp_256r1 =
     ECC_BMODP_SIZE,
     ECC_REDC_SIZE,
     ECC_SECP256R1_INV_ITCH,
-    0, 
+    ECC_SECP256R1_SQRT_ITCH,
     0,
 
     ecc_p,
@@ -285,7 +349,7 @@ const struct ecc_curve _nettle_secp_256r1 =
     ecc_secp256r1_modp,
     USE_REDC ? ecc_secp256r1_redc : ecc_secp256r1_modp,
     ecc_secp256r1_inv,
-    NULL,
+    ecc_secp256r1_sqrt,
     NULL,
   },
   {