+2022-10-07 Niels Möller <nisse@lysator.liu.se>
+
+ * ecc-mod-inv.c: Implement Bernstein-Yang algorithm. It
+ outperforms older inversion code, except except secp192r1 mod p
+ inversion.
+ * eccdata.c (output_modulo): Output iteration count for new
+ modular inversion algorithm, and needed power of two for redc
+ inversion.
+ * ecc-curve25519.c (ecc_curve25519_inv): Deleted.
+ * ecc-curve448.c (ecc_curve448_inv): Deleted.
+ * ecc-secp224r1.c (ecc_secp224r1_inv): Deleted.
+ * ecc-secp256r1.c (ecc_secp256r1_inv): Deleted.
+ * ecc-secp384r1.c (ecc_secp384r1_inv): Deleted.
+ * ecc-secp521r1.c (ecc_secp521r1_inv): Deleted.
+
+2022-10-03 Niels Möller <nisse@lysator.liu.se>
+
+ * ecc-internal.h (struct ecc_modulo): Delete (m+1)/2 constant, and
+ update all curve structs.
+ * eccdata.c (output_curve): Delete corresponding code to generate
+ ecc_pp1h and ecc_qp1h.
+
2022-11-09 Niels Möller <nisse@lysator.liu.se>
From Mamone Tarsha:
x = 0, and we should be fine, since ecc_mod_inv for ecc->p returns 0
in this case. */
ecc_mod_sub (&ecc->p, t0, wp, vp);
- /* Needs a total of 6*size storage. */
+ /* Needs a total of 2*size + scratch for inversion. */
ecc->p.invert (&ecc->p, t1, t0, tp);
ecc_mod_add (&ecc->p, t0, wp, vp);
# include "config.h"
#endif
+#include <assert.h>
#include <string.h>
#include "curve25519.h"
t[0] &= ~7;
t[CURVE25519_SIZE-1] = (t[CURVE25519_SIZE-1] & 0x3f) | 0x40;
- itch = 4*ecc->p.size + ecc->mul_g_itch;
+ assert (ecc->mul_g_itch <= 2*ecc->p.size + ecc->p.invert_itch);
+ itch = 6*ecc->p.size + ecc->p.invert_itch;
+
scratch = gmp_alloc_limbs (itch);
mpn_set_base256_le (x, ecc->p.size, t, CURVE25519_SIZE);
mp_size_t itch;
mp_limb_t *x;
- itch = m->size + ECC_MUL_M_ITCH(m->size);
+ itch = m->size + ECC_MUL_M_ITCH(m->size, m->invert_itch);
x = gmp_alloc_limbs (itch);
/* Note that 255 % GMP_NUMB_BITS == 0 isn't supported, so x always
x = v^2 / u^2 = (V/W)^2 / (U/W)^2 = (V/U)^2
*/
- /* Needs a total of 5*size storage. */
+ /* Needs a total of size + scratch for inversion. */
ecc->p.invert (&ecc->p, t0, up, tp);
ecc_mod_mul (&ecc->p, t0, t0, vp, tp);
ecc_mod_sqr_canonical (&ecc->p, xp, t0, tp);
t[0] &= ~3;
t[CURVE448_SIZE-1] = (t[CURVE448_SIZE-1] & 0x7f) | 0x80;
- itch = 5*ecc->p.size + ecc->mul_g_itch;
+ itch = 4*ecc->p.size + ecc->mul_g_itch;
scratch = gmp_alloc_limbs (itch);
mpn_set_base256_le (x, ecc->p.size, t, CURVE448_SIZE);
mp_size_t itch;
mp_limb_t *x;
- itch = m->size + ECC_MUL_M_ITCH(m->size);
+ itch = m->size + ECC_MUL_M_ITCH(m->size, m->invert_itch);
x = gmp_alloc_limbs (itch);
mpn_set_base256_le (x, m->size, p, CURVE448_SIZE);
#undef tp
}
-/* Scratch as for ecc_mod_pow_252m3 above. */
-#define ECC_25519_INV_ITCH (4*ECC_LIMB_SIZE)
-
-static void
-ecc_curve25519_inv (const struct ecc_modulo *p,
- mp_limb_t *rp, const mp_limb_t *ap,
- mp_limb_t *scratch)
-{
- /* Addition chain
-
- p - 2 = 2^{255} - 21
- = 1 + 2 (1 + 4 (2^{252}-3))
- */
- ecc_mod_pow_252m3 (p, rp, ap, scratch);
- ecc_mod_sqr (p, rp, rp, scratch);
- ecc_mod_sqr (p, rp, rp, scratch);
- ecc_mod_mul (p, rp, ap, rp, scratch);
- ecc_mod_sqr (p, rp, rp, scratch);
- ecc_mod_mul (p, rp, ap, rp, scratch);
-}
-
static int
ecc_curve25519_zero_p (const struct ecc_modulo *p, mp_limb_t *xp)
{
ECC_LIMB_SIZE,
ECC_BMODP_SIZE,
0,
- ECC_25519_INV_ITCH,
+ ECC_MOD_INV_ITCH(ECC_LIMB_SIZE),
0,
ECC_25519_SQRT_RATIO_ITCH,
+ ECC_INVP_COUNT,
+ ECC_BINVP,
ecc_p,
ecc_Bmodp,
ecc_Bmodp_shifted,
ecc_Bm2p,
NULL,
- ecc_pp1h,
+ NULL,
ecc_curve25519_modp,
ecc_curve25519_modp,
- ecc_curve25519_inv,
+ ecc_mod_inv,
NULL,
ecc_curve25519_sqrt_ratio,
},
ECC_MOD_INV_ITCH (ECC_LIMB_SIZE),
0,
0,
+ ECC_INVQ_COUNT,
+ ECC_BINVQ,
ecc_q,
ecc_Bmodq,
ecc_mBmodq_shifted, /* Use q - 2^{252} instead. */
ecc_Bm2q,
NULL,
- ecc_qp1h,
+ NULL,
ecc_curve25519_modq,
ecc_curve25519_modq,
ECC_DUP_TH_ITCH (ECC_LIMB_SIZE),
ECC_MUL_A_EH_ITCH (ECC_LIMB_SIZE),
ECC_MUL_G_EH_ITCH (ECC_LIMB_SIZE),
- ECC_EH_TO_A_ITCH (ECC_LIMB_SIZE, ECC_25519_INV_ITCH),
+ ECC_EH_TO_A_ITCH (ECC_LIMB_SIZE, ECC_MOD_INV_ITCH(ECC_LIMB_SIZE)),
ecc_add_th,
ecc_add_thh,
#undef tp
}
-#define ECC_CURVE448_INV_ITCH (4*ECC_LIMB_SIZE)
-
-static void ecc_curve448_inv (const struct ecc_modulo *p,
- mp_limb_t *rp, const mp_limb_t *ap,
- mp_limb_t *tp)
-{
- ecc_mod_pow_446m224m1 (p, rp, ap, tp);/* a^{2^446-2^222-1} */
- ecc_mod_sqr (p, rp, rp, tp); /* a^{2^447-2^223-2} */
- ecc_mod_sqr (p, rp, rp, tp); /* a^{2^448-2^224-4} */
- ecc_mod_mul (p, rp, ap, rp, tp); /* a^{2^448-2^224-3} */
-}
-
/* To guarantee that inputs to ecc_mod_zero_p are in the required range. */
#if ECC_LIMB_SIZE * GMP_NUMB_BITS != 448
#error Unsupported limb size
ECC_LIMB_SIZE,
ECC_BMODP_SIZE,
0,
- ECC_CURVE448_INV_ITCH,
+ ECC_MOD_INV_ITCH(ECC_LIMB_SIZE),
0,
ECC_CURVE448_SQRT_RATIO_ITCH,
+ ECC_INVP_COUNT,
+ ECC_BINVP,
ecc_p,
ecc_Bmodp,
ecc_Bmodp_shifted,
ecc_Bm2p,
NULL,
- ecc_pp1h,
+ NULL,
ecc_curve448_modp,
ecc_curve448_modp,
- ecc_curve448_inv,
+ ecc_mod_inv,
NULL,
ecc_curve448_sqrt_ratio,
},
ECC_MOD_INV_ITCH (ECC_LIMB_SIZE),
0,
0,
+ ECC_INVQ_COUNT,
+ ECC_BINVQ,
ecc_q,
ecc_Bmodq,
ecc_Bmodq_shifted,
ecc_Bm2q,
NULL,
- ecc_qp1h,
+ NULL,
ecc_mod,
ecc_mod,
ECC_DUP_EH_ITCH (ECC_LIMB_SIZE),
ECC_MUL_A_EH_ITCH (ECC_LIMB_SIZE),
ECC_MUL_G_EH_ITCH (ECC_LIMB_SIZE),
- ECC_EH_TO_A_ITCH (ECC_LIMB_SIZE, ECC_CURVE448_INV_ITCH),
+ ECC_EH_TO_A_ITCH (ECC_LIMB_SIZE, ECC_MOD_INV_ITCH (ECC_LIMB_SIZE)),
ecc_add_eh,
ecc_add_ehh,
ECC_MOD_INV_ITCH (ECC_LIMB_SIZE),
0,
0,
+ ECC_INVP_COUNT,
+ ECC_BINVP,
ecc_p,
ecc_Bmodp,
ecc_Bmodp_shifted,
ecc_Bm2p,
ecc_redc_ppm1,
+ NULL,
- ecc_pp1h,
ecc_gost_gc256b_modp,
ecc_gost_gc256b_modp,
ecc_mod_inv,
ECC_MOD_INV_ITCH (ECC_LIMB_SIZE),
0,
0,
+ ECC_INVQ_COUNT,
+ ECC_BINVQ,
ecc_q,
ecc_Bmodq,
ecc_Bmodq_shifted,
ecc_Bm2q,
NULL,
- ecc_qp1h,
+ NULL,
ecc_gost_gc256b_modq,
ecc_gost_gc256b_modq,
ECC_MOD_INV_ITCH (ECC_LIMB_SIZE),
0,
0,
+ ECC_INVP_COUNT,
+ ECC_BINVP,
ecc_p,
ecc_Bmodp,
ecc_Bmodp_shifted,
ecc_Bm2p,
ecc_redc_ppm1,
+ NULL,
- ecc_pp1h,
ecc_gost_gc512a_modp,
ecc_gost_gc512a_modp,
ecc_mod_inv,
ECC_MOD_INV_ITCH (ECC_LIMB_SIZE),
0,
0,
+ ECC_INVQ_COUNT,
+ ECC_BINVQ,
ecc_q,
ecc_Bmodq,
ecc_Bmodq_shifted,
ecc_Bm2q,
NULL,
- ecc_qp1h,
+ NULL,
ecc_gost_gc512a_modq,
ecc_gost_gc512a_modq,
unsigned short sqrt_itch;
unsigned short sqrt_ratio_itch;
+ /* Needed iterations for the Bernstein-Yang inversion algorithm, for
+ input size bit_size, or bit_size + 1 in case the most significant
+ bit of the most significant modulo limb is zero. */
+ unsigned short invert_count;
+ /* m^{-1} mod 2^{GMP_NUMB_BITS} */
+ mp_limb_t binv;
+
const mp_limb_t *m;
/* B^size mod m. Expected to have at least 32 leading zeros
(equality for secp_256r1). */
/* m +/- 1, for redc, excluding redc_size low limbs. */
const mp_limb_t *redc_mpm1;
- /* (m+1)/2 */
- const mp_limb_t *mp1h;
+
+ const mp_limb_t *invert_power;
ecc_mod_func *mod;
ecc_mod_func *reduce;
mp_limb_t *scratch);
/* Current scratch needs: */
-#define ECC_MOD_INV_ITCH(size) (3*(size))
+#define ECC_MOD_INV_ITCH(size) (5*(size)+4)
#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_MUL_A_EH_ITCH(size) \
(((3 << ECC_MUL_A_EH_WBITS) + 7) * (size))
#endif
-#define ECC_MUL_M_ITCH(size) (8*(size))
+#define ECC_MUL_M_ITCH(size, inv) (4*(size) + (inv))
#define ECC_ECDSA_SIGN_ITCH(size) (11*(size))
#define ECC_GOSTDSA_SIGN_ITCH(size) (11*(size))
#define ECC_MOD_RANDOM_ITCH(size) (size)
/* ecc-mod-inv.c
- Copyright (C) 2013, 2014 Niels Möller
+ Copyright (C) 2013, 2014, 2022 Niels Möller
This file is part of GNU Nettle.
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-internal.h"
-static void
+#if GMP_LIMB_BITS != GMP_NUMB_BITS
+#error Unsupported configuration
+#endif
+
+static mp_limb_t
cnd_neg (int cnd, mp_limb_t *rp, const mp_limb_t *ap, mp_size_t n)
{
mp_limb_t cy = (cnd != 0);
cy = r < cy;
rp[i] = r;
}
+ /* Return cnd except R = A = 0. */
+ return cnd & (1-cy);
}
-/* 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.
+#define BITCNT_BITS (sizeof(mp_bitcnt_t) * CHAR_BIT)
+
+#define STEPS (GMP_NUMB_BITS - 2)
+
+struct matrix {
+ /* Matrix elements interpreted as signed two's complement. Absolute
+ value of elements is at most 2^STEPS. */
+ mp_limb_t a[2][2];
+};
+
+/* Conditionally set (a, b) <-- (b, -a) */
+#define CND_NEG_SWAP_LIMB(cnd, a, b) do {\
+ mp_limb_t __cnd_sum = (a) + (b); \
+ mp_limb_t __cnd_diff = (a) - (b); \
+ (a) -= __cnd_diff & -(cnd); \
+ (b) -= __cnd_sum & -(cnd); \
+ } while (0)
+
+/* Perform K elementary reduction steps on (delta, f, g). Only the
+ least significant K bits of f, g matter. Note that delta has an
+ unsigned type, but is used as two's complement. */
+static mp_bitcnt_t
+steps(struct matrix *M, unsigned k, mp_bitcnt_t delta, mp_limb_t f, mp_limb_t g)
+{
+ mp_limb_t a00, a01, a10, a11;
+ assert (f & 1);
+
+ /* Identity matrix. */
+ a00 = a11 = 1;
+ a01 = a10 = 0;
+
+ /* Preserve invariant (f ; g) = 2^{-i} M (f_orig, g_orig) */
+ for (; k-- > 0; delta++)
+ {
+ mp_limb_t odd = g & 1;
+ mp_limb_t swap = odd & ~(delta >> (BITCNT_BITS-1));
+
+ /* Swap; (f'; g') = (g; -f) = (0,1;-1,0) (g;f) */
+ CND_NEG_SWAP_LIMB(swap, f, g);
+ CND_NEG_SWAP_LIMB(swap, a00, a10);
+ CND_NEG_SWAP_LIMB(swap, a01, a11);
+
+ /* Conditional negation. */
+ delta = (delta ^ - (mp_bitcnt_t) swap) + swap;
+
+ /* Cancel low bit and shift. */
+ g += f & -odd;
+ a10 += a00 & -odd;
+ a11 += a01 & -odd;
+
+ g >>= 1;
+ a00 <<= 1;
+ a01 <<= 1;
+ }
+ M->a[0][0] = a00; M->a[0][1] = a01;
+ M->a[1][0] = a10; M->a[1][1] = a11;
+ return delta;
+}
+
+/* Set R = (u * F + v * G), treating all numbers as two's complement.
+ No overlap allowed. */
+static mp_limb_t
+add_add_mul (mp_ptr rp, const mp_limb_t *fp, const mp_limb_t *gp, mp_size_t n,
+ mp_limb_t u, mp_limb_t v) {
+ mp_limb_t f_sign = fp[n-1] >> (GMP_LIMB_BITS - 1);
+ mp_limb_t g_sign = gp[n-1] >> (GMP_LIMB_BITS - 1);
+ mp_limb_t u_sign = u >> (GMP_LIMB_BITS - 1);
+ mp_limb_t v_sign = v >> (GMP_LIMB_BITS - 1);
+ mp_limb_t hi = mpn_mul_1 (rp, fp, n, u) - ((-f_sign) & u);
+ hi -= ((-u_sign) & fp[n-1]) + mpn_cnd_sub_n (u_sign, rp + 1, rp + 1, fp, n-1);
+
+ hi += mpn_addmul_1 (rp, gp, n, v) - ((-g_sign) & v);
+ hi -= ((-v_sign) & gp[n-1]) + mpn_cnd_sub_n (v_sign, rp + 1, rp + 1, gp, n-1);
+
+ return hi;
+}
+
+/* Update (f'; g') = M (f; g) / 2^{shift}, where all numbers are two's complement. */
+static void
+matrix_vector_mul (const struct matrix *M, unsigned shift,
+ mp_size_t n, mp_limb_t *fp, mp_limb_t *gp, mp_limb_t *tp)
+{
+ mp_limb_t f_hi = add_add_mul (tp, fp, gp, n, M->a[0][0], M->a[0][1]);
+ mp_limb_t g_hi = add_add_mul (tp + n, fp, gp, n, M->a[1][0], M->a[1][1]);
+ mp_limb_t lo = mpn_rshift (fp, tp, n, shift);
+ assert (lo == 0);
+ fp[n-1] += (f_hi << (GMP_LIMB_BITS - shift));
+ lo = mpn_rshift (gp, tp + n, n, shift);
+ assert (lo == 0);
+ gp[n-1] += (g_hi << (GMP_LIMB_BITS - shift));
+}
+
+/* Set R = (u * F + v * G) (mod M), treating u, v as two's complement,
+ but F, G, R unsigned. No overlap allowed. n limb inputs, n+1 limb
+ output.
+
+ Input can be allowed in the range 0 <= F, G < min (2 M, B^n),
+ output always in the range 0 <= R < B M.
+*/
+static void
+add_add_mul_mod (mp_ptr rp, const mp_limb_t *fp, const mp_limb_t *gp,
+ const mp_limb_t *mp, mp_size_t n,
+ mp_limb_t u, mp_limb_t v) {
+ mp_limb_t u_sign = u >> (GMP_LIMB_BITS - 1);
+ mp_limb_t v_sign = v >> (GMP_LIMB_BITS - 1);
+ mp_limb_t r_sign;
+
+ assert ((fp[n-1] >> 1) <= mp[n-1]);
+ assert ((gp[n-1] >> 1) <= mp[n-1]);
+
+ rp[n] = mpn_mul_1 (rp, fp, n, u);
+ mpn_cnd_sub_n (u_sign, rp + 1, rp + 1, fp, n);
+
+ rp[n] += mpn_addmul_1 (rp, gp, n, v);
+ mpn_cnd_sub_n (v_sign, rp + 1, rp + 1, gp, n);
+
+ /* Row sums of the matrix have absolute value <= B/4. With inputs F,
+ G < 2 M, at this point we have
+
+ |R| < B M/2.
+
+ If R < 0, then R + B M is positive negative, adding B M makes the
+ result positive.
+ */
+ r_sign = rp[n] >> (GMP_LIMB_BITS - 1);
+ r_sign -= mpn_cnd_add_n (r_sign, rp + 1, rp + 1, mp, n);
+ assert (r_sign == 0);
+ assert (rp[n] <= mp[n-1]);
+}
+
+/* Input in range 0 <= T < B M, output in range 0 <= R < M. */
+static void
+redc_1 (mp_limb_t *rp, mp_limb_t *tp,
+ const mp_limb_t *mp, mp_size_t n, mp_limb_t m_binv)
+{
+ mp_limb_t cy, hi;
+ /* If we knew that 2M < B^n, we could allow result value in the
+ range 0 <= R < 2M, and use mpn_add_mul_1 without any adjustment
+ step. */
+ hi = mpn_submul_1 (tp, mp, n, m_binv * tp[0]);
+ assert (tp[0] == 0);
+ cy = tp[n] < hi;
+ tp[n] -= hi;
+
+ cy -= mpn_cnd_add_n (cy, rp, tp + 1, mp, n);
+ assert (cy == 0);
+}
+
+/* Update (u'; v') = M (u; v) B^{-1} (mod m), where u, v, m are
+ unsigned n limbs, M has signed elements, and B is the bignum
+ base.
- Needs 3n limbs of scratch space.
+ Output are canonically reduced, 0 <= U, V < M, but
+ inputs are only required to be < 2 M.
*/
+static void
+matrix_vector_mul_mod (const struct matrix *M, const mp_limb_t *mp,
+ mp_limb_t m_binv,
+ mp_size_t n, mp_limb_t *up, mp_limb_t *vp, mp_limb_t *tp)
+{
+ add_add_mul_mod (tp, up, vp, mp, n, M->a[0][0], M->a[0][1]);
+ add_add_mul_mod (tp + n + 1, up, vp, mp, n, M->a[1][0], M->a[1][1]);
+
+ /* Reduce to n limbs, by multiplying with B^-1 (mod m) */
+ redc_1 (up, tp, mp, n, m_binv);
+ redc_1 (vp, tp + n + 1, mp, n, m_binv);
+}
+
+/* Input in range 0 <= T < 2M, output in range 0 <= R < 2M. */
+static void
+redc_bits (mp_limb_t *rp, mp_limb_t *tp, const mp_limb_t *mp,
+ mp_size_t n, mp_limb_t m_binv, mp_bitcnt_t k)
+{
+ mp_limb_t hi;
+
+ for (; k >= GMP_NUMB_BITS; k -= GMP_NUMB_BITS)
+ {
+ hi = mpn_addmul_1 (tp, mp, n, -m_binv * tp[0]);
+ mpn_copyi (tp, tp + 1, n - 1);
+ tp[n-1] = hi;
+ }
+ if (k > 0)
+ {
+ mp_limb_t mask = ((mp_limb_t) 2 << (k-1)) - 1;
+ mp_limb_t q = (-m_binv * tp[0]) & mask;
+ hi = mpn_addmul_1 (tp, mp, n, q);
+ mpn_rshift (rp, tp, n, k);
+ rp[n-1] |= hi << (GMP_NUMB_BITS - k);
+ }
+}
+
+static int
+one_p (const mp_limb_t *p, mp_size_t n)
+{
+ mp_limb_t diff = p[0] ^ 1;
+ mp_size_t i;
+ for (i = 1; i < n; i++)
+ diff |= p[i];
+
+ return diff == 0;
+}
-/* FIXME: Could use mpn_sec_invert (in GMP-6), but with a bit more
- scratch need since it doesn't precompute (m+1)/2. */
+/* 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}.
+ Based on the algorithm in https://eprint.iacr.org/2019/266.pdf */
void
ecc_mod_inv (const struct ecc_modulo *m,
- mp_limb_t *vp, const mp_limb_t *in_ap,
+ mp_limb_t *up, const mp_limb_t *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;
+ mp_bitcnt_t shift, delta, count;
+ mp_limb_t cy;
+ struct matrix M;
+
+ /* Total scratch need 4*(n+1) for fp, gp, tp, n for v, total 5*n+4 */
+#define fp scratch
+#define gp (scratch + n + 1)
+#define vp (scratch + 2*n + 2)
+#define tp (scratch + 3*n + 2)
+
+#define mp (m->m)
+
+ /* Input should satisfy a < 2m, since that is what the iteration
+ count is based on. */
+ assert ((ap[n-1] >> 1) <= mp[n-1]);
- /* Maintain
+ mpn_copyi (fp, mp, n); fp[n] = 0; gp[n] = 0;
- a = u * orig_a (mod m)
- b = v * orig_a (mod m)
+ if (m->invert_power)
+ /* Multiply by suitable pre-computed power to get inverse in redc
+ representation. */
+ ecc_mod_mul (m, gp, ap, m->invert_power, tp);
+ else {
+ mpn_copyi (gp, ap, n);
- and b odd at all times. Initially,
+ shift = GMP_NUMB_BITS*((m->invert_count + STEPS - 1) / STEPS) - m->invert_count;
- a = a_orig, u = 1
- b = m, v = 0
- */
+ /* Premultiply a by 2^{- shift} mod m */
+ redc_bits (gp, gp, mp, n, m->binv, shift);
+ }
- assert (ap != vp);
+ /* Maintain invariant
- 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);
+ a * u = 2^{-count} B^{ceil(count/STEPS)} f (mod m)
+ a * v = 2^{-count} B^{ceil(count/STEPS)} g (mod m)
+ */
+ mpn_zero (up, n);
+ mpn_zero (vp+1, n-1); vp[0] = 1;
- for (i = m->bit_size + GMP_NUMB_BITS * n; i-- > 0; )
+ for (delta = 1, count = m->invert_count; count > STEPS; count -= STEPS)
{
- 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 (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 (cy == 0);
-
- cy = mpn_rshift (ap, ap, n, 1);
- assert (cy == 0);
- cy = mpn_rshift (up, up, n, 1);
- cy = mpn_cnd_add_n (cy, up, up, m->mp1h, n);
- assert (cy == 0);
+ delta = steps (&M, STEPS, delta, fp[0], gp[0]);
+ matrix_vector_mul (&M, STEPS, n+1, fp, gp, tp);
+ matrix_vector_mul_mod (&M, mp, m->binv, n, up, vp, tp);
}
- assert ( (ap[0] | ap[n-1]) == 0);
-#undef ap
-#undef bp
-#undef up
+ delta = steps (&M, count, delta, fp[0], gp[0]);
+ matrix_vector_mul (&M, count, n+1, fp, gp, tp);
+ /* Only compute u, we don't need v. */
+ add_add_mul_mod (tp, up, vp, mp, n, M.a[0][0], M.a[0][1]);
+ redc_1 (up, tp, mp, n, m->binv);
+
+ assert (sec_zero_p (gp, n+1));
+
+ /* Now f = ±1 (if the inverse exists), and a * u = f (mod m) */
+ cy = cnd_neg (fp[n] >> (GMP_LIMB_BITS - 1), up, up, n);
+ /* Make u non-negative */
+ cy -= mpn_cnd_add_n (cy, up, up, mp, n);
+
+ cnd_neg (fp[n] >> (GMP_LIMB_BITS - 1), fp, fp, n + 1);
+ cy = one_p (fp, n+1);
+ /* Set to zero if invert doesn't exist. */
+ while (n > 0)
+ up[--n] &= - cy;
}
ECC_SECP192R1_INV_ITCH,
ECC_SECP192R1_SQRT_ITCH,
0,
+ ECC_INVP_COUNT,
+ ECC_BINVP,
ecc_p,
ecc_Bmodp,
ecc_Bmodp_shifted,
ecc_Bm2p,
ecc_redc_ppm1,
- ecc_pp1h,
+ NULL,
ecc_secp192r1_modp,
ecc_secp192r1_modp,
ECC_MOD_INV_ITCH (ECC_LIMB_SIZE),
0,
0,
+ ECC_INVQ_COUNT,
+ ECC_BINVQ,
ecc_q,
ecc_Bmodq,
ecc_Bmodq_shifted,
ecc_Bm2q,
NULL,
- ecc_qp1h,
+ NULL,
ecc_mod,
ecc_mod,
#undef tp
}
-#define ECC_SECP224R1_INV_ITCH (4*ECC_LIMB_SIZE)
-
-static void
-ecc_secp224r1_inv (const struct ecc_modulo *p,
- mp_limb_t *rp, const mp_limb_t *ap,
- mp_limb_t *scratch)
-{
-#define a96m1 scratch
-#define tp (scratch + ECC_LIMB_SIZE)
-
- /* Compute a^{p - 2}, with
-
- p-2 = 2^{224} - 2^{96} - 1
- = 2^{97}(2^{127} - 1) + 2^{96} - 1
-
- This addition chain needs 97 squarings and one multiply in
- addition to ecc_mod_pow_127m1, for a total of 223 squarings and
- 13 multiplies.
- */
- ecc_mod_pow_127m1 (p, rp, a96m1, ap, tp);
- ecc_mod_pow_2k_mul (p, rp, rp, 97, a96m1, tp); /* a^{2^{224} - 2^{96} - 1 */
-
-#undef a96m1
-#undef tp
-}
-
#define ECC_SECP224R1_SQRT_ITCH (5*ECC_LIMB_SIZE)
static int
ECC_LIMB_SIZE,
ECC_BMODP_SIZE,
-ECC_REDC_SIZE,
- ECC_SECP224R1_INV_ITCH,
+ ECC_MOD_INV_ITCH(ECC_LIMB_SIZE),
ECC_SECP224R1_SQRT_ITCH,
0,
+ ECC_INVP_COUNT,
+ ECC_BINVP,
ecc_p,
ecc_Bmodp,
ecc_Bmodp_shifted,
ecc_Bm2p,
ecc_redc_ppm1,
- ecc_pp1h,
+ USE_REDC ? ecc_invp_redc_power : NULL,
ecc_secp224r1_modp,
USE_REDC ? ecc_secp224r1_redc : ecc_secp224r1_modp,
- ecc_secp224r1_inv,
+ ecc_mod_inv,
ecc_secp224r1_sqrt,
NULL,
},
ECC_MOD_INV_ITCH (ECC_LIMB_SIZE),
0,
0,
+ ECC_INVQ_COUNT,
+ ECC_BINVQ,
ecc_q,
ecc_Bmodq,
ecc_Bmodq_shifted,
ecc_Bm2q,
NULL,
- ecc_qp1h,
+ NULL,
ecc_mod,
ecc_mod,
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_SECP224R1_INV_ITCH),
+ ECC_J_TO_A_ITCH(ECC_LIMB_SIZE, ECC_MOD_INV_ITCH (ECC_LIMB_SIZE)),
ecc_add_jja,
ecc_add_jjj,
#error Unsupported parameters
#endif
-#define ECC_SECP256R1_INV_ITCH (4*ECC_LIMB_SIZE)
-
-static void
-ecc_secp256r1_inv (const struct ecc_modulo *p,
- mp_limb_t *rp, const mp_limb_t *ap,
- mp_limb_t *scratch)
-{
-#define a5m1 scratch
-#define t0 (scratch + ECC_LIMB_SIZE)
-#define a15m1 t0
-#define a32m1 a5m1
-#define tp (scratch + 2*ECC_LIMB_SIZE)
-/*
- Addition chain for p - 2 = 2^{256} - 2^{224} + 2^{192} + 2^{96} - 3
-
- 2^5 - 1 = 1 + 2 (2^4 - 1) = 1 + 2 (2^2+1)(2 + 1) 4 S + 3 M
- 2^{15} - 1 = (2^5 - 1) (1 + 2^5 (1 + 2^5) 10 S + 2 M
- 2^{16} - 1 = 1 + 2 (2^{15} - 1) S + M
- 2^{32} - 1 = (2^{16} + 1) (2^{16} - 1) 16 S + M
- 2^{64} - 2^{32} + 1 = 2^{32} (2^{32} - 1) + 1 32 S + M
- 2^{192} - 2^{160} + 2^{128} + 2^{32} - 1
- = 2^{128} (2^{64} - 2^{32} + 1) + 2^{32} - 1 128 S + M
- 2^{224} - 2^{192} + 2^{160} + 2^{64} - 1
- = 2^{32} (...) + 2^{32} - 1 32 S + M
- 2^{239} - 2^{207} + 2^{175} + 2^{79} - 1
- = 2^{15} (...) + 2^{15} - 1 15 S + M
- 2^{254} - 2^{222} + 2^{190} + 2^{94} - 1
- = 2^{15} (...) + 2^{15} - 1 15 S + M
- p - 2 = 2^2 (...) + 1 2 S M
- ---------------
- 255 S + 13 M
- */
- ecc_mod_sqr (p, rp, ap, tp); /* a^2 */
- ecc_mod_mul (p, rp, rp, ap, tp); /* a^3 */
- ecc_mod_pow_2kp1 (p, t0, rp, 2, tp); /* a^{2^4 - 1} */
- ecc_mod_sqr (p, rp, t0, tp); /* a^{2^5 - 2} */
- ecc_mod_mul (p, a5m1, rp, ap, tp); /* a^{2^5 - 1}, a5m1 */
-
- ecc_mod_pow_2kp1 (p, rp, a5m1, 5, tp); /* a^{2^{10} - 1, a5m1*/
- ecc_mod_pow_2k_mul (p, a15m1, rp, 5, a5m1, tp); /* a^{2^{15} - 1}, a5m1 a15m1 */
- ecc_mod_sqr (p, rp, a15m1, tp); /* a^{2^{16} - 2}, a15m1 */
- ecc_mod_mul (p, rp, rp, ap, tp); /* a^{2^{16} - 1}, a15m1 */
- ecc_mod_pow_2kp1 (p, a32m1, rp, 16, tp); /* a^{2^{32} - 1}, a15m1, a32m1 */
-
- ecc_mod_pow_2k_mul (p, rp, a32m1, 32, ap, tp);/* a^{2^{64} - 2^{32} + 1 */
- ecc_mod_pow_2k_mul (p, rp, rp, 128, a32m1, tp); /* a^{2^{192} - 2^{160} + 2^{128} + 2^{32} - 1} */
- ecc_mod_pow_2k_mul (p, rp, rp, 32, a32m1, tp);/* a^{2^{224} - 2^{192} + 2^{160} + 2^{64} - 1} */
- ecc_mod_pow_2k_mul (p, rp, rp, 15, a15m1, tp);/* a^{2^{239} - 2^{207} + 2^{175} + 2^{79} - 1} */
- ecc_mod_pow_2k_mul (p, rp, rp, 15, a15m1, tp);/* a^{2^{254} - 2^{222} + 2^{190} + 2^{94} - 1} */
- ecc_mod_pow_2k_mul (p, rp, rp, 2, ap, tp); /* a^{2^{256} - 2^{224} + 2^{192} + 2^{96} - 3} */
-
-#undef a5m1
-#undef t0
-#undef a15m1
-#undef a32m1
-#undef tp
-}
-
/* To guarantee that inputs to ecc_mod_zero_p are in the required range. */
#if ECC_LIMB_SIZE * GMP_NUMB_BITS != 256
#error Unsupported limb size
ECC_LIMB_SIZE,
ECC_BMODP_SIZE,
ECC_REDC_SIZE,
- ECC_SECP256R1_INV_ITCH,
+ ECC_MOD_INV_ITCH (ECC_LIMB_SIZE),
ECC_SECP256R1_SQRT_ITCH,
0,
+ ECC_INVP_COUNT,
+ ECC_BINVP,
ecc_p,
ecc_Bmodp,
ecc_Bmodp_shifted,
ecc_Bm2p,
ecc_redc_ppm1,
- ecc_pp1h,
+ USE_REDC ? ecc_invp_redc_power : NULL,
ecc_secp256r1_modp,
USE_REDC ? ecc_secp256r1_redc : ecc_secp256r1_modp,
- ecc_secp256r1_inv,
+ ecc_mod_inv,
ecc_secp256r1_sqrt,
NULL,
},
ECC_MOD_INV_ITCH (ECC_LIMB_SIZE),
0,
0,
+ ECC_INVQ_COUNT,
+ ECC_BINVQ,
ecc_q,
ecc_Bmodq,
ecc_Bmodq_shifted,
ecc_Bm2q,
NULL,
- ecc_qp1h,
-
+ NULL,
ecc_secp256r1_modq,
ecc_secp256r1_modq,
ecc_mod_inv,
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_SECP256R1_INV_ITCH),
+ ECC_J_TO_A_ITCH(ECC_LIMB_SIZE, ECC_MOD_INV_ITCH (ECC_LIMB_SIZE)),
ecc_add_jja,
ecc_add_jjj,
#undef tp
}
-#define ECC_SECP384R1_INV_ITCH (6*ECC_LIMB_SIZE)
-
-static void
-ecc_secp384r1_inv (const struct ecc_modulo *p,
- mp_limb_t *rp, const mp_limb_t *ap,
- mp_limb_t *scratch)
-{
- /*
- Addition chain for
-
- p - 2 = 2^{384} - 2^{128} - 2^{96} + 2^{32} - 3
-
- Start with
-
- a^{2^{288} - 2^{32} - 1}
-
- and then use
-
- 2^{382} - 2^{126} - 2^{94} + 2^{30} - 1
- = 2^{94} (2^{288} - 2^{32} - 1) + 2^{30} - 1
-
- 2^{384} - 2^{128} - 2^{96} + 2^{32} - 3
- = 2^2 (2^{382} - 2^{126} - 2^{94} + 2^{30} - 1) + 1
-
- This addition chain needs 96 additional squarings and 2
- multiplies, for a total of 383 squarings and 14 multiplies.
- */
-
-#define a30m1 scratch
-#define tp (scratch + ECC_LIMB_SIZE)
-
- ecc_mod_pow_288m32m1 (p, rp, a30m1, ap, tp);
- ecc_mod_pow_2k_mul (p, rp, rp, 94, a30m1, tp); /* a^{2^{382} - 2^{126} - 2^{94} + 2^{30} - 1 */
- ecc_mod_pow_2k_mul (p, rp, rp, 2, ap, tp);
-
-#undef a30m1
-#undef tp
-}
-
/* To guarantee that inputs to ecc_mod_zero_p are in the required range. */
#if ECC_LIMB_SIZE * GMP_NUMB_BITS != 384
#error Unsupported limb size
ECC_LIMB_SIZE,
ECC_BMODP_SIZE,
ECC_REDC_SIZE,
- ECC_SECP384R1_INV_ITCH,
+ ECC_MOD_INV_ITCH(ECC_LIMB_SIZE),
ECC_SECP384R1_SQRT_ITCH,
0,
+ ECC_INVP_COUNT,
+ ECC_BINVP,
ecc_p,
ecc_Bmodp,
ecc_Bmodp_shifted,
ecc_Bm2p,
ecc_redc_ppm1,
- ecc_pp1h,
+ NULL,
ecc_secp384r1_modp,
ecc_secp384r1_modp,
- ecc_secp384r1_inv,
+ ecc_mod_inv,
ecc_secp384r1_sqrt,
NULL,
},
ECC_MOD_INV_ITCH (ECC_LIMB_SIZE),
0,
0,
+ ECC_INVQ_COUNT,
+ ECC_BINVQ,
ecc_q,
ecc_Bmodq,
ecc_Bmodq_shifted,
ecc_Bm2q,
NULL,
- ecc_qp1h,
+ NULL,
ecc_mod,
ecc_mod,
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_SECP384R1_INV_ITCH),
+ ECC_J_TO_A_ITCH(ECC_LIMB_SIZE, ECC_MOD_INV_ITCH (ECC_LIMB_SIZE)),
ecc_add_jja,
ecc_add_jjj,
}
#endif
-#define ECC_SECP521R1_INV_ITCH (3*ECC_LIMB_SIZE)
-
-static void
-ecc_secp521r1_inv (const struct ecc_modulo *p,
- mp_limb_t *rp, const mp_limb_t *ap,
- mp_limb_t *scratch)
-{
-#define t0 scratch
-#define tp (scratch + ECC_LIMB_SIZE)
-
- /* Addition chain for p - 2:
-
- 2^{521} - 3
- = 1 + 2^2(2^519 - 1)
- = 1 + 2^2(1 + 2 (2^518 - 1)
- = 1 + 2^2(1 + 2 (2^259 + 1) (1 + 2(2^258 - 1)))
- = 1 + 2^2(1 + 2 (2^259 + 1) (1 + 2(2^129 + 1) (2^129 - 1)))
- = 1 + 2^2(1 + 2 (2^259 + 1) (1 + 2(2^129 + 1) (1 + 2 (2^128 - 1))))
-
- where
-
- 2^{128} - 1 = (2^64 + 1) (2^32+1) (2^16 + 1) (2^8 + 1) (2^4 + 1) (2^2 + 1) (2 + 1)
-
- This addition chain needs 520 squarings and 13 multiplies.
- */
-
- ecc_mod_sqr (p, rp, ap, tp); /* a^2 */
- ecc_mod_mul (p, rp, ap, rp, tp); /* a^3 = a^{2^2 - 1} */
- ecc_mod_pow_2kp1 (p, t0, rp, 2, tp); /* a^15 = a^{2^4 - 1} */
- ecc_mod_pow_2kp1 (p, rp, t0, 4, tp); /* a^{2^8 - 1} */
- ecc_mod_pow_2kp1 (p, t0, rp, 8, tp); /* a^{2^16 - 1} */
- ecc_mod_pow_2kp1 (p, rp, t0, 16, tp); /* a^{2^32 - 1} */
- ecc_mod_pow_2kp1 (p, t0, rp, 32, tp); /* a^{2^64 - 1} */
- ecc_mod_pow_2kp1 (p, rp, t0, 64, tp); /* a^{2^128 - 1} */
- ecc_mod_sqr (p, rp, rp, tp); /* a^{2^129 - 2} */
- ecc_mod_mul (p, rp, rp, ap, tp); /* a^{2^129 - 1} */
- ecc_mod_pow_2kp1 (p, t0, rp, 129, tp);/* a^{2^258 - 1} */
- ecc_mod_sqr (p, rp, t0, tp); /* a^{2^259 - 2} */
- ecc_mod_mul (p, rp, rp, ap, tp); /* a^{2^259 - 1} */
- ecc_mod_pow_2kp1 (p, t0, rp, 259, tp);/* a^{2^518 - 1} */
- ecc_mod_sqr (p, rp, t0, tp); /* a^{2^519 - 2} */
- ecc_mod_mul (p, rp, rp, ap, tp); /* a^{2^519 - 1} */
- ecc_mod_sqr (p, rp, rp, tp); /* a^{2^520 - 2} */
- ecc_mod_sqr (p, rp, rp, tp); /* a^{2^521 - 4} */
- ecc_mod_mul (p, rp, rp, ap, tp); /* a^{2^521 - 3} */
-}
-
#define ECC_SECP521R1_SQRT_ITCH (2*ECC_LIMB_SIZE)
static int
ECC_LIMB_SIZE,
ECC_BMODP_SIZE,
ECC_REDC_SIZE,
- ECC_SECP521R1_INV_ITCH,
+ ECC_MOD_INV_ITCH(ECC_LIMB_SIZE),
ECC_SECP521R1_SQRT_ITCH,
0,
+ ECC_INVP_COUNT,
+ ECC_BINVP,
ecc_p,
ecc_Bmodp,
ecc_Bmodp_shifted,
ecc_Bm2p,
ecc_redc_ppm1,
- ecc_pp1h,
+ NULL,
ecc_secp521r1_modp,
ecc_secp521r1_modp,
- ecc_secp521r1_inv,
+ ecc_mod_inv,
ecc_secp521r1_sqrt,
NULL,
},
ECC_MOD_INV_ITCH (ECC_LIMB_SIZE),
0,
0,
+ ECC_INVQ_COUNT,
+ ECC_BINVQ,
ecc_q,
ecc_Bmodq,
ecc_Bmodq_shifted,
ecc_Bm2q,
NULL,
- ecc_qp1h,
+ NULL,
ecc_mod,
ecc_mod,
ECC_DUP_JJ_ITCH (ECC_LIMB_SIZE),
ECC_MUL_A_ITCH (ECC_LIMB_SIZE),
ECC_MUL_G_ITCH (ECC_LIMB_SIZE),
+#if 0
ECC_J_TO_A_ITCH(ECC_LIMB_SIZE, ECC_SECP521R1_INV_ITCH),
-
+#else
+ ECC_J_TO_A_ITCH(ECC_LIMB_SIZE, ECC_MOD_INV_ITCH(ECC_LIMB_SIZE)),
+#endif
ecc_add_jja,
ecc_add_jjj,
ecc_dup_jj,
abort();
}
+static unsigned
+invert_count (unsigned bit_size) {
+ /* See https://eprint.iacr.org/2019/266.pdf, Theorem 11.2 */
+ return (49 * bit_size + 57) / 17;
+}
+
static void
output_modulo (const char *name, const mpz_t x,
unsigned size, unsigned bits_per_limb)
{
unsigned bit_size;
+ unsigned count;
+ unsigned invert_iterations;
int shift;
char buf[20];
mpz_t t;
/ bits_per_limb));
bit_size = mpz_sizeinbase (x, 2);
-
shift = size * bits_per_limb - bit_size;
assert (shift >= 0);
+
+ /* If most significant bit of top limb is zero, input may be one bit
+ larger than the modulo. */
+ count = invert_count (bit_size + (shift > 0));
+ printf("#define ECC_INV%s_COUNT %u\n", buf, count);
+
+ mpz_set_ui (t, 0);
+ mpz_setbit (t, bits_per_limb);
+ mpz_invert (t, x, t);
+
+ printf ("#define ECC_BINV%s 0x", buf);
+ mpz_out_str (stdout, 16, t);
+ printf("%s\n", bits_per_limb > 32 ? "ULL" : "UL");
+
+ invert_iterations = (count + bits_per_limb - 3) / (bits_per_limb - 2);
+ mpz_set_ui (t, 2);
+ mpz_powm_ui (t, t, bits_per_limb * (invert_iterations + size) - count, x);
+ mpz_invert (t, t, x);
+ snprintf (buf, sizeof(buf), "ecc_inv%s_redc_power", name);
+ output_bignum (buf, t, size, bits_per_limb);
+
if (shift > 0)
{
mpz_set_ui (t, 0);
output_bignum ("ecc_b", ecc->b, limb_size, bits_per_limb);
- mpz_add_ui (t, ecc->p, 1);
- mpz_fdiv_q_2exp (t, t, 1);
- output_bignum ("ecc_pp1h", t, limb_size, bits_per_limb);
-
- mpz_add_ui (t, ecc->q, 1);
- mpz_fdiv_q_2exp (t, t, 1);
- output_bignum ("ecc_qp1h", t, limb_size, bits_per_limb);
-
/* Trailing zeros in p+1 correspond to trailing ones in p. */
redc_limbs = mpz_scan0 (ecc->p, 0) / bits_per_limb;
if (redc_limbs > 0)
mp_limb_t *scratch;
unsigned j;
mpz_t r;
+ mp_bitcnt_t bits;
mpz_init (r);
+ bits = m->bit_size + (m->size * GMP_NUMB_BITS > m->bit_size);
+
a = xalloc_limbs (m->size);
ai = xalloc_limbs (m->size);
ref = xalloc_limbs (m->size);;
fprintf (stderr, " (bad)\n");
abort ();
}
-
+
for (j = 0; j < COUNT; j++)
{
if (j & 1)
- mpz_rrandomb (r, rands, m->size * GMP_NUMB_BITS);
+ mpz_rrandomb (r, rands, bits);
else
- mpz_urandomb (r, rands, m->size * GMP_NUMB_BITS);
+ mpz_urandomb (r, rands, bits);
mpz_limbs_copy (a, r, m->size);
+ if ((a[m->size - 1] >> 1) > m->m[m->size - 1])
+ a[m->size - 1] = 2*m->m[m->size - 1] + 1;
if (!ref_modinv (ref, a, m->m, m->size, use_redc))
{
fprintf (stderr, "%s->invert failed (test %u, bit size %u):\n",
name, j, m->bit_size);
fprintf (stderr, "a = ");
- mpz_out_str (stderr, 16, r);
+ mpn_out_str (stderr, 16, a, m->size);
fprintf (stderr, "\np = ");
mpn_out_str (stderr, 16, m->m, m->size);
fprintf (stderr, "\nt = ");
projects / thirdparty / nettle.git / commitdiff
