Use Edwards curve for curve25519 precomputation.

* eccdata.c: Reorganize curve25519 precomputation to work directly
with the twisted Edwards curve, with new point addition based on a
patch from Daiki Ueno.
* ecc-25519.c (_nettle_curve25519): Update for removed Montgomery
curve constant.

diff --git a/ChangeLog b/ChangeLog
index 3b412639861c20ae3ffd3822e0fb2da426998030..f7d28359516079f7fec6e77ea55aaf1f485593df 100644 (file)
--- a/ChangeLog
+++ b/ChangeLog
@@ -1,5 +1,11 @@
 2017-09-23  Niels Möller  <nisse@lysator.liu.se>
 
+       * eccdata.c: Reorganize curve25519 precomputation to work directly
+       with the twisted Edwards curve, with new point addition based on a
+       patch from Daiki Ueno.
+       * ecc-25519.c (_nettle_curve25519): Update for removed Montgomery
+       curve constant.
+
        * ecc-internal.h (struct ecc_curve): Delete unused pointer
        edwards_root. Update all instances.
        * eccdata.c (output_curve): Don't output it.

diff --git a/ecc-25519.c b/ecc-25519.c
index 3e008d04143daa9cd5ed85d2c21d615031879c2c..c35178c13605ecbd5849671d5fe9bf462eb9edc8 100644 (file)
--- a/ecc-25519.c
+++ b/ecc-25519.c
@@ -349,7 +349,7 @@ const struct ecc_curve _nettle_curve25519 =
   ecc_mul_g_eh,
   ecc_eh_to_a,
 
-  ecc_d, /* Use the Edwards curve constant. */
+  ecc_b, /* Edwards curve constant. */
   ecc_g,
   ecc_unit,
   ecc_table

diff --git a/eccdata.c b/eccdata.c
index b48bdd2d9b1b2c4b8941205bb4e1c234011c646d..b13547419e3ebc135d2768611ed4363ac0537a29 100644 (file)
--- a/eccdata.c
+++ b/eccdata.c
@@ -2,7 +2,9 @@
 
    Generate compile time constant (but machine dependent) tables.
 
-   Copyright (C) 2013, 2014 Niels Möller
+   Copyright (C) 2013, 2014, 2017 Niels Möller
+   Copyright (C) 2017 Daiki Ueno
+   Copyright (C) 2017 Red Hat, Inc.
 
    This file is part of GNU Nettle.
 
@@ -53,8 +55,12 @@ enum ecc_type
   {
     /* y^2 = x^3 - 3x + b (mod p) */
     ECC_TYPE_WEIERSTRASS,
-    /* y^2 = x^3 + b x^2 + x */
-    ECC_TYPE_MONTGOMERY
+#if 0
+    /* x^2 + y^2 = 1 - d x^2 y^2 */
+    ECC_TYPE_EDWARDS,
+#endif
+    /* -x^2 + y^2 = 1 - d x^2 y^2 */
+    ECC_TYPE_TWISTED_EDWARDS,
   };
 
 struct ecc_curve
@@ -73,16 +79,6 @@ struct ecc_curve
   mpz_t q;
   struct ecc_point g;
 
-  /* Non-zero if we want elements represented as point s(u, v) on an
-     equivalent Edwards curve, using
-
-      u = t x / y
-      v = (x-1) / (x+1)
-  */
-  int use_edwards;
-  mpz_t d;
-  mpz_t t;
-
   /* Table for pippenger's algorithm.
      Element
 
@@ -127,9 +123,11 @@ ecc_equal_p (const struct ecc_point *p, const struct ecc_point *q)
 }
 
 static void
-ecc_set_zero (struct ecc_point *r)
+ecc_set_zero (const struct ecc_curve *ecc, struct ecc_point *r)
 {
   r->is_zero = 1;
+  mpz_set_ui (r->x, 0);
+  mpz_set_ui (r->y, ecc->type != ECC_TYPE_WEIERSTRASS);
 }
 
 static void
@@ -140,13 +138,22 @@ ecc_set (struct ecc_point *r, const struct ecc_point *p)
   mpz_set (r->y, p->y);
 }
 
+static void
+ecc_add (const struct ecc_curve *ecc, struct ecc_point *r,
+        const struct ecc_point *p, const struct ecc_point *q);
+
 /* Needs to support in-place operation. */
 static void
 ecc_dup (const struct ecc_curve *ecc,
         struct ecc_point *r, const struct ecc_point *p)
 {
+  if (ecc->type != ECC_TYPE_WEIERSTRASS)
+    {
+      ecc_add (ecc, r, p, p);
+      return;
+    }
   if (ecc_zero_p (p))
-    ecc_set_zero (r);
+    ecc_set_zero (ecc, r);
 
   else
     {
@@ -161,32 +168,18 @@ ecc_dup (const struct ecc_curve *ecc,
       mpz_mul_ui (m, p->y, 2);
       mpz_invert (m, m, ecc->p);
 
-      switch (ecc->type)
-       {
-       case ECC_TYPE_WEIERSTRASS:
-         /* t = 3 (x^2 - 1) * m */
-         mpz_mul (t, p->x, p->x);
-         mpz_mod (t, t, ecc->p);
-         mpz_sub_ui (t, t, 1);
-         mpz_mul_ui (t, t, 3);
-         break;
-       case ECC_TYPE_MONTGOMERY:
-         /* t = (3 x^2 + 2 b x + 1) m = [x(3x+2b)+1] m */
-         mpz_mul_ui (t, ecc->b, 2);
-         mpz_addmul_ui (t, p->x, 3);
-         mpz_mul (t, t, p->x);
-         mpz_mod (t, t, ecc->p);
-         mpz_add_ui (t, t, 1);
-         break;
-       }
+      /* t = 3 (x^2 - 1) * m */
+      mpz_mul (t, p->x, p->x);
+      mpz_mod (t, t, ecc->p);
+      mpz_sub_ui (t, t, 1);
+      mpz_mul_ui (t, t, 3);
+
       mpz_mul (t, t, m);
       mpz_mod (t, t, ecc->p);
 
       /* x' = t^2 - 2 x */
       mpz_mul (x, t, t);
       mpz_submul_ui (x, p->x, 2);
-      if (ecc->type == ECC_TYPE_MONTGOMERY)
-       mpz_sub (x, x, ecc->b);
 
       mpz_mod (x, x, ecc->p);
 
@@ -208,55 +201,114 @@ ecc_dup (const struct ecc_curve *ecc,
 }
 
 static void
-ecc_add (const struct ecc_curve *ecc,
-        struct ecc_point *r, const struct ecc_point *p, const struct ecc_point *q)
+ecc_add (const struct ecc_curve *ecc, struct ecc_point *r,
+        const struct ecc_point *p, const struct ecc_point *q)
 {
-  if (ecc_zero_p (p))
-    ecc_set (r, q);
+  if (ecc->type == ECC_TYPE_WEIERSTRASS)
+    {
+      if (ecc_zero_p (p))
+       ecc_set (r, q);
 
-  else if (ecc_zero_p (q))
-    ecc_set (r, p);
+      else if (ecc_zero_p (q))
+       ecc_set (r, p);
 
-  else if (mpz_cmp (p->x, q->x) == 0)
-    {
-      if (mpz_cmp (p->y, q->y) == 0)
-       ecc_dup (ecc, r, p);
+      else if (mpz_cmp (p->x, q->x) == 0)
+       {
+         if (mpz_cmp (p->y, q->y) == 0)
+           ecc_dup (ecc, r, p);
+         else
+           ecc_set_zero (ecc, r);
+       }
       else
-       ecc_set_zero (r);
+       {
+         mpz_t s, t, x, y;
+         mpz_init (s);
+         mpz_init (t);
+         mpz_init (x);
+         mpz_init (y);
+
+         /* t = (q_y - p_y) / (q_x - p_x) */
+         mpz_sub (t, q->x, p->x);
+         mpz_invert (t, t, ecc->p);
+         mpz_sub (s, q->y, p->y);
+         mpz_mul (t, t, s);
+         mpz_mod (t, t, ecc->p);
+
+         /* x' = t^2 - p_x - q_x */
+         mpz_mul (x, t, t);
+         mpz_sub (x, x, p->x);
+         mpz_sub (x, x, q->x);
+         mpz_mod (x, x, ecc->p);
+
+         /* y' = (x - x') * t - y */
+         mpz_sub (y, p->x, x);
+         mpz_mul (y, y, t);
+         mpz_sub (y, y, p->y);
+         mpz_mod (y, y, ecc->p);
+
+         r->is_zero = 0;
+         mpz_swap (x, r->x);
+         mpz_swap (y, r->y);
+
+         mpz_clear (s);
+         mpz_clear (t);
+         mpz_clear (x);
+         mpz_clear (y);
+       }
     }
   else
     {
+      /* Untwisted:
+        x = (p_x q_y + p_y q_x) / (1 - d p_x p_y q_x q_y)
+        y = (p_y q_y - p_x q_x) / (1 + d p_x p_y q_x q_y)
+
+        Twisted:
+        x = (p_x q_y + p_y q_x) / (1 - d p_x p_y q_x q_y)
+        y = (p_y q_y + p_x q_x) / (1 + d p_x p_y q_x q_y)
+
+        So they differ only by a sign in the expression for the new y
+        coordinate.
+      */
+
       mpz_t s, t, x, y;
       mpz_init (s);
       mpz_init (t);
       mpz_init (x);
       mpz_init (y);
 
-      /* t = (q_y - p_y) / (q_x - p_x) */
-      mpz_sub (t, q->x, p->x);
-      mpz_invert (t, t, ecc->p);
-      mpz_sub (s, q->y, p->y);
-      mpz_mul (t, t, s);
+      /* t = d p_x p_y q_x q_y */
+      mpz_mul (t, ecc->b, p->x);
+      mpz_mod (t, t, ecc->p);
+      mpz_mul (t, t, p->y);
+      mpz_mod (t, t, ecc->p);
+      mpz_mul (t, t, q->x);
+      mpz_mod (t, t, ecc->p);
+      mpz_mul (t, t, q->y);
       mpz_mod (t, t, ecc->p);
 
-      /* x' = t^2 - p_x - q_x */
-      mpz_mul (x, t, t);
-      mpz_sub (x, x, p->x);
-      mpz_sub (x, x, q->x);
-      /* This appears to be the only difference between formulas. */
-      if (ecc->type == ECC_TYPE_MONTGOMERY)
-       mpz_sub (x, x, ecc->b);
+      /* x' = (p_x q_y + q_x p_y) / (1 - t) */
+      mpz_mul (x, p->x, q->y);
+      mpz_mod (x, x, ecc->p);
+      mpz_addmul (x, q->x, p->y);
+      mpz_mod (x, x, ecc->p);
+      mpz_ui_sub (s, 1, t);
+      mpz_invert (s, s, ecc->p);
+      mpz_mul (x, x, s);
       mpz_mod (x, x, ecc->p);
 
-      /* y' = (x - x') * t - y */
-      mpz_sub (y, p->x, x);
-      mpz_mul (y, y, t);
-      mpz_sub (y, y, p->y);
+      /* y' = (p_y q_y - p_x q_x) / (1 + t) */
+      mpz_mul (y, p->y, q->y);
+      mpz_mod (y, y, ecc->p);
+      mpz_addmul (y, p->x, q->x);
+      mpz_mod (y, y, ecc->p);
+      mpz_add_ui (s, t, 1);
+      mpz_invert (s, s, ecc->p);
+      mpz_mul (y, y, s);
       mpz_mod (y, y, ecc->p);
 
-      r->is_zero = 0;
       mpz_swap (x, r->x);
       mpz_swap (y, r->y);
+      r->is_zero = (mpz_cmp_ui (r->x, 0) == 0 && mpz_cmp_ui (r->y, 1) == 0);
 
       mpz_clear (s);
       mpz_clear (t);
@@ -332,16 +384,6 @@ ecc_curve_init_str (struct ecc_curve *ecc, enum ecc_type type,
   ecc->table = NULL;
 
   ecc->ref = NULL;
-
-  mpz_init (ecc->d);
-  mpz_init (ecc->t);
-
-  ecc->use_edwards = (t != NULL);
-  if (ecc->use_edwards)
-    {
-      mpz_set_str (ecc->t, t, 16);
-      mpz_set_str (ecc->d, d, 16);
-    }
 }
 
 static void
@@ -533,43 +575,37 @@ ecc_curve_init (struct ecc_curve *ecc, unsigned bit_size)
 
       break;
     case 255:
-      /* curve25519, y^2 = x^3 + 486662 x^2 + x (mod p), with p = 2^{255} - 19.
-
-        According to http://cr.yp.to/papers.html#newelliptic, this
-        is birationally equivalent to the Edwards curve
-
-          x^2 + y^2 = 1 + (121665/121666) x^2 y^2 (mod p).
+      /* Edwards curve used for eddsa25519 and curve25519,
 
-        And since the constant is not a square, the Edwards formulas
-        should be "complete", with no special cases needed for
-        doubling, neutral element, negatives, etc.
+          -x^2 + y^2 = 1 - (121665/121666) x^2 y^2, with p = 2^{255} - 19.
 
-        Generator is x = 9, with y coordinate
-        14781619447589544791020593568409986887264606134616475288964881837755586237401,
-        according to
+        The generator is
+          x = 0x216936d3cd6e53fec0a4e231fdd6dc5c692cc7609525a7b2c9562d608f25d51a
+           y = 0x6666666666666666666666666666666666666666666666666666666666666658
 
-          x = Mod(9, 2^255-19); sqrt(x^3 + 486662*x^2 + x)
+        Also birationally equivalent to the curve25519 Montgomery curve,
 
-        in PARI/GP. Also, in PARI notation,
-
-          curve25519 = Mod([0, 486662, 0, 1, 0], 2^255-19)
-       */
-      ecc_curve_init_str (ecc, ECC_TYPE_MONTGOMERY,
+          y^2 = x^3 + 486662 x^2 + x (mod p)
+      */
+      ecc_curve_init_str (ecc, ECC_TYPE_TWISTED_EDWARDS,
                          "7fffffffffffffffffffffffffffffff"
                          "ffffffffffffffffffffffffffffffed",
-                         "76d06",
+                         /* (121665/121666) mod p, from PARI/GP
+                            c = Mod(121665, p); c / (c+1)
+                         */
+                         "2dfc9311d490018c7338bf8688861767"
+                         "ff8ff5b2bebe27548a14b235eca6874a",
                          /* Order of the subgroup is 2^252 + q_0, where
                             q_0 = 27742317777372353535851937790883648493,
                             125 bits.
                          */
                          "10000000000000000000000000000000"
                          "14def9dea2f79cd65812631a5cf5d3ed",
-                         "9",
-                         /* y coordinate from PARI/GP
-                            x = Mod(9, 2^255-19); sqrt(x^3 + 486662*x^2 + x)
-                         */
-                         "20ae19a1b8a086b4e01edd2c7748d14c"
-                         "923d4d7e6d7c61b229e9c5a27eced3d9",
+                         /* Generator */
+                         "216936d3cd6e53fec0a4e231fdd6dc5c"
+                         "692cc7609525a7b2c9562d608f25d51a",
+                         "66666666666666666666666666666666"
+                         "66666666666666666666666666666658",
                          /* (121665/121666) mod p, from PARI/GP
                             c = Mod(121665, p); c / (c+1)
                          */
@@ -586,22 +622,21 @@ ecc_curve_init (struct ecc_curve *ecc, unsigned bit_size)
                          );
       ecc->ref = ecc_alloc (3);
       ecc_set_str (&ecc->ref[0], /* 2 g */
-                  "20d342d51873f1b7d9750c687d157114"
-                  "8f3f5ced1e350b5c5cae469cdd684efb",
-                  "13b57e011700e8ae050a00945d2ba2f3"
-                  "77659eb28d8d391ebcd70465c72df563");
+                  "36ab384c9f5a046c3d043b7d1833e7ac"
+                  "080d8e4515d7a45f83c5a14e2843ce0e",
+                  "2260cdf3092329c21da25ee8c9a21f56"
+                  "97390f51643851560e5f46ae6af8a3c9");
       ecc_set_str (&ecc->ref[1], /* 3 g */
-                  "1c12bc1a6d57abe645534d91c21bba64"
-                  "f8824e67621c0859c00a03affb713c12",
-                  "2986855cbe387eaeaceea446532c338c"
-                  "536af570f71ef7cf75c665019c41222b");
+                  "67ae9c4a22928f491ff4ae743edac83a"
+                  "6343981981624886ac62485fd3f8e25c",
+                  "1267b1d177ee69aba126a18e60269ef7"
+                  "9f16ec176724030402c3684878f5b4d4");
 
       ecc_set_str (&ecc->ref[2], /* 4 g */
-                  "79ce98b7e0689d7de7d1d074a15b315f"
-                  "fe1805dfcd5d2a230fee85e4550013ef",
-                  "75af5bf4ebdc75c8fe26873427d275d7"
-                  "3c0fb13da361077a565539f46de1c30");
-
+                  "203da8db56cff1468325d4b87a3520f9"
+                  "1a739ec193ce1547493aa657c4c9f870",
+                  "47d0e827cb1595e1470eb88580d5716c"
+                  "4cf22832ea2f0ff0df38ab61ca32112f");
       break;
 
     default:
@@ -624,7 +659,7 @@ ecc_pippenger_precompute (struct ecc_curve *ecc, unsigned k, unsigned c)
   ecc->table = ecc_alloc (ecc->table_size);
   
   /* Compute the first 2^c entries */
-  ecc_set_zero (&ecc->table[0]);
+  ecc_set_zero (ecc, &ecc->table[0]);
   ecc_set (&ecc->table[1], &ecc->g);
 
   for (j = 2; j < (1U<<c); j <<= 1)
@@ -658,7 +693,7 @@ ecc_mul_pippenger (const struct ecc_curve *ecc,
   mpz_init (n);
   
   mpz_mod (n, n_input, ecc->q);
-  ecc_set_zero (r);
+  ecc_set_zero (ecc, r);
 
   k = ecc->pippenger_k;
   c = ecc->pippenger_c;
@@ -856,38 +891,9 @@ output_point (const char *name, const struct ecc_curve *ecc,
   if (name)
     printf("static const mp_limb_t %s[%u] = {", name, 2*size);
 
-  if (ecc->use_edwards)
-    {
-      if (ecc_zero_p (p))
-       {
-         mpz_set_si (x, 0);
-         mpz_set_si (y, 1);
-       }
-      else if (!mpz_sgn (p->y))
-       {
-         assert (!mpz_sgn (p->x));
-         mpz_set_si (x, 0);
-         mpz_set_si (y, -1);
-       }
-      else
-       {
-         mpz_invert (x, p->y, ecc->p);
-         mpz_mul (x, x, p->x);
-         mpz_mul (x, x, ecc->t);        
-         mpz_mod (x, x, ecc->p);
+  mpz_set (x, p->x);
+  mpz_set (y, p->y);
 
-         mpz_sub_ui (y, p->x, 1);
-         mpz_add_ui (t, p->x, 1);
-         mpz_invert (t, t, ecc->p);
-         mpz_mul (y, y, t);
-         mpz_mod (y, y, ecc->p);
-       }
-    }
-  else
-    {
-      mpz_set (x, p->x);
-      mpz_set (y, p->y);
-    }
   if (use_redc)
     {
       mpz_mul_2exp (x, x, size * bits_per_limb);
@@ -945,8 +951,6 @@ output_curve (const struct ecc_curve *ecc, unsigned bits_per_limb)
 
   output_bignum ("ecc_p", ecc->p, limb_size, bits_per_limb);
   output_bignum ("ecc_b", ecc->b, limb_size, bits_per_limb);
-  if (ecc->use_edwards)
-    output_bignum ("ecc_d", ecc->d, limb_size, bits_per_limb);
   output_bignum ("ecc_q", ecc->q, limb_size, bits_per_limb);
   output_point ("ecc_g", ecc, &ecc->g, 0, limb_size, bits_per_limb);