ecc-dup-eh.c

+/* ecc-dup-eh.c
+
+   Copyright (C) 2014 Niels Möller
+
+   This file is part of GNU Nettle.
+
+   GNU Nettle is free software: you can redistribute it and/or
+   modify it under the terms of either:
+
+     * the GNU Lesser General Public License as published by the Free
+       Software Foundation; either version 3 of the License, or (at your
+       option) any later version.
+
+   or
+
+     * the GNU General Public License as published by the Free
+       Software Foundation; either version 2 of the License, or (at your
+       option) any later version.
+
+   or both in parallel, as here.
+
+   GNU Nettle is distributed in the hope that it will be useful,
+   but WITHOUT ANY WARRANTY; without even the implied warranty of
+   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+   General Public License for more details.
+
+   You should have received copies of the GNU General Public License and
+   the GNU Lesser General Public License along with this program.  If
+   not, see http://www.gnu.org/licenses/.
+*/
+
+#if HAVE_CONFIG_H
+# include "config.h"
+#endif
+
+#include "ecc.h"
+#include "ecc-internal.h"
+
+/* Double a point on an Edwards curve, in homogeneous coordinates */
+void
+ecc_dup_eh (const struct ecc_curve *ecc,
+	    mp_limb_t *r, const mp_limb_t *p,
+	    mp_limb_t *scratch)
+{
+#define x1 p
+#define y1 (p + ecc->p.size)
+#define z1 (p + 2*ecc->p.size)
+
+#define x2 r
+#define y2 (r + ecc->p.size)
+#define z2 (r + 2*ecc->p.size)
+
+  /* Formulas (from djb,
+     http://www.hyperelliptic.org/EFD/g1p/auto-edwards-projective.html#doubling-dbl-2007-bl):
+
+     Computation	Operation	Live variables
+
+     b = (x+y)^2	sqr		b
+     c = x^2		sqr		b, c
+     d = y^2		sqr		b, c, d
+     e = c+d				b, c, d, e
+     h = z^2		sqr		b, c, d, e, h
+     j = e-2*h				b, c, d, e, j
+     x' = (b-e)*j	mul		c, d, e, j
+     y' = e*(c-d)	mul		e, j
+     z' = e*j		mul
+  */
+#define C scratch
+#define D (scratch + 1*ecc->p.size)
+#define B (scratch + 2*ecc->p.size)
+
+#define E C
+
+  ecc_mod_sqr (&ecc->p, C, x1, C);	/* C */
+  ecc_mod_sqr (&ecc->p, D, y1, D);	/* C, D */
+  ecc_mod_add (&ecc->p, B, x1, y1);
+  ecc_mod_sqr (&ecc->p, B, B, x2);	/* C, D, B */
+
+  /* c-d stored at y' */
+  ecc_mod_sub (&ecc->p, y2, C, D);
+  ecc_mod_add (&ecc->p, E, C, D);	/* B, E */
+  /* b-e stored at x' */
+  ecc_mod_sub (&ecc->p, x2, B, E);	/* E */
+
+  /* Use D as scratch for the following multiplies. */
+  ecc_mod_mul (&ecc->p, y2, y2, E, D);
+
+  /* h and j stored at z' */
+  ecc_mod_sqr (&ecc->p, z2, z1, D);
+  ecc_mod_add (&ecc->p, z2, z2, z2);
+  ecc_mod_sub (&ecc->p, z2, E, z2);
+  ecc_mod_mul (&ecc->p, x2, x2, z2, D);
+  ecc_mod_mul (&ecc->p, z2, z2, E, D);
+}

ecc-dup-jj.c

+/* ecc-dup-jj.c
+
+   Copyright (C) 2013 Niels Möller
+
+   This file is part of GNU Nettle.
+
+   GNU Nettle is free software: you can redistribute it and/or
+   modify it under the terms of either:
+
+     * the GNU Lesser General Public License as published by the Free
+       Software Foundation; either version 3 of the License, or (at your
+       option) any later version.
+
+   or
+
+     * the GNU General Public License as published by the Free
+       Software Foundation; either version 2 of the License, or (at your
+       option) any later version.
+
+   or both in parallel, as here.
+
+   GNU Nettle is distributed in the hope that it will be useful,
+   but WITHOUT ANY WARRANTY; without even the implied warranty of
+   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+   General Public License for more details.
+
+   You should have received copies of the GNU General Public License and
+   the GNU Lesser General Public License along with this program.  If
+   not, see http://www.gnu.org/licenses/.
+*/
+
+/* Development of Nettle's ECC support was funded by the .SE Internet Fund. */
+
+#if HAVE_CONFIG_H
+# include "config.h"
+#endif
+
+#include "ecc.h"
+#include "ecc-internal.h"
+
+/* NOTE: Behaviour for corner cases:
+
+   + p = 0  ==>  r = 0, correct!
+*/
+void
+ecc_dup_jj (const struct ecc_curve *ecc,
+	    mp_limb_t *r, const mp_limb_t *p,
+	    mp_limb_t *scratch)
+{
+#define x1 p
+#define y1 (p + ecc->p.size)
+#define z1 (p + 2*ecc->p.size)
+
+#define x2 r
+#define y2 (r + ecc->p.size)
+#define z2 (r + 2*ecc->p.size)
+
+  /* Formulas (from djb,
+     http://www.hyperelliptic.org/EFD/g1p/auto-shortw-jacobian-3.html#doubling-dbl-2001-b):
+
+     Computation			Operation	Live variables
+     delta = z^2			sqr		delta
+     gamma = y^2			sqr		delta, gamma
+     z' = (y+z)^2-gamma-delta		sqr		delta, gamma
+     alpha = 3*(x-delta)*(x+delta)	mul		gamma, beta, alpha
+     beta = x*gamma			mul		gamma, beta, alpha
+     x' = alpha^2-8*beta		sqr		gamma, beta, alpha
+     y' = alpha*(4*beta-x')-8*gamma^2	mul, sqr
+  */
+
+#define gamma scratch
+#define delta (scratch + ecc->p.size)
+#define alpha delta
+
+#define beta (scratch + 2*ecc->p.size)
+#define sum  (scratch + 3*ecc->p.size)
+
+  ecc_mod_sqr (&ecc->p, gamma, y1, gamma);	/* x, y, z, gamma */
+  ecc_mod_sqr (&ecc->p, delta, z1, delta);	/* x, y, z, gamma, delta */
+
+  ecc_mod_add (&ecc->p, sum, z1, y1);		/* x, gamma, delta, s */
+  ecc_mod_sqr (&ecc->p, sum, sum, y2);		/* Can use y-z as scratch */
+  ecc_mod_sub (&ecc->p, z2, sum, delta);	/* x, z, gamma, delta */
+  ecc_mod_sub (&ecc->p, z2, z2, gamma);
+
+  ecc_mod_mul (&ecc->p, beta, x1, gamma, beta);	/* x, z, gamma, delta, beta */
+
+  ecc_mod_add (&ecc->p, sum, x1, delta);	/* x, sum, z', gamma, delta, beta */
+  ecc_mod_sub (&ecc->p, delta, x1, delta);	/* sum, z', gamma, delta, beta */
+  /* This multiplication peaks the storage need; can use x-y for scratch. */
+  ecc_mod_mul (&ecc->p, alpha, sum, delta, x2);	/* z', gamma, alpha, beta */
+  ecc_mod_mul_1 (&ecc->p, alpha, alpha, 3);
+
+  ecc_mod_mul_1 (&ecc->p, y2, beta, 4);
+
+  /* From now on, can use beta as scratch. */
+  ecc_mod_sqr (&ecc->p, x2, alpha, beta);	/* alpha^2 */
+  ecc_mod_submul_1 (&ecc->p, x2, y2, 2);	/*  alpha^2 - 8 beta */
+
+  ecc_mod_sub (&ecc->p, y2, y2, x2);		/* 4 beta - x' */
+  ecc_mod_mul (&ecc->p, y2, y2, alpha, beta);
+  ecc_mod_sqr (&ecc->p, gamma, gamma, beta);
+  ecc_mod_submul_1 (&ecc->p, y2, gamma, 8);
+}

ecc-dup-th.c

+/* ecc-dup-th.c
+
+   Copyright (C) 2014, 2019 Niels Möller
+
+   This file is part of GNU Nettle.
+
+   GNU Nettle is free software: you can redistribute it and/or
+   modify it under the terms of either:
+
+     * the GNU Lesser General Public License as published by the Free
+       Software Foundation; either version 3 of the License, or (at your
+       option) any later version.
+
+   or
+
+     * the GNU General Public License as published by the Free
+       Software Foundation; either version 2 of the License, or (at your
+       option) any later version.
+
+   or both in parallel, as here.
+
+   GNU Nettle is distributed in the hope that it will be useful,
+   but WITHOUT ANY WARRANTY; without even the implied warranty of
+   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+   General Public License for more details.
+
+   You should have received copies of the GNU General Public License and
+   the GNU Lesser General Public License along with this program.  If
+   not, see http://www.gnu.org/licenses/.
+*/
+
+#if HAVE_CONFIG_H
+# include "config.h"
+#endif
+
+#include "ecc.h"
+#include "ecc-internal.h"
+
+/* Double a point on a twisted Edwards curve, in homogeneous coordinates */
+void
+ecc_dup_th (const struct ecc_curve *ecc,
+	    mp_limb_t *r, const mp_limb_t *p,
+	    mp_limb_t *scratch)
+{
+#define x1 p
+#define y1 (p + ecc->p.size)
+#define z1 (p + 2*ecc->p.size)
+
+#define x2 r
+#define y2 (r + ecc->p.size)
+#define z2 (r + 2*ecc->p.size)
+
+  /* Formulas (from djb,
+     http://www.hyperelliptic.org/EFD/g1p/auto-twisted-projective.html#doubling-dbl-2008-bbjlp):
+
+     B = (X1+Y1)^2
+     C = X1^2
+     D = Y1^2
+     (E = a*C = -C)
+     F = E+D
+     H = Z1^2
+     J = F-2*H
+     X3 = (B-C-D)*J
+     Y3 = F*(E-D)
+     Z3 = F*J         (-C+D)*(-C+D - 2Z1^2)
+
+     In the formula for Y3, we have E - D = -(C+D). To avoid explicit
+     negation, negate all of X3, Y3, Z3, and use
+
+     Computation	Operation	Live variables
+
+     B = (X1+Y1)^2	sqr		B
+     C = X1^2		sqr		B, C
+     D = Y1^2		sqr		B, C, D
+     F = -C+D				B, C, D, F
+     H = Z1^2		sqr		B, C, D, F, H
+     J = 2*H - F			B, C, D, F, J
+     X3 = (B-C-D)*J	mul		C, F, J  (Replace C <-- C+D)
+     Y3 = F*(C+D)	mul		F, J
+     Z3 = F*J		mul
+
+     3M+4S
+  */
+
+#define C scratch
+#define D (scratch + 1*ecc->p.size)
+#define B (scratch + 2*ecc->p.size)
+
+#define F C
+
+  ecc_mod_sqr (&ecc->p, C, x1, C);	/* C */
+  ecc_mod_sqr (&ecc->p, D, y1, D);	/* C, D */
+  ecc_mod_add (&ecc->p, B, x1, y1);
+  ecc_mod_sqr (&ecc->p, B, B, x2);	/* C, D, B */
+
+  /* C+D stored at y' */
+  ecc_mod_add (&ecc->p, y2, C, D);
+  /* B - C - C stored at x' */
+  ecc_mod_sub (&ecc->p, x2, B, y2);
+
+  ecc_mod_sub (&ecc->p, F, D, C);	/* F */
+
+  /* Use D as scratch for the following multiplies. */
+  ecc_mod_mul (&ecc->p, y2, y2, F, D);
+
+  /* H and J stored at z' */
+  ecc_mod_sqr (&ecc->p, z2, z1, D);
+  ecc_mod_add (&ecc->p, z2, z2, z2);
+  ecc_mod_sub (&ecc->p, z2, z2, F);
+  ecc_mod_mul (&ecc->p, x2, x2, z2, D);
+  ecc_mod_mul (&ecc->p, z2, z2, F, D);
+}

ecc-ecdsa-sign.c

+/* ecc-ecdsa-sign.c
+
+   Copyright (C) 2013, 2014 Niels Möller
+
+   This file is part of GNU Nettle.
+
+   GNU Nettle is free software: you can redistribute it and/or
+   modify it under the terms of either:
+
+     * the GNU Lesser General Public License as published by the Free
+       Software Foundation; either version 3 of the License, or (at your
+       option) any later version.
+
+   or
+
+     * the GNU General Public License as published by the Free
+       Software Foundation; either version 2 of the License, or (at your
+       option) any later version.
+
+   or both in parallel, as here.
+
+   GNU Nettle is distributed in the hope that it will be useful,
+   but WITHOUT ANY WARRANTY; without even the implied warranty of
+   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+   General Public License for more details.
+
+   You should have received copies of the GNU General Public License and
+   the GNU Lesser General Public License along with this program.  If
+   not, see http://www.gnu.org/licenses/.
+*/
+
+/* Development of Nettle's ECC support was funded by the .SE Internet Fund. */
+
+#if HAVE_CONFIG_H
+# include "config.h"
+#endif
+
+#include <assert.h>
+#include <stdlib.h>
+
+#include "ecdsa.h"
+#include "ecc-internal.h"
+#include "dsa-internal.h"
+
+/* Low-level ECDSA signing */
+
+mp_size_t
+ecc_ecdsa_sign_itch (const struct ecc_curve *ecc)
+{
+  /* Needs 3*ecc->p.size + scratch for ecc_mul_g. */
+  assert (ecc->p.size + ecc->p.invert_itch
+	  <= 3*ecc->p.size + ECC_MUL_G_ITCH (ecc->p.size));
+  return ECC_ECDSA_SIGN_ITCH (ecc->p.size);
+}
+
+/* NOTE: Caller should check if r or s is zero. */
+void
+ecc_ecdsa_sign (const struct ecc_curve *ecc,
+		const mp_limb_t *zp,
+		/* Random nonce, must be invertible mod ecc group
+		   order. */
+		const mp_limb_t *kp,
+		size_t length, const uint8_t *digest,
+		mp_limb_t *rp, mp_limb_t *sp,
+		mp_limb_t *scratch)
+{
+#define P	    scratch
+#define kinv	    scratch
+#define hp	    (scratch  + ecc->p.size) /* NOTE: ecc->p.size + 1 limbs! */
+#define tp	    (scratch + 2*ecc->p.size)
+  /* Procedure, according to RFC 6090, "KT-I". q denotes the group
+     order.
+
+     1. k <-- uniformly random, 0 < k < q
+
+     2. R <-- (r_x, r_y) = k g
+
+     3. s1 <-- r_x mod q
+
+     4. s2 <-- (h + z*s1)/k mod q.
+  */
+
+  ecc_mul_g (ecc, P, kp, P + 3*ecc->p.size);
+  /* x coordinate only, modulo q */
+  ecc_j_to_a (ecc, 2, rp, P, P + 3*ecc->p.size);
+
+  /* Invert k, uses up to 7 * ecc->p.size including scratch (for secp384). */
+  ecc->q.invert (&ecc->q, kinv, kp, tp);
+  
+  /* Process hash digest */
+  _nettle_dsa_hash (hp, ecc->q.bit_size, length, digest);
+
+  ecc_mod_mul (&ecc->q, tp, zp, rp, tp);
+  ecc_mod_add (&ecc->q, hp, hp, tp);
+  ecc_mod_mul_canonical (&ecc->q, sp, hp, kinv, tp);
+
+#undef P
+#undef hp
+#undef kinv
+#undef tp
+}

ecc-ecdsa-verify.c

+/* ecc-ecdsa-verify.c
+
+   Copyright (C) 2013, 2014 Niels Möller
+
+   This file is part of GNU Nettle.
+
+   GNU Nettle is free software: you can redistribute it and/or
+   modify it under the terms of either:
+
+     * the GNU Lesser General Public License as published by the Free
+       Software Foundation; either version 3 of the License, or (at your
+       option) any later version.
+
+   or
+
+     * the GNU General Public License as published by the Free
+       Software Foundation; either version 2 of the License, or (at your
+       option) any later version.
+
+   or both in parallel, as here.
+
+   GNU Nettle is distributed in the hope that it will be useful,
+   but WITHOUT ANY WARRANTY; without even the implied warranty of
+   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+   General Public License for more details.
+
+   You should have received copies of the GNU General Public License and
+   the GNU Lesser General Public License along with this program.  If
+   not, see http://www.gnu.org/licenses/.
+*/
+
+/* Development of Nettle's ECC support was funded by the .SE Internet Fund. */
+
+#if HAVE_CONFIG_H
+# include "config.h"
+#endif
+
+#include <assert.h>
+#include <stdlib.h>
+
+#include "ecdsa.h"
+#include "ecc-internal.h"
+#include "dsa-internal.h"
+
+/* Low-level ECDSA verify */
+
+static int
+ecdsa_in_range (const struct ecc_curve *ecc, const mp_limb_t *xp)
+{
+  return !mpn_zero_p (xp, ecc->p.size)
+    && mpn_cmp (xp, ecc->q.m, ecc->p.size) < 0;
+}
+
+mp_size_t
+ecc_ecdsa_verify_itch (const struct ecc_curve *ecc)
+{
+  /* Largest storage need is for the ecc_mul_a call. */
+  return 5*ecc->p.size + ECC_MUL_A_ITCH (ecc->p.size);
+}
+
+/* FIXME: Use faster primitives, not requiring side-channel silence. */
+int
+ecc_ecdsa_verify (const struct ecc_curve *ecc,
+		  const mp_limb_t *pp, /* Public key */
+		  size_t length, const uint8_t *digest,
+		  const mp_limb_t *rp, const mp_limb_t *sp,
+		  mp_limb_t *scratch)
+{
+  /* Procedure, according to RFC 6090, "KT-I". q denotes the group
+     order.
+
+     1. Check 0 < r, s < q.
+
+     2. s' <-- s^{-1}  (mod q)
+
+     3. u1  <-- h * s' (mod q)
+
+     4. u2  <-- r * s' (mod q)
+
+     5. R = u1 G + u2 Y
+
+     6. Signature is valid if R_x = r (mod q).
+  */
+
+#define P2 scratch
+#define u1 (scratch + 3*ecc->p.size)
+#define u2 (scratch + 4*ecc->p.size)
+
+#define P1 (scratch + 4*ecc->p.size)
+#define sinv (scratch)
+#define hp (scratch + ecc->p.size)
+
+  if (! (ecdsa_in_range (ecc, rp)
+	 && ecdsa_in_range (ecc, sp)))
+    return 0;
+
+  /* FIXME: Micro optimizations: Either simultaneous multiplication.
+     Or convert to projective coordinates (can be done without
+     division, I think), and write an ecc_add_ppp. */
+
+  /* Compute sinv */
+  ecc->q.invert (&ecc->q, sinv, sp, sinv + ecc->p.size);
+
+  /* u1 = h / s, P1 = u1 * G */
+  _nettle_dsa_hash (hp, ecc->q.bit_size, length, digest);
+  ecc_mod_mul_canonical (&ecc->q, u1, hp, sinv, u1);
+
+  /* u2 = r / s, P2 = u2 * Y */
+  ecc_mod_mul_canonical (&ecc->q, u2, rp, sinv, u2);
+
+   /* Total storage: 5*ecc->p.size + ECC_MUL_A_ITCH */
+  ecc_mul_a (ecc, P2, u2, pp, u2 + ecc->p.size);
+
+  /* u = 0 can happen only if h = 0 or h = q, which is extremely
+     unlikely. */
+  if (!mpn_zero_p (u1, ecc->p.size))
+    {
+      /* Total storage: 7*ecc->p.size + ECC_MUL_G_ITCH */
+      ecc_mul_g (ecc, P1, u1, P1 + 3*ecc->p.size);
+
+      /* Total storage: 6*ecc->p.size + ECC_ADD_JJJ_ITCH */
+      if (!ecc_nonsec_add_jjj (ecc, P2, P2, P1, P1 + 3*ecc->p.size))
+	/* Infinity point, not a valid signature. */
+	return 0;
+    }
+  /* x coordinate only, modulo q */
+  ecc_j_to_a (ecc, 2, P1, P2, P1 + 3*ecc->p.size);
+
+  return (mpn_cmp (rp, P1, ecc->p.size) == 0);
+#undef P2
+#undef P1
+#undef sinv
+#undef u2
+#undef hp
+#undef u1
+}

ecc-eh-to-a.c

+/* ecc-eh-to-a.c
+
+   Copyright (C) 2014 Niels Möller
+
+   This file is part of GNU Nettle.
+
+   GNU Nettle is free software: you can redistribute it and/or
+   modify it under the terms of either:
+
+     * the GNU Lesser General Public License as published by the Free
+       Software Foundation; either version 3 of the License, or (at your
+       option) any later version.
+
+   or
+
+     * the GNU General Public License as published by the Free
+       Software Foundation; either version 2 of the License, or (at your
+       option) any later version.
+
+   or both in parallel, as here.
+
+   GNU Nettle is distributed in the hope that it will be useful,
+   but WITHOUT ANY WARRANTY; without even the implied warranty of
+   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+   General Public License for more details.
+
+   You should have received copies of the GNU General Public License and
+   the GNU Lesser General Public License along with this program.  If
+   not, see http://www.gnu.org/licenses/.
+*/
+
+#if HAVE_CONFIG_H
+# include "config.h"
+#endif
+
+#include <assert.h>
+
+#include "ecc.h"
+#include "ecc-internal.h"
+
+/* Convert from homogeneous coordinates on the Edwards curve to affine
+   coordinates. */
+void
+ecc_eh_to_a (const struct ecc_curve *ecc,
+	     int op,
+	     mp_limb_t *r, const mp_limb_t *p,
+	     mp_limb_t *scratch)
+{
+#define izp scratch
+#define tp (scratch + ecc->p.size)
+
+#define xp p
+#define yp (p + ecc->p.size)
+#define zp (p + 2*ecc->p.size)
+
+  assert(op == 0);
+
+  /* Needs size + scratch for the invert call. */
+  ecc->p.invert (&ecc->p, izp, zp, tp);
+
+  ecc_mod_mul_canonical (&ecc->p, r, xp, izp, tp);
+  ecc_mod_mul_canonical (&ecc->p, r + ecc->p.size, yp, izp, tp);
+}

ecc-gost-gc256b.c

+/* ecc-gost-gc256b.c
+
+   Copyright (C) 2016-2020 Dmitry Eremin-Solenikov
+
+   This file is part of GNU Nettle.
+
+   GNU Nettle is free software: you can redistribute it and/or
+   modify it under the terms of either:
+
+     * the GNU Lesser General Public License as published by the Free
+       Software Foundation; either version 3 of the License, or (at your
+       option) any later version.
+
+   or
+
+     * the GNU General Public License as published by the Free
+       Software Foundation; either version 2 of the License, or (at your
+       option) any later version.
+
+   or both in parallel, as here.
+
+   GNU Nettle is distributed in the hope that it will be useful,
+   but WITHOUT ANY WARRANTY; without even the implied warranty of
+   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+   General Public License for more details.
+
+   You should have received copies of the GNU General Public License and
+   the GNU Lesser General Public License along with this program.  If
+   not, see http://www.gnu.org/licenses/.
+*/
+
+#if HAVE_CONFIG_H
+# include "config.h"
+#endif
+
+#include <assert.h>
+
+#include "ecc.h"
+#include "ecc-internal.h"
+
+#define USE_REDC 0
+
+#include "ecc-gost-gc256b.h"
+
+static void
+ecc_gost_gc256b_modp (const struct ecc_modulo *m, mp_limb_t *rp, mp_limb_t *xp)
+{
+  mp_size_t mn = m->size;
+  mp_limb_t hi;
+
+  hi = mpn_addmul_1(xp, xp + mn, mn, 0x269);
+  hi = sec_add_1 (xp, xp, mn, hi * 0x269);
+  hi = sec_add_1 (rp, xp, mn, hi * 0x269);
+  assert(hi == 0);
+}
+
+#define ecc_gost_gc256b_modp ecc_gost_gc256b_modp
+#define ecc_gost_gc256b_modq ecc_mod
+
+const struct ecc_curve _nettle_gost_gc256b =
+{
+  {
+    256,
+    ECC_LIMB_SIZE,
+    ECC_BMODP_SIZE,
+    ECC_REDC_SIZE,
+    ECC_MOD_INV_ITCH (ECC_LIMB_SIZE),
+    0,
+    0,
+
+    ecc_p,
+    ecc_Bmodp,
+    ecc_Bmodp_shifted,
+    ecc_Bm2p,
+    ecc_redc_ppm1,
+
+    ecc_pp1h,
+    ecc_gost_gc256b_modp,
+    ecc_gost_gc256b_modp,
+    ecc_mod_inv,
+    NULL,
+    NULL,
+  },
+  {
+    256,
+    ECC_LIMB_SIZE,
+    ECC_BMODQ_SIZE,
+    0,
+    ECC_MOD_INV_ITCH (ECC_LIMB_SIZE),
+    0,
+    0,
+
+    ecc_q,
+    ecc_Bmodq,
+    ecc_Bmodq_shifted,
+    ecc_Bm2q,
+    NULL,
+    ecc_qp1h,
+
+    ecc_gost_gc256b_modq,
+    ecc_gost_gc256b_modq,
+    ecc_mod_inv,
+    NULL,
+    NULL,
+  },
+
+  USE_REDC,
+  ECC_PIPPENGER_K,
+  ECC_PIPPENGER_C,
+
+  ECC_ADD_JJA_ITCH (ECC_LIMB_SIZE),
+  ECC_ADD_JJJ_ITCH (ECC_LIMB_SIZE),
+  ECC_DUP_JJ_ITCH (ECC_LIMB_SIZE),
+  ECC_MUL_A_ITCH (ECC_LIMB_SIZE),
+  ECC_MUL_G_ITCH (ECC_LIMB_SIZE),
+  ECC_J_TO_A_ITCH (ECC_LIMB_SIZE, ECC_MOD_INV_ITCH(ECC_LIMB_SIZE)),
+
+  ecc_add_jja,
+  ecc_add_jjj,
+  ecc_dup_jj,
+  ecc_mul_a,
+  ecc_mul_g,
+  ecc_j_to_a,
+
+  ecc_b,
+  ecc_unit,
+  ecc_table
+};
+
+const struct ecc_curve *nettle_get_gost_gc256b(void)
+{
+  return &_nettle_gost_gc256b;
+}

ecc-gost-gc512a.c

+/* ecc-gost-gc512a.c
+
+   Copyright (C) 2016-2020 Dmitry Eremin-Solenikov
+
+   This file is part of GNU Nettle.
+
+   GNU Nettle is free software: you can redistribute it and/or
+   modify it under the terms of either:
+
+     * the GNU Lesser General Public License as published by the Free
+       Software Foundation; either version 3 of the License, or (at your
+       option) any later version.
+
+   or
+
+     * the GNU General Public License as published by the Free
+       Software Foundation; either version 2 of the License, or (at your
+       option) any later version.
+
+   or both in parallel, as here.
+
+   GNU Nettle is distributed in the hope that it will be useful,
+   but WITHOUT ANY WARRANTY; without even the implied warranty of
+   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+   General Public License for more details.
+
+   You should have received copies of the GNU General Public License and
+   the GNU Lesser General Public License along with this program.  If
+   not, see http://www.gnu.org/licenses/.
+*/
+
+#if HAVE_CONFIG_H
+# include "config.h"
+#endif
+
+#include <assert.h>
+
+#include "ecc.h"
+#include "ecc-internal.h"
+
+#define USE_REDC 0
+
+#include "ecc-gost-gc512a.h"
+
+static void
+ecc_gost_gc512a_modp (const struct ecc_modulo *m, mp_limb_t *rp, mp_limb_t *xp)
+{
+  mp_size_t mn = m->size;
+  mp_limb_t hi;
+
+  hi = mpn_addmul_1(xp, xp + mn, mn, 0x239);
+  hi = sec_add_1 (xp, xp, mn, hi * 0x239);
+  hi = sec_add_1 (rp, xp, mn, hi * 0x239);
+  assert(hi == 0);
+}
+
+#define ecc_gost_gc512a_modp ecc_gost_gc512a_modp
+#define ecc_gost_gc512a_modq ecc_mod
+
+const struct ecc_curve _nettle_gost_gc512a =
+{
+  {
+    512,
+    ECC_LIMB_SIZE,
+    ECC_BMODP_SIZE,
+    ECC_REDC_SIZE,
+    ECC_MOD_INV_ITCH (ECC_LIMB_SIZE),
+    0,
+    0,
+
+    ecc_p,
+    ecc_Bmodp,
+    ecc_Bmodp_shifted,
+    ecc_Bm2p,
+    ecc_redc_ppm1,
+
+    ecc_pp1h,
+    ecc_gost_gc512a_modp,
+    ecc_gost_gc512a_modp,
+    ecc_mod_inv,
+    NULL,
+    NULL,
+  },
+  {
+    512,
+    ECC_LIMB_SIZE,
+    ECC_BMODQ_SIZE,
+    0,
+    ECC_MOD_INV_ITCH (ECC_LIMB_SIZE),
+    0,
+    0,
+
+    ecc_q,
+    ecc_Bmodq,
+    ecc_Bmodq_shifted,
+    ecc_Bm2q,
+    NULL,
+    ecc_qp1h,
+
+    ecc_gost_gc512a_modq,
+    ecc_gost_gc512a_modq,
+    ecc_mod_inv,
+    NULL,
+    NULL,
+  },
+
+  USE_REDC,
+  ECC_PIPPENGER_K,
+  ECC_PIPPENGER_C,
+
+  ECC_ADD_JJA_ITCH (ECC_LIMB_SIZE),
+  ECC_ADD_JJJ_ITCH (ECC_LIMB_SIZE),
+  ECC_DUP_JJ_ITCH (ECC_LIMB_SIZE),
+  ECC_MUL_A_ITCH (ECC_LIMB_SIZE),
+  ECC_MUL_G_ITCH (ECC_LIMB_SIZE),
+  ECC_J_TO_A_ITCH (ECC_LIMB_SIZE, ECC_MOD_INV_ITCH(ECC_LIMB_SIZE)),
+
+  ecc_add_jja,
+  ecc_add_jjj,
+  ecc_dup_jj,
+  ecc_mul_a,
+  ecc_mul_g,
+  ecc_j_to_a,
+
+  ecc_b,
+  ecc_unit,
+  ecc_table
+};
+
+const struct ecc_curve *nettle_get_gost_gc512a(void)
+{
+  return &_nettle_gost_gc512a;
+}

ecc-gostdsa-sign.c

+/* ecc-gostdsa-sign.c
+
+   Copyright (C) 2015 Dmitry Eremin-Solenikov
+   Copyright (C) 2013, 2014 Niels Möller
+
+   This file is part of GNU Nettle.
+
+   GNU Nettle is free software: you can redistribute it and/or
+   modify it under the terms of either:
+
+     * the GNU Lesser General Public License as published by the Free
+       Software Foundation; either version 3 of the License, or (at your
+       option) any later version.
+
+   or
+
+     * the GNU General Public License as published by the Free
+       Software Foundation; either version 2 of the License, or (at your
+       option) any later version.
+
+   or both in parallel, as here.
+
+   GNU Nettle is distributed in the hope that it will be useful,
+   but WITHOUT ANY WARRANTY; without even the implied warranty of
+   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+   General Public License for more details.
+
+   You should have received copies of the GNU General Public License and
+   the GNU Lesser General Public License along with this program.  If
+   not, see http://www.gnu.org/licenses/.
+*/
+
+#if HAVE_CONFIG_H
+# include "config.h"
+#endif
+
+#include <assert.h>
+#include <stdlib.h>
+
+#include "gostdsa.h"
+#include "dsa-internal.h"
+#include "ecc-internal.h"
+
+/* Low-level GOST DSA signing */
+
+mp_size_t
+ecc_gostdsa_sign_itch (const struct ecc_curve *ecc)
+{
+  /* Needs 3*ecc->p.size + scratch for ecc_mul_g. */
+  return ECC_GOSTDSA_SIGN_ITCH (ecc->p.size);
+}
+
+/* NOTE: Caller should check if r or s is zero. */
+void
+ecc_gostdsa_sign (const struct ecc_curve *ecc,
+		const mp_limb_t *zp,
+		const mp_limb_t *kp,
+		size_t length, const uint8_t *digest,
+		mp_limb_t *rp, mp_limb_t *sp,
+		mp_limb_t *scratch)
+{
+#define P	    scratch
+#define hp	    (scratch + 4*ecc->p.size)
+#define tp	    (scratch + 2*ecc->p.size)
+#define t2p	    scratch
+  /* Procedure, according to GOST 34.10. q denotes the group
+     order.
+
+     1. k <-- uniformly random, 0 < k < q
+
+     2. C <-- (c_x, c_y) = k g
+
+     3. r <-- c_x mod q
+
+     4. s <-- (r*z + k*h) mod q.
+  */
+
+  ecc_mul_g (ecc, P, kp, P + 3*ecc->p.size);
+  /* x coordinate only, modulo q */
+  ecc_j_to_a (ecc, 2, rp, P, P + 3*ecc->p.size);
+
+  /* Process hash digest */
+  _nettle_gostdsa_hash (hp, ecc->q.bit_size, length, digest);
+  if (mpn_zero_p (hp, ecc->p.size))
+    mpn_add_1 (hp, hp, ecc->p.size, 1);
+
+  ecc_mod_mul (&ecc->q, tp, rp, zp, tp);
+  ecc_mod_mul (&ecc->q, t2p, kp, hp, t2p);
+  ecc_mod_add (&ecc->q, sp, tp, t2p);
+
+  /* Also reduce mod ecc->q. It should already be < 2*ecc->q,
+   * so one subtraction should suffice. */
+
+  *scratch = mpn_sub_n (tp, sp, ecc->q.m, ecc->p.size);
+  cnd_copy (is_zero_limb (*scratch), sp, tp, ecc->p.size);
+
+#undef P
+#undef hp
+#undef tp
+#undef t2p
+}

ecc-gostdsa-verify.c

+/* ecc-gostdsa-verify.c
+
+   Copyright (C) 2015 Dmitry Eremin-Solenikov
+   Copyright (C) 2013, 2014 Niels Möller
+
+   This file is part of GNU Nettle.
+
+   GNU Nettle is free software: you can redistribute it and/or
+   modify it under the terms of either:
+
+     * the GNU Lesser General Public License as published by the Free
+       Software Foundation; either version 3 of the License, or (at your
+       option) any later version.
+
+   or
+
+     * the GNU General Public License as published by the Free
+       Software Foundation; either version 2 of the License, or (at your
+       option) any later version.
+
+   or both in parallel, as here.
+
+   GNU Nettle is distributed in the hope that it will be useful,
+   but WITHOUT ANY WARRANTY; without even the implied warranty of
+   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+   General Public License for more details.
+
+   You should have received copies of the GNU General Public License and
+   the GNU Lesser General Public License along with this program.  If
+   not, see http://www.gnu.org/licenses/.
+*/
+
+#if HAVE_CONFIG_H
+# include "config.h"
+#endif
+
+#include <assert.h>
+#include <stdlib.h>
+
+#include "gostdsa.h"
+#include "dsa-internal.h"
+#include "ecc-internal.h"
+
+/* Low-level GOST DSA verify */
+
+static int
+ecdsa_in_range (const struct ecc_curve *ecc, const mp_limb_t *xp)
+{
+  return !mpn_zero_p (xp, ecc->p.size)
+    && mpn_cmp (xp, ecc->q.m, ecc->p.size) < 0;
+}
+
+mp_size_t
+ecc_gostdsa_verify_itch (const struct ecc_curve *ecc)
+{
+  /* Largest storage need is for the ecc_mul_a call. */
+  return 5*ecc->p.size + ECC_MUL_A_ITCH (ecc->p.size);
+}
+
+/* FIXME: Use faster primitives, not requiring side-channel silence. */
+int
+ecc_gostdsa_verify (const struct ecc_curve *ecc,
+		  const mp_limb_t *pp, /* Public key */
+		  size_t length, const uint8_t *digest,
+		  const mp_limb_t *rp, const mp_limb_t *sp,
+		  mp_limb_t *scratch)
+{
+  /* Procedure, according to GOST R 34.10. q denotes the group
+     order.
+
+     1. Check 0 < r, s < q.
+
+     2. v <-- h^{-1}  (mod q)
+
+     3. z1  <-- s * v (mod q)
+
+     4. z2  <-- -r * v (mod q)
+
+     5. R = u1 G + u2 Y
+
+     6. Signature is valid if R_x = r (mod q).
+  */
+
+#define hp (scratch)
+#define vp (scratch + ecc->p.size)
+#define z1 (scratch + 3*ecc->p.size)
+#define z2 (scratch + 4*ecc->p.size)
+
+#define P1 (scratch + 4*ecc->p.size)
+#define P2 (scratch)
+
+
+  if (! (ecdsa_in_range (ecc, rp)
+	 && ecdsa_in_range (ecc, sp)))
+    return 0;
+
+  _nettle_gostdsa_hash (hp, ecc->q.bit_size, length, digest);
+
+  if (mpn_zero_p (hp, ecc->p.size))
+    mpn_add_1 (hp, hp, ecc->p.size, 1);
+
+  /* Compute v */
+  ecc->q.invert (&ecc->q, vp, hp, vp + ecc->p.size);
+
+  /* z1 = s / h, P1 = z1 * G */
+  ecc_mod_mul_canonical (&ecc->q, z1, sp, vp, z1);
+
+  /* z2 = - r / h, P2 = z2 * Y */
+  mpn_sub_n (hp, ecc->q.m, rp, ecc->p.size);
+  ecc_mod_mul_canonical (&ecc->q, z2, hp, vp, z2);
+
+   /* Total storage: 5*ecc->p.size + ECC_MUL_A_ITCH */
+  ecc_mul_a (ecc, P2, z2, pp, z2 + ecc->p.size);
+
+  /* Total storage: 7*ecc->p.size + ECC_MUL_G_ITCH */
+  ecc_mul_g (ecc, P1, z1, P1 + 3*ecc->p.size);
+
+  /* Total storage: 6*ecc->p.size + ECC_ADD_JJJ_ITCH */
+  if (!ecc_nonsec_add_jjj (ecc, P1, P1, P2, P1 + 3*ecc->p.size))
+    return 0;
+
+  /* x coordinate only, modulo q */
+  ecc_j_to_a (ecc, 2, P2, P1, P1 + 3*ecc->p.size);
+
+  return (mpn_cmp (rp, P2, ecc->p.size) == 0);
+#undef P2
+#undef P1
+#undef z2
+#undef z1
+#undef hp
+#undef vp
+}

ecc-internal.h

+/* ecc-internal.h
+
+   Copyright (C) 2013, 2014 Niels Möller
+
+   This file is part of GNU Nettle.
+
+   GNU Nettle is free software: you can redistribute it and/or
+   modify it under the terms of either:
+
+     * the GNU Lesser General Public License as published by the Free
+       Software Foundation; either version 3 of the License, or (at your
+       option) any later version.
+
+   or
+
+     * the GNU General Public License as published by the Free
+       Software Foundation; either version 2 of the License, or (at your
+       option) any later version.
+
+   or both in parallel, as here.
+
+   GNU Nettle is distributed in the hope that it will be useful,
+   but WITHOUT ANY WARRANTY; without even the implied warranty of
+   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+   General Public License for more details.
+
+   You should have received copies of the GNU General Public License and
+   the GNU Lesser General Public License along with this program.  If
+   not, see http://www.gnu.org/licenses/.
+*/
+
+/* Development of Nettle's ECC support was funded by the .SE Internet Fund. */
+
+#ifndef NETTLE_ECC_INTERNAL_H_INCLUDED
+#define NETTLE_ECC_INTERNAL_H_INCLUDED
+
+#include "nettle-types.h"
+#include "bignum.h"
+#include "ecc-curve.h"
+#include "gmp-glue.h"
+
+/* Name mangling */
+#define ecc_pp1_redc _nettle_ecc_pp1_redc
+#define ecc_pm1_redc _nettle_ecc_pm1_redc
+#define ecc_mod_zero_p _nettle_ecc_mod_zero_p
+#define ecc_mod_equal_p _nettle_ecc_mod_equal_p
+#define ecc_mod_add _nettle_ecc_mod_add
+#define ecc_mod_sub _nettle_ecc_mod_sub
+#define ecc_mod_mul_1 _nettle_ecc_mod_mul_1
+#define ecc_mod_addmul_1 _nettle_ecc_mod_addmul_1
+#define ecc_mod_submul_1 _nettle_ecc_mod_submul_1
+#define ecc_mod_mul _nettle_ecc_mod_mul
+#define ecc_mod_sqr _nettle_ecc_mod_sqr
+#define ecc_mod_mul_canonical _nettle_ecc_mod_mul_canonical
+#define ecc_mod_sqr_canonical _nettle_ecc_mod_sqr_canonical
+#define ecc_mod_pow_2k _nettle_ecc_mod_pow_2k
+#define ecc_mod_pow_2k_mul _nettle_ecc_mod_pow_2k_mul
+#define ecc_mod_random _nettle_ecc_mod_random
+#define ecc_mod _nettle_ecc_mod
+#define ecc_mod_inv _nettle_ecc_mod_inv
+#define ecc_a_to_j _nettle_ecc_a_to_j
+#define ecc_j_to_a _nettle_ecc_j_to_a
+#define ecc_eh_to_a _nettle_ecc_eh_to_a
+#define ecc_dup_jj _nettle_ecc_dup_jj
+#define ecc_add_jja _nettle_ecc_add_jja
+#define ecc_add_jjj _nettle_ecc_add_jjj
+#define ecc_nonsec_add_jjj _nettle_ecc_nonsec_add_jjj
+#define ecc_dup_eh _nettle_ecc_dup_eh
+#define ecc_add_eh _nettle_ecc_add_eh
+#define ecc_add_ehh _nettle_ecc_add_ehh
+#define ecc_dup_th _nettle_ecc_dup_th
+#define ecc_add_th _nettle_ecc_add_th
+#define ecc_add_thh _nettle_ecc_add_thh
+#define ecc_mul_g _nettle_ecc_mul_g
+#define ecc_mul_a _nettle_ecc_mul_a
+#define ecc_mul_g_eh _nettle_ecc_mul_g_eh
+#define ecc_mul_a_eh _nettle_ecc_mul_a_eh
+#define ecc_mul_m _nettle_ecc_mul_m
+#define cnd_copy _nettle_cnd_copy
+#define sec_add_1 _nettle_sec_add_1
+#define sec_sub_1 _nettle_sec_sub_1
+#define sec_modinv _nettle_sec_modinv
+#define curve25519_eh_to_x _nettle_curve25519_eh_to_x
+#define curve448_eh_to_x _nettle_curve448_eh_to_x
+
+/* For asserts that are incompatible with sc tests. Currently used
+   only by ECC code. */
+#if WITH_EXTRA_ASSERTS
+# define assert_maybe(x) assert(x)
+#else
+# define assert_maybe(x) ((void)(x))
+#endif
+
+extern const struct ecc_curve _nettle_secp_192r1;
+extern const struct ecc_curve _nettle_secp_224r1;
+extern const struct ecc_curve _nettle_secp_256r1;
+extern const struct ecc_curve _nettle_secp_384r1;
+extern const struct ecc_curve _nettle_secp_521r1;
+
+/* Keep this structure internal for now. It's misnamed (since it's
+   really implementing the equivalent twisted Edwards curve, with
+   different coordinates). And we're not quite ready to provide
+   general ecc operations over an arbitrary type of curve. */
+extern const struct ecc_curve _nettle_curve25519;
+extern const struct ecc_curve _nettle_curve448;
+
+/* GOST curves, visible with underscore prefix for now */
+extern const struct ecc_curve _nettle_gost_gc256b;
+extern const struct ecc_curve _nettle_gost_gc512a;
+
+#define ECC_MAX_SIZE ((521 + GMP_NUMB_BITS - 1) / GMP_NUMB_BITS)
+
+/* Window size for ecc_mul_a. Using 4 bits seems like a good choice,
+   for both Intel x86_64 and ARM Cortex A9. For the larger curves, of
+   384 and 521 bits, we could improve speed by a few percent if we go
+   up to 5 bits, but I don't think that's worth doubling the
+   storage. */
+#define ECC_MUL_A_WBITS 4
+/* And for ecc_mul_a_eh */
+#define ECC_MUL_A_EH_WBITS 4
+
+struct ecc_modulo;
+
+/* Reduces from 2*ecc->size to ecc->size. */
+/* Required to return a result < 2q. This property is inherited by
+   mod_mul and mod_sqr. May clobber input xp. rp may point to the
+   start or the middle of the xp area, but no other overlap is
+   allowed. */
+typedef void ecc_mod_func (const struct ecc_modulo *m, mp_limb_t *rp, mp_limb_t *xp);
+
+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,
+				     const mp_limb_t *up, const mp_limb_t *vp,
+				     mp_limb_t *scratch);
+
+/* Allows in-place operation with r == p, but not r == q */
+typedef void ecc_add_func (const struct ecc_curve *ecc,
+			   mp_limb_t *r,
+			   const mp_limb_t *p, const mp_limb_t *q,
+			   mp_limb_t *scratch);
+
+typedef void ecc_dup_func (const struct ecc_curve *ecc,
+			   mp_limb_t *r, const mp_limb_t *p,
+			   mp_limb_t *scratch);
+
+typedef void ecc_mul_g_func (const struct ecc_curve *ecc, mp_limb_t *r,
+			     const mp_limb_t *np, mp_limb_t *scratch);
+
+typedef void ecc_mul_func (const struct ecc_curve *ecc,
+			   mp_limb_t *r,
+			   const mp_limb_t *np, const mp_limb_t *p,
+			   mp_limb_t *scratch);
+
+typedef void ecc_h_to_a_func (const struct ecc_curve *ecc,
+			      int flags,
+			      mp_limb_t *r, const mp_limb_t *p,
+			      mp_limb_t *scratch);
+
+struct ecc_modulo
+{
+  unsigned short bit_size;
+  unsigned short size;
+  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;
+  /* B^size mod m. Expected to have at least 32 leading zeros
+     (equality for secp_256r1). */
+  const mp_limb_t *B;
+  /* 2^{bit_size} - m. When different from B above, for numbers of
+     interest, usually B has trailing zeros and this is B shifted
+     right. */
+  const mp_limb_t *B_shifted;
+  /* For ecc_mod_sub: B^size - 2m, if that doesn't underflow.
+     Otherwise, same as B */
+  const mp_limb_t *Bm2m;
+
+  /* m +/- 1, for redc, excluding redc_size low limbs. */
+  const mp_limb_t *redc_mpm1;
+  /* (m+1)/2 */
+  const mp_limb_t *mp1h;
+
+  ecc_mod_func *mod;
+  ecc_mod_func *reduce;
+  /* 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;
+};
+
+/* Represents an elliptic curve of the form
+
+     y^2 = x^3 - 3x + b (mod p)
+*/
+struct ecc_curve
+{
+  /* The prime p. */
+  struct ecc_modulo p;
+  /* Group order. Currently, many functions rely on q.size ==
+     p.size. */
+  struct ecc_modulo q;
+
+  unsigned short use_redc;
+  unsigned short pippenger_k;
+  unsigned short pippenger_c;
+
+  unsigned short add_hh_itch;
+  unsigned short add_hhh_itch;
+  unsigned short dup_itch;
+  unsigned short mul_itch;
+  unsigned short mul_g_itch;
+  unsigned short h_to_a_itch;
+
+  ecc_add_func *add_hh;
+  ecc_add_func *add_hhh;
+  ecc_dup_func *dup;
+  ecc_mul_func *mul;
+  ecc_mul_g_func *mul_g;
+  ecc_h_to_a_func *h_to_a;
+
+  /* Curve constant */
+  const mp_limb_t *b;
+
+  /* For redc, same as B mod p, otherwise 1. */
+  const mp_limb_t *unit;
+
+  /* Tables for multiplying by the generator, size determined by k and
+     c. The first 2^c entries are defined by
+
+       T[  j_0 +   j_1 2 +     ... + j_{c-1} 2^{c-1} ]
+         = j_0 g + j_1 2^k g + ... + j_{c-1} 2^{k(c-1)} g
+
+     The following entries differ by powers of 2^{kc},
+
+       T[i] = 2^{kc} T[i-2^c]
+  */
+  const mp_limb_t *pippenger_table;
+};
+
+ecc_mod_func ecc_mod;
+ecc_mod_func ecc_pp1_redc;
+ecc_mod_func ecc_pm1_redc;
+
+ecc_mod_inv_func ecc_mod_inv;
+
+/* Side channel silent. Requires that x < 2m, so checks if x == 0 or x == p */
+int
+ecc_mod_zero_p (const struct ecc_modulo *m, const mp_limb_t *xp);
+
+/* Requires that a < 2m, and ref < m, needs m->size limbs of scratch
+   space. Overlap, a == scratch or ref == scratch, is allowed. */
+int
+ecc_mod_equal_p (const struct ecc_modulo *m, const mp_limb_t *a,
+		 const mp_limb_t *ref, mp_limb_t *scratch);
+
+void
+ecc_mod_add (const struct ecc_modulo *m, mp_limb_t *rp,
+	     const mp_limb_t *ap, const mp_limb_t *bp);
+
+/* If inputs are in the range 0 <= a, b < 2m, then so is the output. */
+void
+ecc_mod_sub (const struct ecc_modulo *m, mp_limb_t *rp,
+	     const mp_limb_t *ap, const mp_limb_t *bp);
+
+void
+ecc_mod_mul_1 (const struct ecc_modulo *m, mp_limb_t *rp,
+	       const mp_limb_t *ap, const mp_limb_t b);
+
+void
+ecc_mod_addmul_1 (const struct ecc_modulo *m, mp_limb_t *rp,
+		  const mp_limb_t *ap, mp_limb_t b);
+void
+ecc_mod_submul_1 (const struct ecc_modulo *m, mp_limb_t *rp,
+		  const mp_limb_t *ap, mp_limb_t b);
+
+/* The mul and sqr function need 2*m->size limbs at tp. rp may overlap
+   ap or bp, and may equal tp or tp + m->size, but no other overlap
+   with tp is allowed. */
+void
+ecc_mod_mul (const struct ecc_modulo *m, mp_limb_t *rp,
+	     const mp_limb_t *ap, const mp_limb_t *bp, mp_limb_t *tp);
+
+void
+ecc_mod_sqr (const struct ecc_modulo *m, mp_limb_t *rp,
+	     const mp_limb_t *ap, mp_limb_t *tp);
+
+/* These mul and sqr functions produce a canonical result, 0 <= R < M.
+   Requirements on input and output areas are similar to the above
+   functions, except that it is *not* allowed to pass rp = tp +
+   m->size.
+ */
+void
+ecc_mod_mul_canonical (const struct ecc_modulo *m, mp_limb_t *rp,
+		       const mp_limb_t *ap, const mp_limb_t *bp, mp_limb_t *tp);
+
+void
+ecc_mod_sqr_canonical (const struct ecc_modulo *m, mp_limb_t *rp,
+		       const mp_limb_t *ap, mp_limb_t *tp);
+
+/* R <-- X^{2^k} mod M. Needs 2*ecc->size limbs of scratch space, same
+   overlap requirements as mul and sqr above. */
+void
+ecc_mod_pow_2k (const struct ecc_modulo *m,
+		mp_limb_t *rp, const mp_limb_t *xp,
+		unsigned k, mp_limb_t *tp);
+
+/* R <-- X^{2^k} Y mod M. Similar requirements as ecc_mod_pow_2k, but
+   rp and yp can't overlap. */
+void
+ecc_mod_pow_2k_mul (const struct ecc_modulo *m,
+		    mp_limb_t *rp, const mp_limb_t *xp,
+		    unsigned k, const mp_limb_t *yp,
+		    mp_limb_t *tp);
+
+/* R <-- X^{2^k + 1}. Here, rp and xp must not overlap. */
+#define ecc_mod_pow_2kp1(m, rp, xp, k, tp) \
+  ecc_mod_pow_2k_mul (m, rp, xp, k, xp, tp)
+
+/* mod q operations. */
+void
+ecc_mod_random (const struct ecc_modulo *m, mp_limb_t *xp,
+		void *ctx, nettle_random_func *random, mp_limb_t *scratch);
+
+/* Converts a point P in affine coordinates into a point R in jacobian
+   coordinates. */
+void
+ecc_a_to_j (const struct ecc_curve *ecc,
+	    mp_limb_t *r, const mp_limb_t *p);
+
+/* Converts a point P in jacobian coordinates into a point R in affine
+   coordinates. If op == 1, produce x coordinate only. If op == 2,
+   produce the x coordinate only, and also reduce it modulo q. */
+void
+ecc_j_to_a (const struct ecc_curve *ecc,
+	    int op,
+	    mp_limb_t *r, const mp_limb_t *p,
+	    mp_limb_t *scratch);
+
+/* Converts a point P in homogeneous coordinates on an Edwards curve
+   to affine coordinates. Meaning of op is the same as for
+   ecc_j_to_a. */
+void
+ecc_eh_to_a (const struct ecc_curve *ecc,
+	     int op,
+	     mp_limb_t *r, const mp_limb_t *p,
+	     mp_limb_t *scratch);
+
+/* Group operations */
+
+/* Point doubling, with jacobian input and output. Corner cases:
+   Correctly sets R = 0 (r_Z = 0) if p = 0 or 2p = 0. */
+void
+ecc_dup_jj (const struct ecc_curve *ecc,
+	    mp_limb_t *r, const mp_limb_t *p,
+	    mp_limb_t *scratch);
+
+/* Point addition, with jacobian output, one jacobian input and one
+   affine input. Corner cases: Fails for the cases
+
+     P = Q != 0                       Duplication of non-zero point
+     P = 0, Q != 0 or P != 0, Q = 0   One input zero
+
+     Correctly gives R = 0 if P = Q = 0 or P = -Q. */
+void
+ecc_add_jja (const struct ecc_curve *ecc,
+	     mp_limb_t *r, const mp_limb_t *p, const mp_limb_t *q,
+	     mp_limb_t *scratch);
+
+/* Point addition with Jacobian input and output. */
+void
+ecc_add_jjj (const struct ecc_curve *ecc,
+	     mp_limb_t *r, const mp_limb_t *p, const mp_limb_t *q,
+	     mp_limb_t *scratch);
+
+/* Variant that handles the checks for the special cases P = ±Q.
+   Returns 1 on success, 0 if result is infinite. Not side-channel
+   silent, so must not be used with secret inputs. */
+int
+ecc_nonsec_add_jjj (const struct ecc_curve *ecc,
+		    mp_limb_t *r, const mp_limb_t *p, const mp_limb_t *q,
+		    mp_limb_t *scratch);
+
+/* Point doubling on a twisted Edwards curve, with homogeneous
+   cooordinates. */
+void
+ecc_dup_eh (const struct ecc_curve *ecc,
+	    mp_limb_t *r, const mp_limb_t *p,
+	    mp_limb_t *scratch);
+
+void
+ecc_add_eh (const struct ecc_curve *ecc,
+	    mp_limb_t *r, const mp_limb_t *p, const mp_limb_t *q,
+	    mp_limb_t *scratch);
+
+void
+ecc_add_ehh (const struct ecc_curve *ecc,
+	     mp_limb_t *r, const mp_limb_t *p, const mp_limb_t *q,
+	     mp_limb_t *scratch);
+
+void
+ecc_dup_th (const struct ecc_curve *ecc,
+	    mp_limb_t *r, const mp_limb_t *p,
+	    mp_limb_t *scratch);
+
+void
+ecc_add_th (const struct ecc_curve *ecc,
+	    mp_limb_t *r, const mp_limb_t *p, const mp_limb_t *q,
+	    mp_limb_t *scratch);
+
+void
+ecc_add_thh (const struct ecc_curve *ecc,
+	     mp_limb_t *r, const mp_limb_t *p, const mp_limb_t *q,
+	     mp_limb_t *scratch);
+
+/* Computes N * the group generator. N is an array of ecc_size()
+   limbs. It must be in the range 0 < N < group order, then R != 0,
+   and the algorithm can work without any intermediate values getting
+   to zero. */
+void
+ecc_mul_g (const struct ecc_curve *ecc, mp_limb_t *r,
+	   const mp_limb_t *np, mp_limb_t *scratch);
+
+/* Computes N * P. The scalar N is the same as for ecc_mul_g. P is a
+   non-zero point on the curve, in affine coordinates. Output R is a
+   non-zero point, in Jacobian coordinates. */
+void
+ecc_mul_a (const struct ecc_curve *ecc,
+	   mp_limb_t *r,
+	   const mp_limb_t *np, const mp_limb_t *p,
+	   mp_limb_t *scratch);
+
+void
+ecc_mul_g_eh (const struct ecc_curve *ecc, mp_limb_t *r,
+	      const mp_limb_t *np, mp_limb_t *scratch);
+
+void
+ecc_mul_a_eh (const struct ecc_curve *ecc,
+	      mp_limb_t *r,
+	      const mp_limb_t *np, const mp_limb_t *p,
+	      mp_limb_t *scratch);
+
+void
+ecc_mul_m (const struct ecc_modulo *m,
+	   mp_limb_t a24,
+	   unsigned bit_low, unsigned bit_high,
+	   mp_limb_t *qx, const uint8_t *n, const mp_limb_t *px,
+	   mp_limb_t *scratch);
+
+/* The cnd argument must be 1 or 0. */
+void
+cnd_copy (int cnd, mp_limb_t *rp, const mp_limb_t *ap, mp_size_t n);
+
+mp_limb_t
+sec_add_1 (mp_limb_t *rp, mp_limb_t *ap, mp_size_t n, mp_limb_t b);
+
+mp_limb_t
+sec_sub_1 (mp_limb_t *rp, mp_limb_t *ap, mp_size_t n, mp_limb_t b);
+
+void
+curve25519_eh_to_x (mp_limb_t *xp, const mp_limb_t *p,
+		    mp_limb_t *scratch);
+
+void
+curve448_eh_to_x (mp_limb_t *xp, const mp_limb_t *p,
+		  mp_limb_t *scratch);
+
+/* Current scratch needs: */
+#define ECC_MOD_INV_ITCH(size) (3*(size))
+#define ECC_J_TO_A_ITCH(size, inv) ((size)+(inv))
+#define ECC_EH_TO_A_ITCH(size, inv) ((size)+(inv))
+#define ECC_DUP_JJ_ITCH(size) (4*(size))
+#define ECC_DUP_EH_ITCH(size) (3*(size))
+#define ECC_DUP_TH_ITCH(size) (3*(size))
+#define ECC_ADD_JJA_ITCH(size) (5*(size))
+#define ECC_ADD_JJJ_ITCH(size) (5*(size))
+#define ECC_ADD_EH_ITCH(size) (4*(size))
+#define ECC_ADD_EHH_ITCH(size) (4*(size))
+#define ECC_ADD_TH_ITCH(size) (4*(size))
+#define ECC_ADD_THH_ITCH(size) (4*(size))
+#define ECC_MUL_G_ITCH(size) (8*(size))
+#define ECC_MUL_G_EH_ITCH(size) (7*(size))
+#if ECC_MUL_A_WBITS == 0
+#define ECC_MUL_A_ITCH(size) (11*(size))
+#else
+#define ECC_MUL_A_ITCH(size) \
+  (((3 << ECC_MUL_A_WBITS) + 8) * (size))
+#endif
+#if ECC_MUL_A_EH_WBITS == 0
+#define ECC_MUL_A_EH_ITCH(size) (10*(size))
+#else
+#define ECC_MUL_A_EH_ITCH(size) \
+  (((3 << ECC_MUL_A_EH_WBITS) + 7) * (size))
+#endif
+#define ECC_MUL_M_ITCH(size) (8*(size))
+#define ECC_ECDSA_SIGN_ITCH(size) (11*(size))
+#define ECC_GOSTDSA_SIGN_ITCH(size) (11*(size))
+#define ECC_MOD_RANDOM_ITCH(size) (size)
+#define ECC_HASH_ITCH(size) (1+(size))
+
+#endif /* NETTLE_ECC_INTERNAL_H_INCLUDED */

ecc-j-to-a.c

+/* ecc-j-to-a.c
+
+   Copyright (C) 2013 Niels Möller
+
+   This file is part of GNU Nettle.
+
+   GNU Nettle is free software: you can redistribute it and/or
+   modify it under the terms of either:
+
+     * the GNU Lesser General Public License as published by the Free
+       Software Foundation; either version 3 of the License, or (at your
+       option) any later version.
+
+   or
+
+     * the GNU General Public License as published by the Free
+       Software Foundation; either version 2 of the License, or (at your
+       option) any later version.
+
+   or both in parallel, as here.
+
+   GNU Nettle is distributed in the hope that it will be useful,
+   but WITHOUT ANY WARRANTY; without even the implied warranty of
+   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+   General Public License for more details.
+
+   You should have received copies of the GNU General Public License and
+   the GNU Lesser General Public License along with this program.  If
+   not, see http://www.gnu.org/licenses/.
+*/
+
+/* Development of Nettle's ECC support was funded by the .SE Internet Fund. */
+
+#if HAVE_CONFIG_H
+# include "config.h"
+#endif
+
+#include "ecc.h"
+#include "ecc-internal.h"
+
+void
+ecc_j_to_a (const struct ecc_curve *ecc,
+	    int op,
+	    mp_limb_t *r, const mp_limb_t *p,
+	    mp_limb_t *scratch)
+{
+#define izp   scratch
+#define iz2p (scratch + ecc->p.size)
+#define iz3p (scratch + 2*ecc->p.size)
+#define tp    scratch
+
+  ecc->p.invert (&ecc->p, izp, p+2*ecc->p.size, izp + ecc->p.size);
+  ecc_mod_sqr (&ecc->p, iz2p, izp, iz2p);
+
+  if (ecc->use_redc)
+    {
+      /* Divide this common factor by B, instead of applying redc to
+	 both x and y outputs. */
+      mpn_zero (iz2p + ecc->p.size, ecc->p.size);
+      ecc->p.reduce (&ecc->p, iz2p, iz2p);
+    }
+
+  /* r_x <-- x / z^2 */
+  ecc_mod_mul_canonical (&ecc->p, r, iz2p, p, iz3p);
+  if (op)
+    {
+      /* Skip y coordinate */
+      if (op > 1)
+	{
+	  mp_limb_t cy;
+	  /* Also reduce the x coordinate mod ecc->q. It should
+	     already be < 2*ecc->q, so one subtraction should
+	     suffice. */
+	  cy = mpn_sub_n (scratch, r, ecc->q.m, ecc->p.size);
+	  cnd_copy (1 - cy, r, scratch, ecc->p.size);
+	}
+      return;
+    }
+  ecc_mod_mul (&ecc->p, iz3p, iz2p, izp, iz3p);
+  ecc_mod_mul_canonical (&ecc->p, r + ecc->p.size, iz3p, p + ecc->p.size, tp);
+
+#undef izp
+#undef iz2p
+#undef iz3p
+#undef tp
+}

ecc-mod-arith.c

+/* ecc-mod-arith.c
+
+   Copyright (C) 2013, 2014 Niels Möller
+
+   This file is part of GNU Nettle.
+
+   GNU Nettle is free software: you can redistribute it and/or
+   modify it under the terms of either:
+
+     * the GNU Lesser General Public License as published by the Free
+       Software Foundation; either version 3 of the License, or (at your
+       option) any later version.
+
+   or
+
+     * the GNU General Public License as published by the Free
+       Software Foundation; either version 2 of the License, or (at your
+       option) any later version.
+
+   or both in parallel, as here.
+
+   GNU Nettle is distributed in the hope that it will be useful,
+   but WITHOUT ANY WARRANTY; without even the implied warranty of
+   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+   General Public License for more details.
+
+   You should have received copies of the GNU General Public License and
+   the GNU Lesser General Public License along with this program.  If
+   not, see http://www.gnu.org/licenses/.
+*/
+
+/* Development of Nettle's ECC support was funded by the .SE Internet Fund. */
+
+#if HAVE_CONFIG_H
+# include "config.h"
+#endif
+
+#include <assert.h>
+
+#include "ecc-internal.h"
+
+/* Routines for modp arithmetic. All values are ecc->size limbs, but
+   not necessarily < p. */
+
+int
+ecc_mod_zero_p (const struct ecc_modulo *m, const mp_limb_t *xp_in)
+{
+  volatile mp_limb_t is_non_zero, is_not_p;
+  const volatile mp_limb_t *xp;
+  mp_size_t i;
+
+  for (xp = xp_in, i = 0, is_non_zero = is_not_p = 0; i < m->size; i++)
+    {
+      is_non_zero |= xp[i];
+      is_not_p |= (xp[i] ^ m->m[i]);
+    }
+
+  return is_zero_limb (is_non_zero) | is_zero_limb (is_not_p);
+}
+
+int
+ecc_mod_equal_p (const struct ecc_modulo *m, const mp_limb_t *a,
+		 const mp_limb_t *ref, mp_limb_t *scratch)
+{
+  mp_limb_t cy;
+  cy = mpn_sub_n (scratch, a, ref, m->size);
+  /* If cy > 0, i.e., a < ref, then they can't be equal mod m. */
+  return (1 - cy) & ecc_mod_zero_p (m, scratch);
+}
+
+void
+ecc_mod_add (const struct ecc_modulo *m, mp_limb_t *rp,
+	     const mp_limb_t *ap, const mp_limb_t *bp)
+{
+  mp_limb_t cy;
+  cy = mpn_add_n (rp, ap, bp, m->size);
+  cy = mpn_cnd_add_n (cy, rp, rp, m->B, m->size);
+  cy = mpn_cnd_add_n (cy, rp, rp, m->B, m->size);
+  assert_maybe (cy == 0);
+}
+
+void
+ecc_mod_sub (const struct ecc_modulo *m, mp_limb_t *rp,
+	     const mp_limb_t *ap, const mp_limb_t *bp)
+{
+  mp_limb_t cy;
+  cy = mpn_sub_n (rp, ap, bp, m->size);
+  /* The adjustments for this function work differently depending on
+     the value of the most significant bit of m.
+
+     If m has a most significant bit of zero, then the first
+     adjustment step conditionally adds 2m. If in addition, inputs are
+     in the 0 <= a,b < 2m range, then the first adjustment guarantees
+     that result is in that same range. The second adjustment step is
+     needed only if b > 2m, it then ensures output is correct modulo
+     m, but nothing more.
+
+     If m has a most significant bit of one, Bm2m and B are the same,
+     and this function works analogously to ecc_mod_add.
+   */
+  cy = mpn_cnd_sub_n (cy, rp, rp, m->Bm2m, m->size);
+  cy = mpn_cnd_sub_n (cy, rp, rp, m->B, m->size);
+  assert_maybe (cy == 0);
+}
+
+void
+ecc_mod_mul_1 (const struct ecc_modulo *m, mp_limb_t *rp,
+	       const mp_limb_t *ap, mp_limb_t b)
+{
+  mp_limb_t hi;
+
+  assert (b <= 0xffffffff);
+  hi = mpn_mul_1 (rp, ap, m->size, b);
+  hi = mpn_addmul_1 (rp, m->B, m->size, hi);
+  assert_maybe (hi <= 1);
+  hi = mpn_cnd_add_n (hi, rp, rp, m->B, m->size);
+  /* Sufficient if b < B^size / p */
+  assert_maybe (hi == 0);
+}
+
+void
+ecc_mod_addmul_1 (const struct ecc_modulo *m, mp_limb_t *rp,
+		  const mp_limb_t *ap, mp_limb_t b)
+{
+  mp_limb_t hi;
+
+  assert (b <= 0xffffffff);
+  hi = mpn_addmul_1 (rp, ap, m->size, b);
+  hi = mpn_addmul_1 (rp, m->B, m->size, hi);
+  assert_maybe (hi <= 1);
+  hi = mpn_cnd_add_n (hi, rp, rp, m->B, m->size);
+  /* Sufficient roughly if b < B^size / p */
+  assert_maybe (hi == 0);
+}
+  
+void
+ecc_mod_submul_1 (const struct ecc_modulo *m, mp_limb_t *rp,
+		  const mp_limb_t *ap, mp_limb_t b)
+{
+  mp_limb_t hi;
+
+  assert (b <= 0xffffffff);
+  hi = mpn_submul_1 (rp, ap, m->size, b);
+  hi = mpn_submul_1 (rp, m->B, m->size, hi);
+  assert_maybe (hi <= 1);
+  hi = mpn_cnd_sub_n (hi, rp, rp, m->B, m->size);
+  /* Sufficient roughly if b < B^size / p */
+  assert_maybe (hi == 0);
+}
+
+void
+ecc_mod_mul (const struct ecc_modulo *m, mp_limb_t *rp,
+	     const mp_limb_t *ap, const mp_limb_t *bp, mp_limb_t *tp)
+{
+  mpn_mul_n (tp, ap, bp, m->size);
+  m->reduce (m, rp, tp);
+}
+
+void
+ecc_mod_sqr (const struct ecc_modulo *m, mp_limb_t *rp,
+	     const mp_limb_t *ap, mp_limb_t *tp)
+{
+  mpn_sqr (tp, ap, m->size);
+  m->reduce (m, rp, tp);
+}
+
+void
+ecc_mod_mul_canonical (const struct ecc_modulo *m, mp_limb_t *rp,
+		       const mp_limb_t *ap, const mp_limb_t *bp, mp_limb_t *tp)
+{
+  mp_limb_t cy;
+  mpn_mul_n (tp, ap, bp, m->size);
+  m->reduce (m, tp + m->size, tp);
+
+  cy = mpn_sub_n (rp, tp + m->size, m->m, m->size);
+  cnd_copy (cy, rp, tp + m->size, m->size);
+}
+
+void
+ecc_mod_sqr_canonical (const struct ecc_modulo *m, mp_limb_t *rp,
+		       const mp_limb_t *ap, mp_limb_t *tp)
+{
+  mp_limb_t cy;
+  mpn_sqr (tp, ap, m->size);
+  m->reduce (m, tp + m->size, tp);
+
+  cy = mpn_sub_n (rp, tp + m->size, m->m, m->size);
+  cnd_copy (cy, rp, tp + m->size, m->size);
+}
+
+void
+ecc_mod_pow_2k (const struct ecc_modulo *m,
+		mp_limb_t *rp, const mp_limb_t *xp,
+		unsigned k, mp_limb_t *tp)
+{
+  ecc_mod_sqr (m, rp, xp, tp);
+  while (--k > 0)
+    ecc_mod_sqr (m, rp, rp, tp);
+}
+
+void
+ecc_mod_pow_2k_mul (const struct ecc_modulo *m,
+		    mp_limb_t *rp, const mp_limb_t *xp,
+		    unsigned k, const mp_limb_t *yp,
+		    mp_limb_t *tp)
+{
+  ecc_mod_pow_2k (m, rp, xp, k, tp);
+  ecc_mod_mul (m, rp, rp, yp, tp);
+}

ecc-mod-inv.c

+/* ecc-mod-inv.c
+
+   Copyright (C) 2013, 2014 Niels Möller
+
+   This file is part of GNU Nettle.
+
+   GNU Nettle is free software: you can redistribute it and/or
+   modify it under the terms of either:
+
+     * the GNU Lesser General Public License as published by the Free
+       Software Foundation; either version 3 of the License, or (at your
+       option) any later version.
+
+   or
+
+     * the GNU General Public License as published by the Free
+       Software Foundation; either version 2 of the License, or (at your
+       option) any later version.
+
+   or both in parallel, as here.
+
+   GNU Nettle is distributed in the hope that it will be useful,
+   but WITHOUT ANY WARRANTY; without even the implied warranty of
+   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+   General Public License for more details.
+
+   You should have received copies of the GNU General Public License and
+   the GNU Lesser General Public License along with this program.  If
+   not, see http://www.gnu.org/licenses/.
+*/
+
+/* Development of Nettle's ECC support was funded by the .SE Internet Fund. */
+
+#if HAVE_CONFIG_H
+# include "config.h"
+#endif
+
+#include <assert.h>
+
+#include "ecc-internal.h"
+
+static void
+cnd_neg (int cnd, mp_limb_t *rp, const mp_limb_t *ap, mp_size_t n)
+{
+  mp_limb_t cy = (cnd != 0);
+  mp_limb_t mask = -cy;
+  mp_size_t i;
+
+  for (i = 0; i < n; i++)
+    {
+      mp_limb_t r = (ap[i] ^ mask) + cy;
+      cy = r < cy;
+      rp[i] = r;
+    }
+}
+
+/* Compute a^{-1} mod m, with running time depending only on the size.
+   Returns zero if a == 0 (mod m), to be consistent with a^{phi(m)-1}.
+   Also needs (m+1)/2, and m must be odd.
+
+   Needs 3n limbs of scratch space.
+*/
+
+/* FIXME: Could use mpn_sec_invert (in GMP-6), but with a bit more
+   scratch need since it doesn't precompute (m+1)/2. */
+void
+ecc_mod_inv (const struct ecc_modulo *m,
+	     mp_limb_t *vp, const mp_limb_t *in_ap,
+	     mp_limb_t *scratch)
+{
+#define ap scratch
+#define bp (scratch + n)
+#define up (scratch + 2*n)
+
+  mp_size_t n = m->size;
+  /* Avoid the mp_bitcnt_t type for compatibility with older GMP
+     versions. */  
+  unsigned i;
+
+  /* Maintain
+
+       a = u * orig_a (mod m)
+       b = v * orig_a (mod m)
+
+     and b odd at all times. Initially,
+
+       a = a_orig, u = 1
+       b = m,      v = 0
+     */
+
+  assert (ap != vp);
+
+  up[0] = 1;
+  mpn_zero (up+1, n - 1);
+  mpn_copyi (bp, m->m, n);
+  mpn_zero (vp, n);
+  mpn_copyi (ap, in_ap, n);
+
+  for (i = m->bit_size + GMP_NUMB_BITS * n; i-- > 0; )
+    {
+      mp_limb_t odd, swap, cy;
+      
+      /* Always maintain b odd. The logic of the iteration is as
+	 follows. For a, b:
+
+	   odd = a & 1
+	   a -= odd * b
+	   if (underflow from a-b)
+	     {
+	       b += a, assigns old a
+	       a = B^n-a
+	     }
+	   
+	   a /= 2
+
+	 For u, v:
+
+	   if (underflow from a - b)
+	     swap u, v
+	   u -= odd * v
+	   if (underflow from u - v)
+	     u += m
+
+	   u /= 2
+	   if (a one bit was shifted out)
+	     u += (m+1)/2
+
+	 As long as a > 0, the quantity
+
+	   (bitsize of a) + (bitsize of b)
+
+	 is reduced by at least one bit per iteration, hence after
+         (bit_size of orig_a) + (bit_size of m) - 1 iterations we
+         surely have a = 0. Then b = gcd(orig_a, m) and if b = 1 then
+         also v = orig_a^{-1} (mod m)
+      */
+
+      assert_maybe (bp[0] & 1);
+      odd = ap[0] & 1;
+
+      swap = mpn_cnd_sub_n (odd, ap, ap, bp, n);
+      mpn_cnd_add_n (swap, bp, bp, ap, n);
+      cnd_neg (swap, ap, ap, n);
+
+      mpn_cnd_swap (swap, up, vp, n);
+      cy = mpn_cnd_sub_n (odd, up, up, vp, n);
+      cy -= mpn_cnd_add_n (cy, up, up, m->m, n);
+      assert_maybe (cy == 0);
+
+      cy = mpn_rshift (ap, ap, n, 1);
+      assert_maybe (cy == 0);
+      cy = mpn_rshift (up, up, n, 1);
+      cy = mpn_cnd_add_n (cy, up, up, m->mp1h, n);
+      assert_maybe (cy == 0);
+    }
+  assert_maybe ( (ap[0] | ap[n-1]) == 0);
+#undef ap
+#undef bp
+#undef up
+}

ecc-mod.c

+/* ecc-mod.c
+
+   Copyright (C) 2013 Niels Möller
+
+   This file is part of GNU Nettle.
+
+   GNU Nettle is free software: you can redistribute it and/or
+   modify it under the terms of either:
+
+     * the GNU Lesser General Public License as published by the Free
+       Software Foundation; either version 3 of the License, or (at your
+       option) any later version.
+
+   or
+
+     * the GNU General Public License as published by the Free
+       Software Foundation; either version 2 of the License, or (at your
+       option) any later version.
+
+   or both in parallel, as here.
+
+   GNU Nettle is distributed in the hope that it will be useful,
+   but WITHOUT ANY WARRANTY; without even the implied warranty of
+   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+   General Public License for more details.
+
+   You should have received copies of the GNU General Public License and
+   the GNU Lesser General Public License along with this program.  If
+   not, see http://www.gnu.org/licenses/.
+*/
+
+/* Development of Nettle's ECC support was funded by the .SE Internet Fund. */
+
+#if HAVE_CONFIG_H
+# include "config.h"
+#endif
+
+#include <assert.h>
+
+#include "ecc-internal.h"
+
+/* Computes r <-- x mod m, input 2*m->size, output m->size. It's
+ * allowed to have rp == xp or rp == xp + m->size, but no other kind
+ * of overlap is allowed. */
+void
+ecc_mod (const struct ecc_modulo *m, mp_limb_t *rp, mp_limb_t *xp)
+{
+  mp_limb_t hi;
+  mp_size_t mn = m->size;
+  mp_size_t bn = m->B_size;
+  mp_size_t sn = mn - bn;
+  mp_size_t rn = 2*mn;
+  mp_size_t i;
+  unsigned shift;
+
+  assert (bn < mn);
+
+  /* FIXME: Could use mpn_addmul_2. */
+  /* Eliminate sn limbs at a time */
+  if (m->B[bn-1] < ((mp_limb_t) 1 << (GMP_NUMB_BITS - 1)))
+    {
+      /* Multiply sn + 1 limbs at a time, so we get a mn+1 limb
+	 product. Then we can absorb the carry in the high limb */
+      while (rn > 2 * mn - bn)
+	{
+	  rn -= sn;
+
+	  for (i = 0; i <= sn; i++)
+	    xp[rn+i-1] = mpn_addmul_1 (xp + rn - mn - 1 + i, m->B, bn, xp[rn+i-1]);
+	  xp[rn-1] = xp[rn+sn-1]
+	    + mpn_add_n (xp + rn - sn - 1, xp + rn - sn - 1, xp + rn - 1, sn);
+	}
+    }
+  else
+    {
+      while (rn > 2 * mn - bn)
+	{
+	  rn -= sn;
+
+	  for (i = 0; i < sn; i++)
+	    xp[rn+i] = mpn_addmul_1 (xp + rn - mn + i, m->B, bn, xp[rn+i]);
+				     
+	  hi = mpn_add_n (xp + rn - sn, xp + rn - sn, xp + rn, sn);
+	  hi = mpn_cnd_add_n (hi, xp + rn - mn, xp + rn - mn, m->B, mn);
+	  assert_maybe (hi == 0);
+	}
+    }
+
+  assert (rn > mn);
+  rn -= mn;
+  assert (rn <= sn);
+
+  for (i = 0; i < rn; i++)
+    xp[mn+i] = mpn_addmul_1 (xp + i, m->B, bn, xp[mn+i]);
+
+  hi = mpn_add_n (xp + bn, xp + bn, xp + mn, rn);
+  if (rn < sn)
+    hi = sec_add_1 (xp + bn + rn, xp + bn + rn, sn - rn, hi);
+
+  shift = m->size * GMP_NUMB_BITS - m->bit_size;
+  if (shift > 0)
+    {
+      /* Combine hi with top bits, add in */
+      hi = (hi << shift) | (xp[mn-1] >> (GMP_NUMB_BITS - shift));
+      xp[mn-1] = (xp[mn-1] & (((mp_limb_t) 1 << (GMP_NUMB_BITS - shift)) - 1))
+	+ mpn_addmul_1 (xp, m->B_shifted, mn-1, hi);
+      /* FIXME: Can this copying be eliminated? */
+      if (rp != xp)
+	mpn_copyi (rp, xp, mn);
+    }
+  else
+    {
+      hi = mpn_cnd_add_n (hi, rp, xp, m->B, mn);
+      assert_maybe (hi == 0);
+    }
+}

ecc-mul-a-eh.c

+/* ecc-mul-a-eh.c
+
+   Copyright (C) 2013, 2014 Niels Möller
+
+   This file is part of GNU Nettle.
+
+   GNU Nettle is free software: you can redistribute it and/or
+   modify it under the terms of either:
+
+     * the GNU Lesser General Public License as published by the Free
+       Software Foundation; either version 3 of the License, or (at your
+       option) any later version.
+
+   or
+
+     * the GNU General Public License as published by the Free
+       Software Foundation; either version 2 of the License, or (at your
+       option) any later version.
+
+   or both in parallel, as here.
+
+   GNU Nettle is distributed in the hope that it will be useful,
+   but WITHOUT ANY WARRANTY; without even the implied warranty of
+   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+   General Public License for more details.
+
+   You should have received copies of the GNU General Public License and
+   the GNU Lesser General Public License along with this program.  If
+   not, see http://www.gnu.org/licenses/.
+*/
+
+#if HAVE_CONFIG_H
+# include "config.h"
+#endif
+
+#include <assert.h>
+
+#include "ecc.h"
+#include "ecc-internal.h"
+
+/* Binary algorithm needs 6*ecc->p.size + scratch for ecc_add_ehh,
+   total 10 ecc->p.size
+
+   Window algorithm needs (3<<w) * ecc->p.size for the table,
+   3*ecc->p.size for a temporary point, and scratch for
+   ecc_add_ehh. */
+
+#if ECC_MUL_A_EH_WBITS == 0
+void
+ecc_mul_a_eh (const struct ecc_curve *ecc,
+	      mp_limb_t *r,
+	      const mp_limb_t *np, const mp_limb_t *p,
+	      mp_limb_t *scratch)
+{
+#define pe scratch
+#define tp (scratch + 3*ecc->p.size)
+#define scratch_out (scratch + 6*ecc->p.size)
+
+  unsigned i;
+
+  ecc_a_to_j (ecc, pe, p);
+
+  /* x = 0, y = 1, z = 1 */
+  mpn_zero (r, 3*ecc->p.size);
+  r[ecc->p.size] = r[2*ecc->p.size] = 1;
+  
+  for (i = ecc->p.size; i-- > 0; )
+    {
+      mp_limb_t w = np[i] << (GMP_LIMB_BITS - GMP_NUMB_BITS);
+      unsigned j;
+
+      for (j = 0; j < GMP_NUMB_BITS; j++, w <<= 1)
+	{
+	  int bit;
+	  ecc->dup (ecc, r, r, scratch_out);
+	  ecc->add_hh (ecc, tp, r, pe, scratch_out);
+
+	  bit = w >> (GMP_LIMB_BITS - 1);
+	  /* If we had a one-bit, use the sum. */
+	  cnd_copy (bit, r, tp, 3*ecc->p.size);
+	}
+    }
+}
+#else /* ECC_MUL_A_EH_WBITS > 1 */
+
+#define TABLE_SIZE (1U << ECC_MUL_A_EH_WBITS)
+#define TABLE_MASK (TABLE_SIZE - 1)
+
+#define TABLE(j) (table + (j) * 3*ecc->p.size)
+
+static void
+table_init (const struct ecc_curve *ecc,
+	    mp_limb_t *table, unsigned bits,
+	    const mp_limb_t *p,
+	    mp_limb_t *scratch)
+{
+  unsigned size = 1 << bits;
+  unsigned j;
+
+  mpn_zero (TABLE(0), 3*ecc->p.size);
+  TABLE(0)[ecc->p.size] = TABLE(0)[2*ecc->p.size] = 1;
+
+  ecc_a_to_j (ecc, TABLE(1), p);
+
+  for (j = 2; j < size; j += 2)
+    {
+      ecc->dup (ecc, TABLE(j), TABLE(j/2), scratch);
+      ecc->add_hh (ecc, TABLE(j+1), TABLE(j), TABLE(1), scratch);
+    }
+}
+
+void
+ecc_mul_a_eh (const struct ecc_curve *ecc,
+	      mp_limb_t *r,
+	      const mp_limb_t *np, const mp_limb_t *p,
+	      mp_limb_t *scratch)
+{
+#define tp scratch
+#define table (scratch + 3*ecc->p.size)
+  mp_limb_t *scratch_out = table + (3*ecc->p.size << ECC_MUL_A_EH_WBITS);
+
+  /* Avoid the mp_bitcnt_t type for compatibility with older GMP
+     versions. */
+  unsigned blocks = (ecc->p.bit_size + ECC_MUL_A_EH_WBITS - 1) / ECC_MUL_A_EH_WBITS;
+  unsigned bit_index = (blocks-1) * ECC_MUL_A_EH_WBITS;
+
+  mp_size_t limb_index = bit_index / GMP_NUMB_BITS;
+  unsigned shift = bit_index % GMP_NUMB_BITS;
+  mp_limb_t w, bits;
+
+  table_init (ecc, table, ECC_MUL_A_EH_WBITS, p, scratch_out);
+
+  w = np[limb_index];
+  bits = w >> shift;
+  if (limb_index < ecc->p.size - 1)
+    bits |= np[limb_index + 1] << (GMP_NUMB_BITS - shift);
+
+  assert (bits < TABLE_SIZE);
+
+  mpn_sec_tabselect (r, table, 3*ecc->p.size, TABLE_SIZE, bits);
+
+  for (;;)
+    {
+      unsigned j;
+      if (shift >= ECC_MUL_A_EH_WBITS)
+	{
+	  shift -= ECC_MUL_A_EH_WBITS;
+	  bits = w >> shift;
+	}
+      else
+	{
+	  if (limb_index == 0)
+	    {
+	      assert (shift == 0);
+	      break;
+	    }
+	  bits = w << (ECC_MUL_A_EH_WBITS - shift);
+	  w = np[--limb_index];
+	  shift = shift + GMP_NUMB_BITS - ECC_MUL_A_EH_WBITS;
+	  bits |= w >> shift;
+	}
+      for (j = 0; j < ECC_MUL_A_EH_WBITS; j++)
+	ecc->dup (ecc, r, r, scratch_out);
+
+      bits &= TABLE_MASK;
+      mpn_sec_tabselect (tp, table, 3*ecc->p.size, TABLE_SIZE, bits);
+      ecc->add_hhh (ecc, r, r, tp, scratch_out);
+    }
+#undef table
+#undef tp
+}
+
+#endif /* ECC_MUL_A_EH_WBITS > 1 */

ecc-mul-a.c

+/* ecc-mul-a.c
+
+   Copyright (C) 2013 Niels Möller
+
+   This file is part of GNU Nettle.
+
+   GNU Nettle is free software: you can redistribute it and/or
+   modify it under the terms of either:
+
+     * the GNU Lesser General Public License as published by the Free
+       Software Foundation; either version 3 of the License, or (at your
+       option) any later version.
+
+   or
+
+     * the GNU General Public License as published by the Free
+       Software Foundation; either version 2 of the License, or (at your
+       option) any later version.
+
+   or both in parallel, as here.
+
+   GNU Nettle is distributed in the hope that it will be useful,
+   but WITHOUT ANY WARRANTY; without even the implied warranty of
+   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+   General Public License for more details.
+
+   You should have received copies of the GNU General Public License and
+   the GNU Lesser General Public License along with this program.  If
+   not, see http://www.gnu.org/licenses/.
+*/
+
+/* Development of Nettle's ECC support was funded by the .SE Internet Fund. */
+
+#if HAVE_CONFIG_H
+# include "config.h"
+#endif
+
+#include <assert.h>
+
+#include "ecc.h"
+#include "ecc-internal.h"
+#include "nettle-internal.h"
+
+/* Binary algorithm needs 6*ecc->p.size + scratch for ecc_add_jja.
+   Current total is 12 ecc->p.size, at most 864 bytes.
+
+   Window algorithm needs (3<<w) * ecc->p.size for the table,
+   3*ecc->p.size for a temporary point, and scratch for
+   ecc_add_jjj. */
+
+#if ECC_MUL_A_WBITS == 0
+void
+ecc_mul_a (const struct ecc_curve *ecc,
+	   mp_limb_t *r,
+	   const mp_limb_t *np, const mp_limb_t *p,
+	   mp_limb_t *scratch)
+{
+#define tp scratch
+#define pj (scratch + 3*ecc->p.size)
+#define scratch_out (scratch + 6*ecc->p.size)
+
+  int is_zero;
+
+  unsigned i;
+
+  ecc_a_to_j (ecc, pj, p);
+  mpn_zero (r, 3*ecc->p.size);
+  
+  for (i = ecc->p.size, is_zero = 1; i-- > 0; )
+    {
+      mp_limb_t w = np[i] << (GMP_LIMB_BITS - GMP_NUMB_BITS);
+      unsigned j;
+
+      for (j = 0; j < GMP_NUMB_BITS; j++, w <<= 1)
+	{
+	  int bit;
+
+	  ecc_dup_jj (ecc, r, r, scratch_out);
+	  ecc_add_jja (ecc, tp, r, pj, scratch_out);
+
+	  bit = w >> (GMP_LIMB_BITS - 1);
+	  /* If is_zero is set, r is the zero point,
+	     and ecc_add_jja produced garbage. */
+	  cnd_copy (is_zero, tp, pj, 3*ecc->p.size);
+	  is_zero &= 1 - bit;
+	  /* If we had a one-bit, use the sum. */
+	  cnd_copy (bit, r, tp, 3*ecc->p.size);
+	}
+    }
+}
+#else /* ECC_MUL_A_WBITS > 1 */
+
+#define TABLE_SIZE (1U << ECC_MUL_A_WBITS)
+#define TABLE_MASK (TABLE_SIZE - 1)
+
+#define TABLE(j) (table + (j) * 3*ecc->p.size)
+
+static void
+table_init (const struct ecc_curve *ecc,
+	    mp_limb_t *table, unsigned bits,
+	    const mp_limb_t *p,
+	    mp_limb_t *scratch)
+{
+  unsigned size = 1 << bits;
+  unsigned j;
+
+  mpn_zero (TABLE(0), 3*ecc->p.size);
+  ecc_a_to_j (ecc, TABLE(1), p);
+
+  for (j = 2; j < size; j += 2)
+    {
+      ecc_dup_jj (ecc, TABLE(j), TABLE(j/2), scratch);
+      ecc_add_jja (ecc, TABLE(j+1), TABLE(j), TABLE(1), scratch);
+    }  
+}
+
+void
+ecc_mul_a (const struct ecc_curve *ecc,
+	   mp_limb_t *r,
+	   const mp_limb_t *np, const mp_limb_t *p,
+	   mp_limb_t *scratch)
+{
+#define tp scratch
+#define table (scratch + 3*ecc->p.size)
+  mp_limb_t *scratch_out = table + (3*ecc->p.size << ECC_MUL_A_WBITS);
+  int is_zero = 0;
+
+  /* Avoid the mp_bitcnt_t type for compatibility with older GMP
+     versions. */
+  unsigned blocks = (ecc->p.bit_size + ECC_MUL_A_WBITS - 1) / ECC_MUL_A_WBITS;
+  unsigned bit_index = (blocks-1) * ECC_MUL_A_WBITS;
+
+  mp_size_t limb_index = bit_index / GMP_NUMB_BITS;
+  unsigned shift = bit_index % GMP_NUMB_BITS;
+  mp_limb_t w, bits;
+
+  table_init (ecc, table, ECC_MUL_A_WBITS, p, scratch_out);
+
+  w = np[limb_index];
+  bits = w >> shift;
+  if (limb_index < ecc->p.size - 1)
+    bits |= np[limb_index + 1] << (GMP_NUMB_BITS - shift);
+
+  assert (bits < TABLE_SIZE);
+
+  mpn_sec_tabselect (r, table, 3*ecc->p.size, TABLE_SIZE, bits);
+  is_zero = IS_ZERO_SMALL (bits);
+
+  for (;;)
+    {
+      int bits_is_zero;
+      unsigned j;
+      if (shift >= ECC_MUL_A_WBITS)
+	{
+	  shift -= ECC_MUL_A_WBITS;
+	  bits = w >> shift;
+	}
+      else
+	{
+	  if (limb_index == 0)
+	    {
+	      assert (shift == 0);
+	      break;
+	    }
+	  bits = w << (ECC_MUL_A_WBITS - shift);
+	  w = np[--limb_index];
+	  shift = shift + GMP_NUMB_BITS - ECC_MUL_A_WBITS;
+	  bits |= w >> shift;
+	}
+      for (j = 0; j < ECC_MUL_A_WBITS; j++)
+	ecc_dup_jj (ecc, r, r, scratch_out);
+
+      bits &= TABLE_MASK;
+      mpn_sec_tabselect (tp, table, 3*ecc->p.size, TABLE_SIZE, bits);
+      cnd_copy (is_zero, r, tp, 3*ecc->p.size);
+      ecc_add_jjj (ecc, tp, tp, r, scratch_out);
+      bits_is_zero = IS_ZERO_SMALL (bits);
+
+      /* Use the sum when valid. ecc_add_jja produced garbage if
+	 is_zero or bits_is_zero. */
+      cnd_copy (1 - (bits_is_zero | is_zero), r, tp, 3*ecc->p.size);
+      is_zero &= bits_is_zero;
+    }
+#undef table
+#undef tp
+}
+
+#endif /* ECC_MUL_A_WBITS > 1 */

ecc-mul-g-eh.c

+/* ecc-mul-g-eh.c
+
+   Copyright (C) 2013, 2014 Niels Möller
+
+   This file is part of GNU Nettle.
+
+   GNU Nettle is free software: you can redistribute it and/or
+   modify it under the terms of either:
+
+     * the GNU Lesser General Public License as published by the Free
+       Software Foundation; either version 3 of the License, or (at your
+       option) any later version.
+
+   or
+
+     * the GNU General Public License as published by the Free
+       Software Foundation; either version 2 of the License, or (at your
+       option) any later version.
+
+   or both in parallel, as here.
+
+   GNU Nettle is distributed in the hope that it will be useful,
+   but WITHOUT ANY WARRANTY; without even the implied warranty of
+   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+   General Public License for more details.
+
+   You should have received copies of the GNU General Public License and
+   the GNU Lesser General Public License along with this program.  If
+   not, see http://www.gnu.org/licenses/.
+*/
+
+/* Development of Nettle's ECC support was funded by the .SE Internet Fund. */
+
+#if HAVE_CONFIG_H
+# include "config.h"
+#endif
+
+#include <assert.h>
+
+#include "ecc.h"
+#include "ecc-internal.h"
+
+void
+ecc_mul_g_eh (const struct ecc_curve *ecc, mp_limb_t *r,
+	      const mp_limb_t *np, mp_limb_t *scratch)
+{
+  /* Scratch need determined by the ecc_add_eh call. Current total is
+     7 * ecc->p.size, at most 392 bytes (for curve448). */
+#define tp scratch
+#define scratch_out (scratch + 3*ecc->p.size)
+
+  unsigned k, c;
+  unsigned i, j;
+  unsigned bit_rows;
+
+  k = ecc->pippenger_k;
+  c = ecc->pippenger_c;
+
+  bit_rows = (ecc->p.bit_size + k - 1) / k;
+
+  /* x = 0, y = 1, z = 1 */
+  mpn_zero (r, 3*ecc->p.size);
+  r[ecc->p.size] = r[2*ecc->p.size] = 1;
+
+  for (i = k; i-- > 0; )
+    {
+      ecc->dup (ecc, r, r, scratch);
+      for (j = 0; j * c < bit_rows; j++)
+	{
+	  unsigned bits;
+	  /* Avoid the mp_bitcnt_t type for compatibility with older GMP
+	     versions. */
+	  unsigned bit_index;
+	  
+	  /* Extract c bits from n, stride k, starting at i + kcj,
+	     ending at i + k (cj + c - 1)*/
+	  for (bits = 0, bit_index = i + k*(c*j+c); bit_index > i + k*c*j; )
+	    {
+	      mp_size_t limb_index;
+	      unsigned shift;
+	      
+	      bit_index -= k;
+
+	      limb_index = bit_index / GMP_NUMB_BITS;
+	      if (limb_index >= ecc->p.size)
+		continue;
+
+	      shift = bit_index % GMP_NUMB_BITS;
+	      bits = (bits << 1) | ((np[limb_index] >> shift) & 1);
+	    }
+	  mpn_sec_tabselect (tp,
+			     (ecc->pippenger_table
+			      + (2*ecc->p.size * (mp_size_t) j << c)),
+			      2*ecc->p.size, 1<<c, bits);
+
+	  ecc->add_hh (ecc, r, r, tp, scratch_out);
+	}
+    }
+#undef tp
+#undef scratch_out
+}

ecc-mul-g.c

+/* ecc-mul-g.c
+
+   Copyright (C) 2013 Niels Möller
+
+   This file is part of GNU Nettle.
+
+   GNU Nettle is free software: you can redistribute it and/or
+   modify it under the terms of either:
+
+     * the GNU Lesser General Public License as published by the Free
+       Software Foundation; either version 3 of the License, or (at your
+       option) any later version.
+
+   or
+
+     * the GNU General Public License as published by the Free
+       Software Foundation; either version 2 of the License, or (at your
+       option) any later version.
+
+   or both in parallel, as here.
+
+   GNU Nettle is distributed in the hope that it will be useful,
+   but WITHOUT ANY WARRANTY; without even the implied warranty of
+   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+   General Public License for more details.
+
+   You should have received copies of the GNU General Public License and
+   the GNU Lesser General Public License along with this program.  If
+   not, see http://www.gnu.org/licenses/.
+*/
+
+/* Development of Nettle's ECC support was funded by the .SE Internet Fund. */
+
+#if HAVE_CONFIG_H
+# include "config.h"
+#endif
+
+#include <assert.h>
+
+#include "ecc.h"
+#include "ecc-internal.h"
+#include "nettle-internal.h"
+
+void
+ecc_mul_g (const struct ecc_curve *ecc, mp_limb_t *r,
+	   const mp_limb_t *np, mp_limb_t *scratch)
+{
+  /* Scratch need determined by the ecc_add_jja call. Current total is
+     8 * ecc->p.size, at most 576 bytes. */
+#define tp scratch
+#define scratch_out (scratch + 3*ecc->p.size)
+
+  unsigned k, c;
+  unsigned i, j;
+  unsigned bit_rows;
+
+  int is_zero;
+
+  k = ecc->pippenger_k;
+  c = ecc->pippenger_c;
+
+  bit_rows = (ecc->p.bit_size + k - 1) / k;
+  
+  mpn_zero (r, 3*ecc->p.size);
+  
+  for (i = k, is_zero = 1; i-- > 0; )
+    {
+      ecc_dup_jj (ecc, r, r, scratch);
+      for (j = 0; j * c < bit_rows; j++)
+	{
+	  unsigned bits;
+	  /* Avoid the mp_bitcnt_t type for compatibility with older GMP
+	     versions. */
+	  unsigned bit_index;
+	  int bits_is_zero;
+
+	  /* Extract c bits from n, stride k, starting at i + kcj,
+	     ending at i + k (cj + c - 1)*/
+	  for (bits = 0, bit_index = i + k*(c*j+c); bit_index > i + k*c*j; )
+	    {
+	      mp_size_t limb_index;
+	      unsigned shift;
+	      
+	      bit_index -= k;
+
+	      limb_index = bit_index / GMP_NUMB_BITS;
+	      if (limb_index >= ecc->p.size)
+		continue;
+
+	      shift = bit_index % GMP_NUMB_BITS;
+	      bits = (bits << 1) | ((np[limb_index] >> shift) & 1);
+	    }
+	  mpn_sec_tabselect (tp,
+			     (ecc->pippenger_table
+			      + (2*ecc->p.size * (mp_size_t) j << c)),
+			      2*ecc->p.size, 1<<c, bits);
+	  cnd_copy (is_zero, r, tp, 2*ecc->p.size);
+	  cnd_copy (is_zero, r + 2*ecc->p.size, ecc->unit, ecc->p.size);
+	  
+	  ecc_add_jja (ecc, tp, r, tp, scratch_out);
+	  bits_is_zero = IS_ZERO_SMALL (bits);
+
+	  /* Use the sum when valid. ecc_add_jja produced garbage if
+	     is_zero or bits_is_zero. */
+	  cnd_copy (1 - (bits_is_zero | is_zero), r, tp, 3*ecc->p.size);
+	  is_zero &= bits_is_zero;
+	}
+    }
+#undef tp
+#undef scratch_out
+}

ecc-mul-m.c

+/* ecc-mul-m.c
+
+   Point multiplication using Montgomery curve representation.
+
+   Copyright (C) 2014 Niels Möller
+
+   This file is part of GNU Nettle.
+
+   GNU Nettle is free software: you can redistribute it and/or
+   modify it under the terms of either:
+
+     * the GNU Lesser General Public License as published by the Free
+       Software Foundation; either version 3 of the License, or (at your
+       option) any later version.
+
+   or
+
+     * the GNU General Public License as published by the Free
+       Software Foundation; either version 2 of the License, or (at your
+       option) any later version.
+
+   or both in parallel, as here.
+
+   GNU Nettle is distributed in the hope that it will be useful,
+   but WITHOUT ANY WARRANTY; without even the implied warranty of
+   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+   General Public License for more details.
+
+   You should have received copies of the GNU General Public License and
+   the GNU Lesser General Public License along with this program.  If
+   not, see http://www.gnu.org/licenses/.
+*/
+
+#if HAVE_CONFIG_H
+# include "config.h"
+#endif
+
+#include <assert.h>
+
+#include "ecc.h"
+#include "ecc-internal.h"
+
+void
+ecc_mul_m (const struct ecc_modulo *m,
+	   mp_limb_t a24,
+	   unsigned bit_low, unsigned bit_high,
+	   mp_limb_t *qx, const uint8_t *n, const mp_limb_t *px,
+	   mp_limb_t *scratch)
+{
+  unsigned i;
+  mp_limb_t swap;
+
+#define x2 (scratch)
+#define z2 (scratch + m->size)
+#define x3 (scratch + 2*m->size)
+#define z3 (scratch + 3*m->size)
+
+  /* Formulas from RFC 7748:
+
+       A = x_2 + z_2
+       AA = A^2
+       B = x_2 - z_2
+       BB = B^2
+       E = AA - BB
+       C = x_3 + z_3
+       D = x_3 - z_3
+       DA = D * A
+       CB = C * B
+       x_3 = (DA + CB)^2
+       z_3 = x_1 * (DA - CB)^2
+       x_2 = AA * BB
+       z_2 = E * (AA + a24 * E)
+
+     For pure doubling, we use:
+
+       A = x_2 + z_2
+       AA = A^2
+       B = x_2 - z_2
+       BB = B^2
+       E = AA - BB
+       x3 = AA * BB
+       z3 =  E * (AA + a24 * E)
+  */
+
+#define A (scratch + 4*m->size)
+#define AA A
+#define D (scratch + 5*m->size)
+#define DA D
+
+#define tp (scratch + 6*m->size)
+
+  /* For the doubling formulas. */
+#define B D
+#define BB D
+#define E D
+
+  /* Initialize, x2 = px, z2 = 1 */
+  mpn_copyi (x2, px, m->size);
+  z2[0] = 1;
+  mpn_zero (z2+1, m->size - 1);
+
+  /* Get x3, z3 from doubling. Since most significant bit is forced to 1. */
+  ecc_mod_add (m, A, x2, z2);
+  ecc_mod_sub (m, B, x2, z2);
+  ecc_mod_sqr (m, AA, A, tp);
+  ecc_mod_sqr (m, BB, B, tp);
+  ecc_mod_mul (m, x3, AA, BB, tp);
+  ecc_mod_sub (m, E, AA, BB);
+  ecc_mod_addmul_1 (m, AA, E, a24);
+  ecc_mod_mul (m, z3, E, AA, tp);
+
+  for (i = bit_high, swap = 0; i >= bit_low; i--)
+    {
+      mp_limb_t bit = (n[i/8] >> (i & 7)) & 1;
+
+      mpn_cnd_swap (swap ^ bit, x2, x3, 2*m->size);
+      swap = bit;
+
+      ecc_mod_add (m, A, x2, z2);
+      ecc_mod_sub (m, D, x3, z3);
+      ecc_mod_mul (m, DA, D, A, tp);
+      ecc_mod_sqr (m, AA, A, tp);
+
+      /* Store B, BB and E at z2 */
+      ecc_mod_sub (m, z2, x2, z2);	/* B */
+      /* Store C and CB at z3 */
+      ecc_mod_add (m, z3, x3, z3);	/* C */
+      ecc_mod_mul (m, z3, z3, z2, tp);	/* CB */
+      ecc_mod_sqr (m, z2, z2, tp);	/* BB */
+
+      /* Finish x2 */
+      ecc_mod_mul (m, x2, AA, z2, tp);
+
+      ecc_mod_sub (m, z2, AA, z2);	/* E */
+
+      /* Finish z2 */
+      ecc_mod_addmul_1 (m, AA, z2, a24);
+      ecc_mod_mul (m, z2, z2, AA, tp);
+
+      /* Finish x3 */
+      ecc_mod_add (m, x3, DA, z3);
+      ecc_mod_sqr (m, x3, x3, tp);
+
+      /* Finish z3 */
+      ecc_mod_sub (m, z3, DA, z3);	/* DA - CB */
+      ecc_mod_sqr (m, z3, z3, tp);
+      ecc_mod_mul (m, z3, z3, px, tp);
+    }
+  mpn_cnd_swap (swap, x2, x3, 2*m->size);
+
+  /* Do the low zero bits, just duplicating x2 */
+  for (i = 0; i < bit_low; i++)
+    {
+      ecc_mod_add (m, A, x2, z2);
+      ecc_mod_sub (m, B, x2, z2);
+      ecc_mod_sqr (m, AA, A, tp);
+      ecc_mod_sqr (m, BB, B, tp);
+      ecc_mod_mul (m, x2, AA, BB, tp);
+      ecc_mod_sub (m, E, AA, BB);
+      ecc_mod_addmul_1 (m, AA, E, a24);
+      ecc_mod_mul (m, z2, E, AA, tp);
+    }
+  assert (m->invert_itch <= 7 * m->size);
+  m->invert (m, x3, z2, z3 + m->size);
+  ecc_mod_mul_canonical (m, qx, x2, x3, z3);
+}