Implement secp192r1 square root, based on patch by Wim Lewis.

12 files changed, 201 insertions, 6 deletions

diff --git a/ChangeLog b/ChangeLog
index 3b411504..30aa255f 100644
--- a/ChangeLog
+++ b/ChangeLog
@@ -1,3 +1,16 @@
+2021-11-08  Niels Möller  <nisse@lysator.liu.se>
+
+	Square root functions, based on patch by Wim Lewis.
+	* ecc-internal.h (ecc_mod_sqrt_func): New typedef.
+	(struct ecc_modulo): Add sqrt function pointer and sqrt_itch.
+	Update all curve definitions.
+	* ecc-secp192r1.c (ECC_SECP192R1_SQRT_ITCH): New constant.
+	(ecc_secp192r1_sqrt): New function.
+
+	* testsuite/ecc-sqrt-test.c (test_sqrt): New function.
+	(test_sqrt_ratio): Renamed function (was test_modulo).
+	(test_main): Test sqrt function, for curves that define it.
+
 2021-11-07  Niels Möller  <nisse@lysator.liu.se>
 
 	* ecc-internal.h (struct ecc_modulo): Renamed sqrt_itch to
diff --git a/ecc-curve25519.c b/ecc-curve25519.c
index e461c197..56abcf23 100644
--- a/ecc-curve25519.c
+++ b/ecc-curve25519.c
@@ -260,6 +260,7 @@ const struct ecc_curve _nettle_curve25519 =
     ECC_BMODP_SIZE,
     0,
     ECC_25519_INV_ITCH,
+    0,
     ECC_25519_SQRT_RATIO_ITCH,
 
     ecc_p,
@@ -271,6 +272,7 @@ const struct ecc_curve _nettle_curve25519 =
     ecc_curve25519_modp,
     ecc_curve25519_modp,
     ecc_curve25519_inv,
+    NULL,
     ecc_curve25519_sqrt_ratio,
   },
   {
@@ -280,6 +282,7 @@ const struct ecc_curve _nettle_curve25519 =
     0,
     ECC_MOD_INV_ITCH (ECC_LIMB_SIZE),
     0,
+    0,
 
     ecc_q,
     ecc_Bmodq,  
@@ -291,6 +294,7 @@ const struct ecc_curve _nettle_curve25519 =
     ecc_curve25519_modq,
     ecc_mod_inv,
     NULL,
+    NULL,
   },
 
   0, /* No redc */
diff --git a/ecc-curve448.c b/ecc-curve448.c
index 67d197eb..1bd4e11f 100644
--- a/ecc-curve448.c
+++ b/ecc-curve448.c
@@ -214,6 +214,7 @@ const struct ecc_curve _nettle_curve448 =
     ECC_BMODP_SIZE,
     0,
     ECC_CURVE448_INV_ITCH,
+    0,
     ECC_CURVE448_SQRT_RATIO_ITCH,
 
     ecc_p,
@@ -225,6 +226,7 @@ const struct ecc_curve _nettle_curve448 =
     ecc_curve448_modp,
     ecc_curve448_modp,
     ecc_curve448_inv,
+    NULL,
     ecc_curve448_sqrt_ratio,
   },
   {
@@ -234,6 +236,7 @@ const struct ecc_curve _nettle_curve448 =
     0,
     ECC_MOD_INV_ITCH (ECC_LIMB_SIZE),
     0,
+    0,
 
     ecc_q,
     ecc_Bmodq,
@@ -245,6 +248,7 @@ const struct ecc_curve _nettle_curve448 =
     ecc_mod,
     ecc_mod_inv,
     NULL,
+    NULL,
   },
 
   0, /* No redc */
diff --git a/ecc-gost-gc256b.c b/ecc-gost-gc256b.c
index 988368e9..0cf753e4 100644
--- a/ecc-gost-gc256b.c
+++ b/ecc-gost-gc256b.c
@@ -66,6 +66,7 @@ const struct ecc_curve _nettle_gost_gc256b =
     ECC_REDC_SIZE,
     ECC_MOD_INV_ITCH (ECC_LIMB_SIZE),
     0,
+    0,
 
     ecc_p,
     ecc_Bmodp,
@@ -77,6 +78,7 @@ const struct ecc_curve _nettle_gost_gc256b =
     ecc_gost_gc256b_modp,
     ecc_mod_inv,
     NULL,
+    NULL,
   },
   {
     256,
@@ -85,6 +87,7 @@ const struct ecc_curve _nettle_gost_gc256b =
     0,
     ECC_MOD_INV_ITCH (ECC_LIMB_SIZE),
     0,
+    0,
 
     ecc_q,
     ecc_Bmodq,
@@ -96,6 +99,7 @@ const struct ecc_curve _nettle_gost_gc256b =
     ecc_gost_gc256b_modq,
     ecc_mod_inv,
     NULL,
+    NULL,
   },
 
   USE_REDC,
diff --git a/ecc-gost-gc512a.c b/ecc-gost-gc512a.c
index 0b9864ef..338ed001 100644
--- a/ecc-gost-gc512a.c
+++ b/ecc-gost-gc512a.c
@@ -66,6 +66,7 @@ const struct ecc_curve _nettle_gost_gc512a =
     ECC_REDC_SIZE,
     ECC_MOD_INV_ITCH (ECC_LIMB_SIZE),
     0,
+    0,
 
     ecc_p,
     ecc_Bmodp,
@@ -77,6 +78,7 @@ const struct ecc_curve _nettle_gost_gc512a =
     ecc_gost_gc512a_modp,
     ecc_mod_inv,
     NULL,
+    NULL,
   },
   {
     512,
@@ -85,6 +87,7 @@ const struct ecc_curve _nettle_gost_gc512a =
     0,
     ECC_MOD_INV_ITCH (ECC_LIMB_SIZE),
     0,
+    0,
 
     ecc_q,
     ecc_Bmodq,
@@ -96,6 +99,7 @@ const struct ecc_curve _nettle_gost_gc512a =
     ecc_gost_gc512a_modq,
     ecc_mod_inv,
     NULL,
+    NULL,
   },
 
   USE_REDC,
diff --git a/ecc-internal.h b/ecc-internal.h
index a3e19331..277f5a02 100644
--- a/ecc-internal.h
+++ b/ecc-internal.h
@@ -125,6 +125,11 @@ typedef void ecc_mod_inv_func (const struct ecc_modulo *m,
 			       mp_limb_t *vp, const mp_limb_t *ap,
 			       mp_limb_t *scratch);
 
+/* Computes the square root of ap mod p. No overlap between input and output. */
+typedef int ecc_mod_sqrt_func (const struct ecc_modulo *m,
+			       mp_limb_t *vp, const mp_limb_t *ap,
+			       mp_limb_t *scratch);
+
 /* Computes the square root of (u/v) (mod p). */
 typedef int ecc_mod_sqrt_ratio_func (const struct ecc_modulo *m,
 				     mp_limb_t *rp,
@@ -161,6 +166,7 @@ struct ecc_modulo
   unsigned short B_size;
   unsigned short redc_size;
   unsigned short invert_itch;
+  unsigned short sqrt_itch;
   unsigned short sqrt_ratio_itch;
 
   const mp_limb_t *m;
@@ -179,6 +185,7 @@ struct ecc_modulo
   /* For moduli where we use redc, the invert and sqrt functions work
      with inputs and outputs in redc form. */
   ecc_mod_inv_func *invert;
+  ecc_mod_sqrt_func *sqrt;
   ecc_mod_sqrt_ratio_func *sqrt_ratio;
 };
 
diff --git a/ecc-secp192r1.c b/ecc-secp192r1.c
index abdad0c9..391ba528 100644
--- a/ecc-secp192r1.c
+++ b/ecc-secp192r1.c
@@ -2,7 +2,8 @@
 
    Compile time constant (but machine dependent) tables.
 
-   Copyright (C) 2013, 2014 Niels Möller
+   Copyright (C) 2013, 2014, 2019, 2021 Niels Möller
+   Copyright (C) 2019 Wim Lewis
 
    This file is part of GNU Nettle.
 
@@ -180,6 +181,57 @@ ecc_secp192r1_inv (const struct ecc_modulo *p,
   ecc_mod_sqr (p, rp, rp, tp);		/* a^{2^191 - 2^63 - 2} */
   ecc_mod_sqr (p, rp, rp, tp);		/* a^{2^192 - 2^64 - 4} */
   ecc_mod_mul (p, rp, rp, ap, tp);
+
+#undef a62m1
+#undef t0
+#undef tp
+}
+
+/* To guarantee that inputs to ecc_mod_zero_p are in the required range. */
+#if ECC_LIMB_SIZE * GMP_NUMB_BITS != 192
+#error Unsupported limb size
+#endif
+
+#define ECC_SECP192R1_SQRT_ITCH (3*ECC_LIMB_SIZE)
+
+static int
+ecc_secp192r1_sqrt (const struct ecc_modulo *p,
+		    mp_limb_t *rp,
+		    const mp_limb_t *cp,
+		    mp_limb_t *scratch)
+{
+  /* This computes the square root modulo p192 using the identity:
+
+     sqrt(c) = c^(2^190 - 2^62)  (mod P-192)
+
+     which can be seen as a special case of Tonelli-Shanks with e=1.
+  */
+
+  /* We need one t0 (and use clobbering rp) and scratch space for mul and sqr. */
+
+#define t0 scratch
+#define tp (scratch + ECC_LIMB_SIZE)
+
+  ecc_mod_sqr(p, rp, cp, tp);		/* c^2 */
+  ecc_mod_mul(p, rp, rp, cp, tp);	/* c^3 */
+  ecc_mod_pow_2kp1(p, t0, rp, 2, tp);	/* c^(2^4 - 1) */
+  ecc_mod_pow_2kp1(p, rp, t0, 4, tp);	/* c^(2^8 - 1) */
+  ecc_mod_pow_2kp1(p, t0, rp, 8,  tp);	/* c^(2^16 - 1) */
+  ecc_mod_pow_2kp1(p, rp, t0, 16, tp);	/* c^(2^32 - 1) */
+  ecc_mod_pow_2kp1(p, t0, rp, 32, tp);	/* c^(2^64 - 1) */
+  ecc_mod_pow_2kp1(p, rp, t0, 64, tp);	/* c^(2^128 - 1) */
+
+  ecc_mod_pow_2k    (p, rp, rp,     62, tp);   /* c^(2^190 - 2^62) */
+
+  /* Check that input was a square, R^2 = C, for non-squares we'd get
+     R^2 = -C. */
+  ecc_mod_sqr(p, t0, rp, tp);
+  ecc_mod_sub(p, t0, t0, cp);
+
+  return ecc_mod_zero_p (p, t0);
+
+#undef t0
+#undef tp
 }
 
 const struct ecc_curve _nettle_secp_192r1 =
@@ -190,6 +242,7 @@ const struct ecc_curve _nettle_secp_192r1 =
     ECC_BMODP_SIZE,
     ECC_REDC_SIZE,
     ECC_SECP192R1_INV_ITCH,
+    ECC_SECP192R1_SQRT_ITCH,
     0,
 
     ecc_p,
@@ -201,6 +254,7 @@ const struct ecc_curve _nettle_secp_192r1 =
     ecc_secp192r1_modp,
     ecc_secp192r1_modp,
     ecc_secp192r1_inv,
+    ecc_secp192r1_sqrt,
     NULL,
   },
   {
@@ -210,6 +264,7 @@ const struct ecc_curve _nettle_secp_192r1 =
     0,
     ECC_MOD_INV_ITCH (ECC_LIMB_SIZE),
     0,
+    0,
 
     ecc_q,
     ecc_Bmodq,
@@ -221,6 +276,7 @@ const struct ecc_curve _nettle_secp_192r1 =
     ecc_mod,
     ecc_mod_inv,
     NULL,
+    NULL,
   },
   
   USE_REDC,
diff --git a/ecc-secp224r1.c b/ecc-secp224r1.c
index 91c337d4..767eccf7 100644
--- a/ecc-secp224r1.c
+++ b/ecc-secp224r1.c
@@ -121,6 +121,7 @@ const struct ecc_curve _nettle_secp_224r1 =
     -ECC_REDC_SIZE,
     ECC_SECP224R1_INV_ITCH,
     0,
+    0,
 
     ecc_p,
     ecc_Bmodp,
@@ -132,6 +133,7 @@ const struct ecc_curve _nettle_secp_224r1 =
     USE_REDC ? ecc_secp224r1_redc : ecc_secp224r1_modp,
     ecc_secp224r1_inv,
     NULL,
+    NULL,
   },
   {
     224,
@@ -140,6 +142,7 @@ const struct ecc_curve _nettle_secp_224r1 =
     0,
     ECC_MOD_INV_ITCH (ECC_LIMB_SIZE),
     0,
+    0,
 
     ecc_q,
     ecc_Bmodq,
@@ -151,6 +154,7 @@ const struct ecc_curve _nettle_secp_224r1 =
     ecc_mod,
     ecc_mod_inv,
     NULL,
+    NULL,
   },
   
   USE_REDC,
diff --git a/ecc-secp256r1.c b/ecc-secp256r1.c
index 625c9f7a..3bdbdc21 100644
--- a/ecc-secp256r1.c
+++ b/ecc-secp256r1.c
@@ -273,6 +273,7 @@ const struct ecc_curve _nettle_secp_256r1 =
     ECC_BMODP_SIZE,
     ECC_REDC_SIZE,
     ECC_SECP256R1_INV_ITCH,
+    0, 
     0,
 
     ecc_p,
@@ -285,6 +286,7 @@ const struct ecc_curve _nettle_secp_256r1 =
     USE_REDC ? ecc_secp256r1_redc : ecc_secp256r1_modp,
     ecc_secp256r1_inv,
     NULL,
+    NULL,
   },
   {
     256,
@@ -293,6 +295,7 @@ const struct ecc_curve _nettle_secp_256r1 =
     0,
     ECC_MOD_INV_ITCH (ECC_LIMB_SIZE),
     0,
+    0,
 
     ecc_q,
     ecc_Bmodq,
@@ -304,6 +307,7 @@ const struct ecc_curve _nettle_secp_256r1 =
     ecc_secp256r1_modq,
     ecc_mod_inv,
     NULL,
+    NULL,
   },
 
   USE_REDC,
diff --git a/ecc-secp384r1.c b/ecc-secp384r1.c
index 2c00523b..ae76390f 100644
--- a/ecc-secp384r1.c
+++ b/ecc-secp384r1.c
@@ -213,6 +213,7 @@ const struct ecc_curve _nettle_secp_384r1 =
     ECC_REDC_SIZE,
     ECC_SECP384R1_INV_ITCH,
     0,
+    0,
 
     ecc_p,
     ecc_Bmodp,
@@ -224,6 +225,7 @@ const struct ecc_curve _nettle_secp_384r1 =
     ecc_secp384r1_modp,
     ecc_secp384r1_inv,
     NULL,
+    NULL,
   },
   {
     384,
@@ -232,6 +234,7 @@ const struct ecc_curve _nettle_secp_384r1 =
     0,
     ECC_MOD_INV_ITCH (ECC_LIMB_SIZE),
     0,
+    0,
 
     ecc_q,
     ecc_Bmodq,
@@ -243,6 +246,7 @@ const struct ecc_curve _nettle_secp_384r1 =
     ecc_mod,
     ecc_mod_inv,
     NULL,
+    NULL,
   },
 
   USE_REDC,
diff --git a/ecc-secp521r1.c b/ecc-secp521r1.c
index 7eda3732..1b09e7b1 100644
--- a/ecc-secp521r1.c
+++ b/ecc-secp521r1.c
@@ -130,6 +130,7 @@ const struct ecc_curve _nettle_secp_521r1 =
     ECC_REDC_SIZE,
     ECC_SECP521R1_INV_ITCH,
     0,
+    0,
 
     ecc_p,
     ecc_Bmodp,
@@ -141,6 +142,7 @@ const struct ecc_curve _nettle_secp_521r1 =
     ecc_secp521r1_modp,
     ecc_secp521r1_inv,
     NULL,
+    NULL,
   },
   {
     521,
@@ -149,6 +151,7 @@ const struct ecc_curve _nettle_secp_521r1 =
     0,
     ECC_MOD_INV_ITCH (ECC_LIMB_SIZE),
     0,
+    0,
 
     ecc_q,
     ecc_Bmodq,
@@ -160,6 +163,7 @@ const struct ecc_curve _nettle_secp_521r1 =
     ecc_mod,
     ecc_mod_inv,
     NULL,
+    NULL,
   },
   
   USE_REDC,
diff --git a/testsuite/ecc-sqrt-test.c b/testsuite/ecc-sqrt-test.c
index 026d7f7a..69e08aa4 100644
--- a/testsuite/ecc-sqrt-test.c
+++ b/testsuite/ecc-sqrt-test.c
@@ -66,7 +66,92 @@ mpz_ui_kronecker (mp_limb_t ul, const mpz_t p)
 #endif /* NETTLE_USE_MINI_GMP */
 
 static void
-test_modulo (gmp_randstate_t rands, const struct ecc_modulo *m)
+test_sqrt (gmp_randstate_t rands, const struct ecc_modulo *m, int use_redc)
+{
+  mpz_t u;
+  mpz_t p;
+  mpz_t r;
+  mpz_t t;
+
+  unsigned z, i;
+  mp_limb_t *up;
+  mp_limb_t *rp;
+  mp_limb_t *scratch;
+
+  mpz_init (u);
+  mpz_init (t);
+
+  mpz_roinit_n (p, m->m, m->size);
+
+  up = xalloc_limbs (m->size);
+  rp = xalloc_limbs (2*m->size);
+  scratch = xalloc_limbs (m->sqrt_itch);
+
+  /* Find a non-square */
+  for (z = 2; mpz_ui_kronecker (z, p) != -1; z++)
+    ;
+
+  if (verbose)
+    fprintf(stderr, "test_sqrt on %d-bit modulo. Non square: %d\n", m->bit_size, z);
+
+  for (i = 0; i < COUNT; i++)
+    {
+      if (i & 1)
+	mpz_rrandomb (u, rands, m->bit_size);
+      else
+	mpz_urandomb (u, rands, m->bit_size);
+
+      if (use_redc)
+	{
+	  /* We get non-redc sqrt if we reduce u before calling m->sqrt */
+	  mpz_limbs_copy (scratch, u, m->size);
+	  mpn_zero (scratch + m->size, m->size);
+	  m->reduce (m, up, scratch);
+	}
+      else
+	{
+	  mpz_limbs_copy (up, u, m->size);
+	}
+      if (!m->sqrt (m, rp, up, scratch))
+	{
+	  mpz_mul_ui (u, u, z);
+	  mpz_mod (u, u, p);
+	  ecc_mod_mul_1 (m, up, up, z);
+
+	  if (!m->sqrt (m, rp, up, scratch))
+	    {
+	      fprintf (stderr, "m->sqrt returned failure, bit_size = %d\n"
+		       "u = 0x",
+		       m->bit_size);
+	      mpz_out_str (stderr, 16, u);
+	      fprintf (stderr, "\n");
+	      abort ();
+	    }
+	}
+      /* Check that r^2 = u */
+      mpz_roinit_n (r, rp, m->size);
+      mpz_mul (t, r, r);
+      if (!mpz_congruent_p (t, u, p))
+	{
+	  fprintf (stderr, "m->sqrt gave incorrect result, bit_size = %d\n"
+		   "u = 0x",
+		   m->bit_size);
+	  mpz_out_str (stderr, 16, u);
+	  fprintf (stderr, "\nr = 0x");
+	  mpz_out_str (stderr, 16, r);
+	  fprintf (stderr, "\n");
+	  abort ();
+	}
+    }
+  mpz_clear (u);
+  mpz_clear (t);
+  free (up);
+  free (rp);
+  free (scratch);
+}
+
+static void
+test_sqrt_ratio (gmp_randstate_t rands, const struct ecc_modulo *m)
 {
   mpz_t u;
   mpz_t v;
@@ -96,7 +181,7 @@ test_modulo (gmp_randstate_t rands, const struct ecc_modulo *m)
     ;
 
   if (verbose)
-    fprintf(stderr, "Non square: %d\n", z);
+    fprintf(stderr, "test_sqrt_ratio on %d-bit modulo. Non square: %d\n", m->bit_size, z);
 
   for (i = 0; i < COUNT; i++)
     {
@@ -119,7 +204,7 @@ test_modulo (gmp_randstate_t rands, const struct ecc_modulo *m)
 	  mpz_limbs_copy (up, u, m->size);
 	  if (!m->sqrt_ratio (m, rp, up, vp, scratch))
 	    {
-	      fprintf (stderr, "m->sqrt returned failure, bit_size = %d\n"
+	      fprintf (stderr, "m->sqrt_ratio returned failure, bit_size = %d\n"
 		       "u = 0x",
 		       m->bit_size);
 	      mpz_out_str (stderr, 16, u);
@@ -135,7 +220,7 @@ test_modulo (gmp_randstate_t rands, const struct ecc_modulo *m)
       mpz_mul (t, t, v);
       if (!mpz_congruent_p (t, u, p))
 	{
-	  fprintf (stderr, "m->sqrt gave incorrect result, bit_size = %d\n"
+	  fprintf (stderr, "m->sqrt_ratio gave incorrect result, bit_size = %d\n"
 		   "u = 0x",
 		   m->bit_size);
 	  mpz_out_str (stderr, 16, u);
@@ -165,8 +250,10 @@ test_main (void)
   gmp_randinit_default (rands);
   for (i = 0; ecc_curves[i]; i++)
     {
+      if (ecc_curves[i]->p.sqrt)
+	test_sqrt (rands, &ecc_curves[i]->p, ecc_curves[i]->use_redc);
       if (ecc_curves[i]->p.sqrt_ratio)
-	test_modulo (rands, &ecc_curves[i]->p);
+	test_sqrt_ratio (rands, &ecc_curves[i]->p);
     }
   gmp_randclear (rands);
 }