diff options
Diffstat (limited to 'chip/g/dcrypto/bn.c')
-rw-r--r-- | chip/g/dcrypto/bn.c | 1244 |
1 files changed, 1244 insertions, 0 deletions
diff --git a/chip/g/dcrypto/bn.c b/chip/g/dcrypto/bn.c new file mode 100644 index 0000000000..94aafa1799 --- /dev/null +++ b/chip/g/dcrypto/bn.c @@ -0,0 +1,1244 @@ +/* Copyright 2015 The Chromium OS Authors. All rights reserved. + * Use of this source code is governed by a BSD-style license that can be + * found in the LICENSE file. + */ + +#ifdef PRINT_PRIMES +#include "console.h" +#endif + +#include "dcrypto.h" +#include "internal.h" + +#include "trng.h" + +#include "cryptoc/util.h" + +#include <assert.h> + +#ifdef CONFIG_WATCHDOG +extern void watchdog_reload(void); +#else +static inline void watchdog_reload(void) { } +#endif + +void bn_init(struct LITE_BIGNUM *b, void *buf, size_t len) +{ + DCRYPTO_bn_wrap(b, buf, len); + always_memset(buf, 0x00, len); +} + +void DCRYPTO_bn_wrap(struct LITE_BIGNUM *b, void *buf, size_t len) +{ + /* Only word-multiple sized buffers accepted. */ + assert((len & 0x3) == 0); + b->dmax = len / LITE_BN_BYTES; + b->d = (struct access_helper *) buf; +} + +int bn_eq(const struct LITE_BIGNUM *a, const struct LITE_BIGNUM *b) +{ + int i; + uint32_t top = 0; + + for (i = a->dmax - 1; i > b->dmax - 1; --i) + top |= BN_DIGIT(a, i); + if (top) + return 0; + + for (i = b->dmax - 1; i > a->dmax - 1; --i) + top |= BN_DIGIT(b, i); + if (top) + return 0; + + for (i = MIN(a->dmax, b->dmax) - 1; i >= 0; --i) + if (BN_DIGIT(a, i) != BN_DIGIT(b, i)) + return 0; + + return 1; +} + +static void bn_copy(struct LITE_BIGNUM *dst, const struct LITE_BIGNUM *src) +{ + dst->dmax = src->dmax; + memcpy(dst->d, src->d, bn_size(dst)); +} + +int bn_check_topbit(const struct LITE_BIGNUM *N) +{ + return BN_DIGIT(N, N->dmax - 1) >> 31; +} + +/* a[n]. */ +int bn_is_bit_set(const struct LITE_BIGNUM *a, int n) +{ + int i, j; + + if (n < 0) + return 0; + + i = n / LITE_BN_BITS2; + j = n % LITE_BN_BITS2; + if (a->dmax <= i) + return 0; + + return (BN_DIGIT(a, i) >> j) & 1; +} + +static int bn_set_bit(const struct LITE_BIGNUM *a, int n) +{ + int i, j; + + if (n < 0) + return 0; + + i = n / LITE_BN_BITS2; + j = n % LITE_BN_BITS2; + if (a->dmax <= i) + return 0; + + BN_DIGIT(a, i) |= 1 << j; + return 1; +} + +/* a[] >= b[]. */ +/* TODO(ngm): constant time. */ +static int bn_gte(const struct LITE_BIGNUM *a, const struct LITE_BIGNUM *b) +{ + int i; + uint32_t top = 0; + + for (i = a->dmax - 1; i > b->dmax - 1; --i) + top |= BN_DIGIT(a, i); + if (top) + return 1; + + for (i = b->dmax - 1; i > a->dmax - 1; --i) + top |= BN_DIGIT(b, i); + if (top) + return 0; + + for (i = MIN(a->dmax, b->dmax) - 1; + BN_DIGIT(a, i) == BN_DIGIT(b, i) && i > 0; --i) + ; + return BN_DIGIT(a, i) >= BN_DIGIT(b, i); +} + +/* c[] = c[] - a[], assumes c > a. */ +uint32_t bn_sub(struct LITE_BIGNUM *c, const struct LITE_BIGNUM *a) +{ + int64_t A = 0; + int i; + + for (i = 0; i < a->dmax; i++) { + A += (uint64_t) BN_DIGIT(c, i) - BN_DIGIT(a, i); + BN_DIGIT(c, i) = (uint32_t) A; + A >>= 32; + } + + for (; A && i < c->dmax; i++) { + A += (uint64_t) BN_DIGIT(c, i); + BN_DIGIT(c, i) = (uint32_t) A; + A >>= 32; + } + + return (uint32_t) A; /* 0 or -1. */ +} + +/* c[] = c[] - a[], negative numbers in 2's complement representation. */ +/* Returns borrow bit. */ +static uint32_t bn_signed_sub(struct LITE_BIGNUM *c, int *c_neg, + const struct LITE_BIGNUM *a, int a_neg) +{ + uint32_t carry = 0; + uint64_t A = 1; + int i; + + for (i = 0; i < a->dmax; ++i) { + A += (uint64_t) BN_DIGIT(c, i) + ~BN_DIGIT(a, i); + BN_DIGIT(c, i) = (uint32_t) A; + A >>= 32; + } + + for (; i < c->dmax; ++i) { + A += (uint64_t) BN_DIGIT(c, i) + 0xFFFFFFFF; + BN_DIGIT(c, i) = (uint32_t) A; + A >>= 32; + } + + A &= 0x01; + carry = (!*c_neg && a_neg && A) || (*c_neg && !a_neg && !A); + *c_neg = carry ? *c_neg : (*c_neg + !a_neg + A) & 0x01; + return carry; +} + +/* c[] = c[] + a[]. */ +uint32_t bn_add(struct LITE_BIGNUM *c, const struct LITE_BIGNUM *a) +{ + uint64_t A = 0; + int i; + + for (i = 0; i < a->dmax; ++i) { + A += (uint64_t) BN_DIGIT(c, i) + BN_DIGIT(a, i); + BN_DIGIT(c, i) = (uint32_t) A; + A >>= 32; + } + + for (; A && i < c->dmax; ++i) { + A += (uint64_t) BN_DIGIT(c, i); + BN_DIGIT(c, i) = (uint32_t) A; + A >>= 32; + } + + return (uint32_t) A; /* 0 or 1. */ +} + +/* c[] = c[] + a[], negative numbers in 2's complement representation. */ +/* Returns carry bit. */ +static uint32_t bn_signed_add(struct LITE_BIGNUM *c, int *c_neg, + const struct LITE_BIGNUM *a, int a_neg) +{ + uint32_t A = bn_add(c, a); + uint32_t carry; + + carry = (!*c_neg && !a_neg && A) || (*c_neg && a_neg && !A); + *c_neg = carry ? *c_neg : (*c_neg + a_neg + A) & 0x01; + return carry; +} + +/* r[] <<= 1. */ +static uint32_t bn_lshift(struct LITE_BIGNUM *r) +{ + int i; + uint32_t w; + uint32_t carry = 0; + + for (i = 0; i < r->dmax; i++) { + w = (BN_DIGIT(r, i) << 1) | carry; + carry = BN_DIGIT(r, i) >> 31; + BN_DIGIT(r, i) = w; + } + return carry; +} + +/* r[] >>= 1. Handles 2's complement negative numbers. */ +static void bn_rshift(struct LITE_BIGNUM *r, uint32_t carry, uint32_t neg) +{ + int i; + uint32_t ones = ~0; + uint32_t highbit = (!carry && neg) || (carry && !neg); + + for (i = 0; i < r->dmax - 1; ++i) { + uint32_t accu; + + ones &= BN_DIGIT(r, i); + accu = (BN_DIGIT(r, i) >> 1); + accu |= (BN_DIGIT(r, i + 1) << (LITE_BN_BITS2 - 1)); + BN_DIGIT(r, i) = accu; + } + ones &= BN_DIGIT(r, i); + BN_DIGIT(r, i) = (BN_DIGIT(r, i) >> 1) | + (highbit << (LITE_BN_BITS2 - 1)); + + if (ones == ~0 && highbit && neg) + memset(r->d, 0x00, bn_size(r)); /* -1 >> 1 = 0. */ +} + +/* Montgomery c[] += a * b[] / R % N. */ +/* TODO(ngm): constant time. */ +static void bn_mont_mul_add(struct LITE_BIGNUM *c, const uint32_t a, + const struct LITE_BIGNUM *b, const uint32_t nprime, + const struct LITE_BIGNUM *N) +{ + uint32_t A, B, d0; + int i; + + { + register uint64_t tmp; + + tmp = BN_DIGIT(c, 0) + (uint64_t) a * BN_DIGIT(b, 0); + A = tmp >> 32; + d0 = (uint32_t) tmp * (uint32_t) nprime; + tmp = (uint32_t)tmp + (uint64_t) d0 * BN_DIGIT(N, 0); + B = tmp >> 32; + } + + for (i = 0; i < N->dmax - 1;) { + register uint64_t tmp; + + tmp = A + (uint64_t) a * BN_DIGIT(b, i + 1) + + BN_DIGIT(c, i + 1); + A = tmp >> 32; + tmp = B + (uint64_t) d0 * BN_DIGIT(N, i + 1) + (uint32_t) tmp; + BN_DIGIT(c, i) = (uint32_t) tmp; + B = tmp >> 32; + ++i; + } + + { + uint64_t tmp = (uint64_t) A + B; + + BN_DIGIT(c, i) = (uint32_t) tmp; + A = tmp >> 32; /* 0 or 1. */ + if (A) + bn_sub(c, N); + } +} + +/* Montgomery c[] = a[] * b[] / R % N. */ +static void bn_mont_mul(struct LITE_BIGNUM *c, const struct LITE_BIGNUM *a, + const struct LITE_BIGNUM *b, const uint32_t nprime, + const struct LITE_BIGNUM *N) +{ + int i; + + for (i = 0; i < N->dmax; i++) + BN_DIGIT(c, i) = 0; + + bn_mont_mul_add(c, a ? BN_DIGIT(a, 0) : 1, b, nprime, N); + for (i = 1; i < N->dmax; i++) + bn_mont_mul_add(c, a ? BN_DIGIT(a, i) : 0, b, nprime, N); +} + +/* Mongomery R * R % N, R = 1 << (1 + log2N). */ +/* TODO(ngm): constant time. */ +static void bn_compute_RR(struct LITE_BIGNUM *RR, const struct LITE_BIGNUM *N) +{ + int i; + + bn_sub(RR, N); /* R - N = R % N since R < 2N */ + + /* Repeat 2 * R % N, log2(R) times. */ + for (i = 0; i < N->dmax * LITE_BN_BITS2; i++) { + if (bn_lshift(RR)) + assert(bn_sub(RR, N) == -1); + if (bn_gte(RR, N)) + bn_sub(RR, N); + } +} + +/* Montgomery nprime = -1 / n0 % (2 ^ 32). */ +static uint32_t bn_compute_nprime(const uint32_t n0) +{ + int i; + uint32_t ninv = 1; + + /* Repeated Hensel lifting. */ + for (i = 0; i < 5; i++) + ninv *= 2 - (n0 * ninv); + + return ~ninv + 1; /* Two's complement. */ +} + +/* TODO(ngm): this implementation not timing or side-channel safe by + * any measure. */ +static void bn_modexp_internal(struct LITE_BIGNUM *output, + const struct LITE_BIGNUM *input, + const struct LITE_BIGNUM *exp, + const struct LITE_BIGNUM *N) +{ + int i; + uint32_t nprime; + uint32_t RR_buf[RSA_MAX_WORDS]; + uint32_t acc_buf[RSA_MAX_WORDS]; + uint32_t aR_buf[RSA_MAX_WORDS]; + + struct LITE_BIGNUM RR; + struct LITE_BIGNUM acc; + struct LITE_BIGNUM aR; + + bn_init(&RR, RR_buf, bn_size(N)); + bn_init(&acc, acc_buf, bn_size(N)); + bn_init(&aR, aR_buf, bn_size(N)); + + nprime = bn_compute_nprime(BN_DIGIT(N, 0)); + bn_compute_RR(&RR, N); + bn_mont_mul(&acc, NULL, &RR, nprime, N); /* R = 1 * RR / R % N */ + bn_mont_mul(&aR, input, &RR, nprime, N); /* aR = a * RR / R % N */ + + /* TODO(ngm): burn stack space and use windowing. */ + for (i = exp->dmax * LITE_BN_BITS2 - 1; i >= 0; i--) { + bn_mont_mul(output, &acc, &acc, nprime, N); + if (bn_is_bit_set(exp, i)) { + bn_mont_mul(&acc, output, &aR, nprime, N); + } else { + struct LITE_BIGNUM tmp = *output; + + *output = acc; + acc = tmp; + } + /* Poke the watchdog. + * TODO(ngm): may be unnecessary with + * a faster implementation. + */ + watchdog_reload(); + } + + bn_mont_mul(output, NULL, &acc, nprime, N); /* Convert out. */ + /* Copy to output buffer if necessary. */ + if (acc.d != (struct access_helper *) acc_buf) { + memcpy(acc.d, acc_buf, bn_size(output)); + *output = acc; + } + + /* TODO(ngm): constant time. */ + if (bn_sub(output, N)) + bn_add(output, N); /* Final reduce. */ + output->dmax = N->dmax; + + always_memset(RR_buf, 0, sizeof(RR_buf)); + always_memset(acc_buf, 0, sizeof(acc_buf)); + always_memset(aR_buf, 0, sizeof(aR_buf)); +} + +/* output = input ^ exp % N */ +int bn_modexp(struct LITE_BIGNUM *output, const struct LITE_BIGNUM *input, + const struct LITE_BIGNUM *exp, const struct LITE_BIGNUM *N) +{ +#ifndef CR50_NO_BN_ASM + if ((bn_bits(N) & 255) == 0) { + /* Use hardware support for standard key sizes. */ + return dcrypto_modexp(output, input, exp, N); + } +#endif + bn_modexp_internal(output, input, exp, N); + return 1; +} + +/* output = input ^ exp % N */ +int bn_modexp_word(struct LITE_BIGNUM *output, const struct LITE_BIGNUM *input, + uint32_t exp, const struct LITE_BIGNUM *N) +{ +#ifndef CR50_NO_BN_ASM + if ((bn_bits(N) & 255) == 0) { + /* Use hardware support for standard key sizes. */ + return dcrypto_modexp_word(output, input, exp, N); + } +#endif + { + struct LITE_BIGNUM pubexp; + + DCRYPTO_bn_wrap(&pubexp, &exp, sizeof(exp)); + bn_modexp_internal(output, input, &pubexp, N); + return 1; + } +} + +/* output = input ^ exp % N */ +int bn_modexp_blinded(struct LITE_BIGNUM *output, + const struct LITE_BIGNUM *input, + const struct LITE_BIGNUM *exp, + const struct LITE_BIGNUM *N, + uint32_t pubexp) +{ +#ifndef CR50_NO_BN_ASM + if ((bn_bits(N) & 255) == 0) { + /* Use hardware support for standard key sizes. */ + return dcrypto_modexp_blinded(output, input, exp, N, pubexp); + } +#endif + bn_modexp_internal(output, input, exp, N); + return 1; +} + +/* c[] += a * b[] */ +static uint32_t bn_mul_add(struct LITE_BIGNUM *c, uint32_t a, + const struct LITE_BIGNUM *b, uint32_t offset) +{ + int i; + uint64_t carry = 0; + + for (i = 0; i < b->dmax; i++) { + carry += BN_DIGIT(c, offset + i) + + (uint64_t) BN_DIGIT(b, i) * a; + BN_DIGIT(c, offset + i) = (uint32_t) carry; + carry >>= 32; + } + + return carry; +} + +/* c[] = a[] * b[] */ +void DCRYPTO_bn_mul(struct LITE_BIGNUM *c, const struct LITE_BIGNUM *a, + const struct LITE_BIGNUM *b) +{ + int i; + uint32_t carry = 0; + + memset(c->d, 0, bn_size(c)); + for (i = 0; i < a->dmax; i++) { + BN_DIGIT(c, i + b->dmax - 1) = carry; + carry = bn_mul_add(c, BN_DIGIT(a, i), b, i); + } + + BN_DIGIT(c, i + b->dmax - 1) = carry; +} + +/* c[] = a[] * b[] */ +static void bn_mul_ex(struct LITE_BIGNUM *c, + const struct LITE_BIGNUM *a, int a_len, + const struct LITE_BIGNUM *b) +{ + int i; + uint32_t carry = 0; + + memset(c->d, 0, bn_size(c)); + for (i = 0; i < a_len; i++) { + BN_DIGIT(c, i + b->dmax - 1) = carry; + carry = bn_mul_add(c, BN_DIGIT(a, i), b, i); + } + + BN_DIGIT(c, i + b->dmax - 1) = carry; +} + +static int bn_div_word_ex(struct LITE_BIGNUM *q, + struct LITE_BIGNUM *r, + const struct LITE_BIGNUM *u, int m, + uint32_t div) +{ + uint32_t rem = 0; + int i; + + for (i = m - 1; i >= 0; --i) { + uint64_t tmp = ((uint64_t)rem << 32) + BN_DIGIT(u, i); + uint32_t qd = tmp / div; + + BN_DIGIT(q, i) = qd; + rem = tmp - (uint64_t)qd * div; + } + + if (r != NULL) + BN_DIGIT(r, 0) = rem; + + return 1; +} + +/* + * Knuth's long division. + * + * Returns 0 on error. + * |u| >= |v| + * v[n-1] must not be 0 + * r gets |v| digits written to. + * q gets |u| - |v| + 1 digits written to. + */ +static int bn_div_ex(struct LITE_BIGNUM *q, + struct LITE_BIGNUM *r, + const struct LITE_BIGNUM *u, int m, + const struct LITE_BIGNUM *v, int n) +{ + uint32_t vtop; + int s, i, j; + uint32_t vn[RSA_MAX_WORDS]; /* Normalized v */ + uint32_t un[RSA_MAX_WORDS + 1]; /* Normalized u */ + + if (m < n || n <= 0) + return 0; + + vtop = BN_DIGIT(v, n - 1); + + if (vtop == 0) + return 0; + + if (n == 1) + return bn_div_word_ex(q, r, u, m, vtop); + + /* Compute shift factor to make v have high bit set */ + s = 0; + while ((vtop & 0x80000000) == 0) { + s = s + 1; + vtop = vtop << 1; + } + + /* Normalize u and v into un and vn. + * Note un always gains a leading digit + */ + if (s != 0) { + for (i = n - 1; i > 0; i--) + vn[i] = (BN_DIGIT(v, i) << s) | + (BN_DIGIT(v, i - 1) >> (32 - s)); + vn[0] = BN_DIGIT(v, 0) << s; + + un[m] = BN_DIGIT(u, m - 1) >> (32 - s); + for (i = m - 1; i > 0; i--) + un[i] = (BN_DIGIT(u, i) << s) | + (BN_DIGIT(u, i - 1) >> (32 - s)); + un[0] = BN_DIGIT(u, 0) << s; + } else { + for (i = 0; i < n; ++i) + vn[i] = BN_DIGIT(v, i); + for (i = 0; i < m; ++i) + un[i] = BN_DIGIT(u, i); + un[m] = 0; + } + + /* Main loop, reducing un digit by digit */ + for (j = m - n; j >= 0; j--) { + uint32_t qd; + int64_t t, k; + + /* Estimate quotient digit */ + if (un[j + n] == vn[n - 1]) { + /* Maxed out */ + qd = 0xFFFFFFFF; + } else { + /* Fine tune estimate */ + uint64_t rhat = ((uint64_t)un[j + n] << 32) + + un[j + n - 1]; + + qd = rhat / vn[n - 1]; + rhat = rhat - (uint64_t)qd * vn[n - 1]; + while ((rhat >> 32) == 0 && + (uint64_t)qd * vn[n - 2] > + (rhat << 32) + un[j + n - 2]) { + qd = qd - 1; + rhat = rhat + vn[n - 1]; + } + } + + /* Multiply and subtract */ + k = 0; + for (i = 0; i < n; i++) { + uint64_t p = (uint64_t)qd * vn[i]; + + t = un[i + j] - k - (p & 0xFFFFFFFF); + un[i + j] = t; + k = (p >> 32) - (t >> 32); + } + t = un[j + n] - k; + un[j + n] = t; + + /* If borrowed, add one back and adjust estimate */ + if (t < 0) { + k = 0; + qd = qd - 1; + for (i = 0; i < n; i++) { + t = (uint64_t)un[i + j] + vn[i] + k; + un[i + j] = t; + k = t >> 32; + } + un[j + n] = un[j + n] + k; + } + + BN_DIGIT(q, j) = qd; + } + + if (r != NULL) { + /* Denormalize un into r */ + if (s != 0) { + for (i = 0; i < n - 1; i++) + BN_DIGIT(r, i) = (un[i] >> s) | + (un[i + 1] << (32 - s)); + BN_DIGIT(r, n - 1) = un[n - 1] >> s; + } else { + for (i = 0; i < n; i++) + BN_DIGIT(r, i) = un[i]; + } + } + + return 1; +} + +static void bn_set_bn(struct LITE_BIGNUM *d, const struct LITE_BIGNUM *src, + size_t n) +{ + size_t i = 0; + + for (; i < n && i < d->dmax; ++i) + BN_DIGIT(d, i) = BN_DIGIT(src, i); + for (; i < d->dmax; ++i) + BN_DIGIT(d, i) = 0; +} + +static size_t bn_digits(const struct LITE_BIGNUM *a) +{ + size_t n = a->dmax - 1; + + while (BN_DIGIT(a, n) == 0 && n) + --n; + return n + 1; +} + +int DCRYPTO_bn_div(struct LITE_BIGNUM *quotient, + struct LITE_BIGNUM *remainder, + const struct LITE_BIGNUM *src, + const struct LITE_BIGNUM *divisor) +{ + int src_len = bn_digits(src); + int div_len = bn_digits(divisor); + int i, result; + + if (src_len < div_len) + return 0; + + result = bn_div_ex(quotient, remainder, + src, src_len, + divisor, div_len); + + if (!result) + return 0; + + /* 0-pad the destinations. */ + for (i = src_len - div_len + 1; i < quotient->dmax; ++i) + BN_DIGIT(quotient, i) = 0; + if (remainder) { + for (i = div_len; i < remainder->dmax; ++i) + BN_DIGIT(remainder, i) = 0; + } + + return result; +} + +/* + * Extended Euclid modular inverse. + * + * https://en.wikipedia.org/wiki/Extended_Euclidean_algorithm + * #Computing_multiplicative_inverses_in_modular_structures: + + * function inverse(a, n) + * t := 0; newt := 1; + * r := n; newr := a; + * while newr ≠ 0 + * quotient := r div newr + * (t, newt) := (newt, t - quotient * newt) + * (r, newr) := (newr, r - quotient * newr) + * if r > 1 then return "a is not invertible" + * if t < 0 then t := t + n + * return t + */ +int bn_modinv_vartime(struct LITE_BIGNUM *dst, const struct LITE_BIGNUM *src, + const struct LITE_BIGNUM *mod) +{ + uint32_t R_buf[RSA_MAX_WORDS]; + uint32_t nR_buf[RSA_MAX_WORDS]; + uint32_t Q_buf[RSA_MAX_WORDS]; + + uint32_t nT_buf[RSA_MAX_WORDS + 1]; /* Can go negative, hence +1 */ + uint32_t T_buf[RSA_MAX_WORDS + 1]; /* Can go negative */ + uint32_t tmp_buf[2 * RSA_MAX_WORDS + 1]; /* needs to hold Q*nT */ + + struct LITE_BIGNUM R; + struct LITE_BIGNUM nR; + struct LITE_BIGNUM Q; + struct LITE_BIGNUM T; + struct LITE_BIGNUM nT; + struct LITE_BIGNUM tmp; + + struct LITE_BIGNUM *pT = &T; + struct LITE_BIGNUM *pnT = &nT; + struct LITE_BIGNUM *pR = &R; + struct LITE_BIGNUM *pnR = &nR; + struct LITE_BIGNUM *bnswap; + + int t_neg = 0; + int nt_neg = 0; + int iswap; + + size_t r_len, nr_len; + + bn_init(&R, R_buf, bn_size(mod)); + bn_init(&nR, nR_buf, bn_size(mod)); + bn_init(&Q, Q_buf, bn_size(mod)); + bn_init(&T, T_buf, bn_size(mod) + sizeof(uint32_t)); + bn_init(&nT, nT_buf, bn_size(mod) + sizeof(uint32_t)); + bn_init(&tmp, tmp_buf, bn_size(mod) + sizeof(uint32_t)); + + r_len = bn_digits(mod); + nr_len = bn_digits(src); + + BN_DIGIT(&nT, 0) = 1; /* T = 0, nT = 1 */ + bn_set_bn(&R, mod, r_len); /* R = n */ + bn_set_bn(&nR, src, nr_len); /* nR = input */ + + /* Trim nR */ + while (nr_len && BN_DIGIT(&nR, nr_len - 1) == 0) + --nr_len; + + while (nr_len) { + size_t q_len = r_len - nr_len + 1; + + /* (r, nr) = (nr, r % nr), q = r / nr */ + if (!bn_div_ex(&Q, pR, pR, r_len, pnR, nr_len)) + return 0; + + /* swap R and nR */ + r_len = nr_len; + bnswap = pR; pR = pnR; pnR = bnswap; + + /* trim nR and Q */ + while (nr_len && BN_DIGIT(pnR, nr_len - 1) == 0) + --nr_len; + while (q_len && BN_DIGIT(&Q, q_len - 1) == 0) + --q_len; + + Q.dmax = q_len; + + /* compute t - q*nt */ + if (q_len == 1 && BN_DIGIT(&Q, 0) <= 2) { + /* Doing few direct subs is faster than mul + sub */ + uint32_t n = BN_DIGIT(&Q, 0); + + while (n--) + bn_signed_sub(pT, &t_neg, pnT, nt_neg); + } else { + /* Call bn_mul_ex with smallest operand first */ + if (nt_neg) { + /* Negative numbers use all digits, + * thus pnT is large + */ + bn_mul_ex(&tmp, &Q, q_len, pnT); + } else { + int nt_len = bn_digits(pnT); + + if (q_len < nt_len) + bn_mul_ex(&tmp, &Q, q_len, pnT); + else + bn_mul_ex(&tmp, pnT, nt_len, &Q); + } + bn_signed_sub(pT, &t_neg, &tmp, nt_neg); + } + + /* swap T and nT */ + bnswap = pT; pT = pnT; pnT = bnswap; + iswap = t_neg; t_neg = nt_neg; nt_neg = iswap; + } + + if (r_len != 1 || BN_DIGIT(pR, 0) != 1) { + /* gcd not 1; no direct inverse */ + return 0; + } + + if (t_neg) + bn_signed_add(pT, &t_neg, mod, 0); + + bn_set_bn(dst, pT, bn_digits(pT)); + + return 1; +} + +#define PRIME1 3 + +/* + * The array below is an encoding of the first 4096 primes, starting with + * PRIME1. Using 4096 of the first primes results in at least 5% improvement + * in running time over using the first 2048. + * + * Most byte entries in the array contain two sequential differentials between + * two adjacent prime numbers, each differential halved (as the difference is + * always even) and packed into 4 bits. + * + * If a halved differential value exceeds 0xf (and as such does not fit into 4 + * bits), a zero is placed in the array followed by the value literal (no + * halving). + * + * If out of two consecutive differencials only the second one exceeds 0xf, + * the first one still is put into the array in its own byte prepended by a + * zero. + */ +const uint8_t PRIME_DELTAS[] = { + 1, 18, 18, 18, 49, 50, 18, 51, 19, 33, 50, 52, + 33, 33, 39, 35, 21, 19, 50, 51, 21, 18, 22, 98, + 18, 49, 83, 51, 19, 33, 87, 33, 39, 53, 18, 52, + 51, 35, 66, 69, 21, 19, 35, 66, 18, 100, 36, 35, + 97, 147, 83, 49, 53, 51, 19, 50, 22, 81, 35, 49, + 98, 52, 84, 84, 51, 36, 50, 66, 117, 97, 81, 33, + 87, 33, 39, 33, 42, 36, 84, 35, 55, 35, 52, 54, + 35, 21, 19, 81, 81, 57, 33, 35, 52, 51, 177, 84, + 83, 52, 98, 51, 19, 101, 145, 35, 19, 33, 38, 19, + 0, 34, 51, 73, 87, 33, 35, 66, 19, 101, 18, 18, + 54, 100, 99, 35, 66, 66, 114, 49, 35, 19, 90, 50, + 28, 33, 86, 21, 67, 51, 147, 33, 101, 100, 135, 50, + 18, 21, 99, 57, 24, 27, 52, 50, 18, 67, 81, 87, + 83, 97, 33, 86, 24, 19, 33, 84, 156, 35, 72, 18, + 72, 18, 67, 50, 97, 179, 19, 35, 115, 33, 50, 54, + 51, 114, 54, 67, 45, 149, 66, 49, 59, 97, 132, 38, + 117, 18, 67, 50, 18, 52, 33, 53, 21, 66, 117, 97, + 50, 24, 114, 52, 50, 148, 83, 52, 86, 114, 51, 30, + 21, 66, 114, 70, 54, 35, 165, 24, 210, 22, 50, 99, + 66, 75, 18, 22, 225, 51, 50, 49, 98, 97, 81, 129, + 131, 168, 66, 18, 27, 70, 53, 18, 49, 53, 22, 81, + 87, 50, 52, 51, 134, 18, 115, 36, 84, 51, 179, 21, + 114, 57, 21, 114, 21, 114, 73, 35, 18, 49, 98, 171, + 97, 35, 49, 59, 19, 131, 97, 54, 129, 35, 114, 25, + 197, 49, 81, 81, 83, 21, 21, 52, 245, 21, 67, 89, + 54, 97, 147, 35, 57, 21, 115, 33, 44, 22, 56, 67, + 57, 129, 35, 19, 53, 54, 105, 19, 41, 76, 33, 35, + 22, 39, 245, 54, 115, 86, 18, 52, 53, 18, 115, 50, + 49, 81, 134, 73, 35, 97, 51, 62, 55, 36, 84, 105, + 33, 44, 99, 24, 51, 117, 114, 243, 51, 67, 33, 99, + 33, 59, 49, 41, 18, 97, 50, 211, 50, 69, 0, 32, + 129, 50, 18, 21, 115, 36, 83, 162, 19, 242, 69, 51, + 67, 98, 49, 50, 49, 81, 131, 162, 103, 227, 162, 148, + 50, 55, 51, 81, 86, 69, 21, 70, 92, 18, 67, 36, + 149, 51, 19, 86, 21, 51, 52, 53, 49, 51, 53, 76, + 59, 25, 36, 95, 73, 33, 83, 19, 41, 70, 152, 49, + 99, 81, 81, 53, 114, 193, 129, 81, 90, 33, 36, 131, + 49, 104, 66, 63, 21, 19, 35, 52, 50, 99, 70, 39, + 101, 195, 99, 27, 73, 83, 114, 19, 84, 50, 63, 117, + 22, 81, 129, 156, 147, 137, 49, 146, 49, 84, 83, 52, + 35, 21, 22, 35, 49, 98, 121, 35, 162, 67, 36, 39, + 50, 118, 33, 242, 195, 54, 103, 50, 18, 147, 100, 50, + 97, 111, 129, 59, 115, 86, 49, 36, 83, 60, 115, 36, + 105, 81, 81, 35, 163, 39, 33, 39, 54, 197, 52, 81, + 242, 49, 98, 115, 0, 34, 100, 53, 18, 165, 72, 21, + 114, 22, 56, 52, 36, 35, 67, 54, 50, 51, 73, 42, + 38, 21, 49, 86, 18, 163, 243, 36, 86, 49, 225, 50, + 24, 97, 53, 76, 99, 147, 39, 50, 100, 54, 35, 99, + 97, 138, 33, 89, 66, 114, 19, 179, 115, 53, 49, 81, + 33, 177, 35, 54, 55, 86, 52, 0, 4, 0, 36, 118, + 50, 49, 99, 104, 21, 75, 22, 50, 57, 22, 50, 100, + 54, 35, 99, 22, 98, 115, 131, 21, 73, 0, 6, 0, + 34, 30, 27, 49, 86, 19, 36, 179, 21, 66, 52, 38, + 150, 162, 51, 66, 24, 97, 84, 81, 35, 118, 180, 225, + 42, 33, 39, 86, 22, 129, 228, 180, 35, 55, 36, 99, + 50, 162, 145, 99, 35, 121, 84, 0, 10, 0, 32, 53, + 51, 19, 131, 22, 62, 21, 72, 52, 53, 202, 81, 81, + 98, 58, 33, 105, 81, 81, 42, 141, 36, 50, 99, 70, + 99, 36, 177, 135, 83, 102, 115, 42, 38, 49, 51, 132, + 177, 228, 50, 162, 108, 162, 69, 24, 22, 0, 12, 0, + 34, 18, 54, 51, 67, 33, 60, 42, 83, 55, 35, 49, + 99, 81, 83, 162, 210, 19, 177, 194, 49, 35, 195, 66, + 0, 2, 0, 34, 52, 134, 21, 21, 52, 36, 107, 55, + 45, 33, 101, 66, 70, 39, 56, 52, 35, 52, 53, 97, + 51, 132, 51, 101, 19, 146, 51, 54, 148, 53, 73, 39, + 57, 84, 86, 19, 102, 0, 36, 35, 66, 49, 41, 99, + 67, 50, 145, 33, 194, 51, 127, 50, 54, 58, 36, 36, + 51, 47, 21, 100, 84, 195, 98, 114, 49, 231, 129, 99, + 42, 83, 51, 69, 103, 87, 135, 87, 56, 52, 56, 165, + 19, 33, 38, 21, 19, 179, 18, 148, 84, 177, 89, 114, + 18, 145, 35, 69, 31, 47, 21, 25, 41, 55, 81, 42, + 0, 36, 50, 55, 42, 87, 179, 31, 101, 145, 39, 59, + 145, 99, 36, 36, 53, 22, 149, 120, 114, 51, 19, 33, + 225, 227, 18, 55, 38, 120, 114, 52, 50, 51, 52, 36, + 39, 132, 50, 100, 129, 84, 35, 211, 84, 35, 103, 242, + 123, 70, 35, 69, 55, 83, 21, 102, 115, 57, 83, 73, + 35, 19, 81, 84, 51, 81, 149, 22, 35, 69, 103, 98, + 69, 51, 162, 120, 117, 69, 97, 147, 101, 97, 33, 99, + 36, 0, 4, 0, 44, 33, 33, 86, 51, 114, 51, 52, + 0, 6, 0, 36, 146, 49, 99, 51, 39, 182, 25, 83, + 220, 33, 33, 39, 35, 52, 134, 0, 2, 0, 42, 33, + 44, 51, 25, 39, 62, 151, 53, 97, 54, 243, 35, 55, + 33, 194, 51, 213, 147, 67, 63, 38, 97, 129, 50, 105, + 19, 45, 99, 98, 204, 99, 22, 228, 35, 97, 147, 35, + 58, 129, 51, 149, 49, 36, 51, 200, 52, 83, 123, 72, + 49, 98, 27, 73, 0, 34, 19, 146, 51, 69, 73, 50, + 18, 72, 22, 99, 146, 51, 49, 54, 90, 105, 35, 24, + 21, 114, 241, 86, 28, 56, 69, 22, 179, 24, 165, 22, + 105, 86, 49, 81, 53, 145, 99, 35, 28, 225, 33, 81, + 134, 75, 19, 33, 83, 166, 84, 99, 51, 41, 18, 105, + 22, 50, 24, 102, 114, 73, 38, 115, 50, 67, 42, 101, + 114, 24, 22, 242, 60, 172, 84, 101, 99, 102, 52, 135, + 50, 0, 6, 0, 36, 165, 246, 18, 30, 103, 59, 66, + 147, 121, 35, 19, 0, 34, 145, 131, 145, 194, 19, 99, + 101, 67, 134, 69, 0, 14, 0, 40, 49, 50, 103, 33, + 33, 36, 53, 51, 19, 51, 99, 197, 21, 54, 51, 115, + 0, 6, 0, 52, 163, 81, 84, 86, 97, 50, 120, 70, + 59, 21, 67, 177, 179, 69, 102, 21, 54, 18, 117, 19, + 146, 100, 150, 51, 35, 55, 33, 102, 35, 153, 97, 134, + 73, 93, 35, 67, 50, 21, 162, 52, 42, 81, 0, 34, + 18, 193, 102, 83, 22, 243, 104, 97, 185, 103, 81, 102, + 33, 35, 97, 137, 0, 2, 0, 40, 72, 52, 81, 41, + 69, 70, 41, 25, 81, 33, 36, 225, 59, 99, 121, 35, + 67, 53, 66, 25, 83, 171, 67, 242, 18, 147, 241, 36, + 50, 54, 0, 14, 0, 34, 115, 33, 50, 114, 19, 225, + 35, 69, 21, 21, 18, 241, 102, 89, 103, 81, 99, 83, + 118, 39, 41, 21, 66, 69, 105, 148, 57, 135, 51, 87, + 35, 22, 98, 51, 97, 129, 99, 39, 50, 22, 146, 0, + 36, 150, 97, 33, 36, 98, 0, 36, 57, 22, 83, 108, + 67, 56, 97, 149, 165, 19, 146, 0, 2, 0, 40, 49, + 129, 36, 149, 99, 21, 66, 54, 21, 148, 50, 162, 0, + 6, 0, 36, 49, 83, 195, 120, 57, 21, 165, 67, 35, + 21, 22, 33, 36, 83, 105, 118, 132, 56, 66, 19, 156, + 149, 97, 39, 83, 51, 150, 30, 151, 134, 124, 107, 49, + 84, 33, 39, 99, 35, 114, 18, 243, 19, 81, 251, 18, + 52, 51, 134, 99, 66, 28, 98, 52, 51, 81, 54, 231, + 50, 100, 54, 35, 115, 101, 51, 67, 50, 18, 70, 39, + 149, 24, 58, 53, 66, 0, 30, 0, 36, 100, 182, 19, + 104, 51, 25, 45, 36, 149, 69, 55, 42, 185, 100, 230, + 51, 67, 108, 135, 39, 99, 86, 163, 36, 150, 149, 18, + 165, 114, 49, 92, 145, 42, 135, 87, 50, 58, 53, 49, + 99, 245, 67, 35, 0, 8, 0, 40, 18, 22, 146, 52, + 83, 153, 22, 132, 50, 51, 0, 2, 0, 52, 114, 168, + 18, 54, 19, 102, 50, 117, 51, 117, 120, 67, 98, 75, + 49, 155, 49, 147, 135, 83, 97, 50, 73, 104, 18, 114, + 70, 111, 132, 33, 59, 100, 83, 51, 115, 149, 97, 81, + 45, 38, 66, 148, 87, 131, 52, 83, 67, 101, 165, 66, + 109, 146, 105, 63, 52, 59, 97, 35, 49, 81, 35, 49, + 59, 147, 150, 70, 53, 97, 129, 81, 89, 58, 33, 59, + 51, 147, 118, 129, 51, 39, 98, 25, 0, 16, 0, 36, + 99, 126, 22, 54, 50, 24, 244, 195, 245, 25, 35, 100, + 177, 59, 145, 81, 95, 30, 55, 131, 168, 19, 0, 4, + 0, 32, 33, 35, 22, 35, 54, 19, 35, 67, 42, 0, + 4, 0, 32, 84, 129, 177, 35, 67, 135, 41, 66, 163, + 102, 53, 21, 22, 230, 145, 149, 69, 0, 48, 18, 52, + 81, 95, 0, 2, 0, 36, 53, 49, 146, 52, 135, 131, + 114, 162, 49, 86, 19, 99, 50, 97, 50, 99, 66, 19, + 149, 52, 99, 177, 54, 146, 115, 42, 56, 66, 75, 70, + 51, 134, 159, 66, 18, 61, 39, 203, 49, 53, 55, 51, + 101, 49, 101, 100, 153, 83, 72, 51, 72, 162, 21, 21, + 99, 67, 90, 89, 210, 63, 18, 67, 102, 146, 75, 49, + 0, 12, 0, 34, 57, 99, 30, 120, 114, 118, 35, 49, + 0, 36, 35, 166, 195, 177, 137, 102, 145, 51, 50, 55, + 33, 180, 99, 83, 70, 150, 53, 27, 115, 50, 147, 171, + 22, 194, 153, 27, 18, 100, 101, 114, 25, 0, 16, 0, + 38, 51, 54, 83, 100, 50, 55, 243, 84, 179, 70, 81, + 81, 53, 21, 105, 163, 36, 179, 63, 55, 54, 99, 81, + 95, 24, 66, 19, 146, 19, 45, 36, 53, 18, 52, 35, + 246, 19, 50, 171, 66, 18, 0, 72, 66, 75, 18, 117, + 18, 163, 89, 58, 131, 67, 42, 107, 18, 22, 89, 27, + 57, 241, 87, 84, 0, 16, 0, 50, 53, 69, 99, 145, + 179, 18, 52, 51, 89, 27, 24, 117, 49, 101, 162, 115, + 0, 4, 0, 36, 18, 54, 18, 118, 50, 49, 50, 165, + 21, 54, 28, 102, 51, 44, 18, 193, 50, 52, 131, 21, + 103, 0, 6, 0, 34, 55, 50, 31, 180, 35, 66, 30, + 19, 45, 155, 19, 131, 24, 97, 98, 51, 117, 52, 98, + 145, 84, 131, 63, 21, 145, 84, 36, 108, 0, 40, 22, + 83, 97, 98, 18, 57, 118, 50, 127, 36, 84, 53, 148, + 39, 131, 66, 49, 81, 98, 18, 52, 35, 0, 32, 197, + 73, 81, 53, 18, 147, 97, 129, 179, 52, 146, 150, 67, + 42, 63, 182, 19, 146, 0, 62, 33, 99, 81, 102, 225, + 39, 179, 19, 53, 114, 21, 52, 87, 83, 22, 185, 69, + 150, 22, 38, 21, 19, 147, 0, 6, 0, 34, 49, 98, + 57, 145, 131, 52, 53, 148, 84, 81, 41, 214, 177, 33, + 179, 55, 131, 165, 97, 0, 18, 0, 42, 44, 19, 86, + 19, 84, 35, 102, 66, 54, 250, 60, 53, 97, 90, 51, + 38, 117, 150, 67, 98, 117, 22, 248, 22, 50, 18, 61, + 41, 18, 55, 0, 54, 0, 6, 0, 52, 24, 51, 109, + 33, 59, 49, 102, 53, 145, 102, 89, 99, 67, 83, 66, + 18, 172, 51, 87, 81, 179, 117, 210, 148, 102, 86, 52, + 131, 67, 59, 21, 165, 0, 6, 0, 44, 147, 81, 35, + 114, 210, 22, 84, 36, 98, 100, 180, 53, 147, 52, 54, + 36, 149, 99, 97, 50, 24, 102, 117, 115, 86, 22, 50, + 49, 98, 211, 147, 83, 25, 84, 45, 90, 56, 166, 84, + 81, 131, 165, 162, 241, 36, 129, 146, 19, 89, 103, 147, + 138, 50, 67, 35, 100, 81, 99, 33, 53, 24, 103, 83, + 67, 225, 57, 0, 30, 0, 34, 24, 97, 152, 52, 84, + 84, 0, 10, 0, 44, 51, 42, 33, 39, 228, 56, 127, + 63, 39, 83, 52, 41, 99, 27, 100, 54, 39, 35, 18, + 154, 56, 0, 38, 129, 35, 0, 2, 0, 40, 0, 42, + 114, 49, 197, 49, 149, 97, 129, 56, 52, 33, 83, 69, + 25, 132, 105, 99, 101, 51, +}; + +static uint32_t bn_mod_word16(const struct LITE_BIGNUM *p, uint16_t word) +{ + int i; + uint32_t rem = 0; + + for (i = p->dmax - 1; i >= 0; i--) { + rem = ((rem << 16) | + ((BN_DIGIT(p, i) >> 16) & 0xFFFFUL)) % word; + rem = ((rem << 16) | (BN_DIGIT(p, i) & 0xFFFFUL)) % word; + } + + return rem; +} + +static uint32_t bn_mod_f4(const struct LITE_BIGNUM *d) +{ + int i = bn_size(d) - 1; + const uint8_t *p = (const uint8_t *) (d->d); + uint32_t rem = 0; + + for (; i >= 0; --i) { + uint32_t q = RSA_F4 * (rem >> 8); + + if (rem < q) + q -= RSA_F4; + rem <<= 8; + rem |= p[i]; + rem -= q; + } + + if (rem >= RSA_F4) + rem -= RSA_F4; + + return rem; +} + +#define bn_is_even(b) !bn_is_bit_set((b), 0) +/* From HAC Fact 4.48 (ii), the following number of + * rounds suffice for ~2^145 confidence. Each additional + * round provides about another k/100 bits of confidence. */ +#define ROUNDS_1024 7 +#define ROUNDS_512 15 +#define ROUNDS_384 22 + +/* Miller-Rabin from HAC, algorithm 4.24. */ +static int bn_probable_prime(const struct LITE_BIGNUM *p) +{ + int j; + int s = 0; + + uint32_t ONE_buf = 1; + uint8_t r_buf[RSA_MAX_BYTES / 2]; + uint8_t A_buf[RSA_MAX_BYTES / 2]; + uint8_t y_buf[RSA_MAX_BYTES / 2]; + + struct LITE_BIGNUM ONE; + struct LITE_BIGNUM r; + struct LITE_BIGNUM A; + struct LITE_BIGNUM y; + + const int rounds = bn_bits(p) >= 1024 ? ROUNDS_1024 : + bn_bits(p) >= 512 ? ROUNDS_512 : + ROUNDS_384; + + /* Failsafe: update rounds table above to support smaller primes. */ + if (bn_bits(p) < 384) + return 0; + + if (bn_size(p) > sizeof(r_buf)) + return 0; + + DCRYPTO_bn_wrap(&ONE, &ONE_buf, sizeof(ONE_buf)); + DCRYPTO_bn_wrap(&r, r_buf, bn_size(p)); + bn_copy(&r, p); + + /* r * (2 ^ s) = p - 1 */ + bn_sub(&r, &ONE); + while (bn_is_even(&r)) { + bn_rshift(&r, 0, 0); + s++; + } + + DCRYPTO_bn_wrap(&A, A_buf, bn_size(p)); + DCRYPTO_bn_wrap(&y, y_buf, bn_size(p)); + for (j = 0; j < rounds; j++) { + int i; + + /* pick random A, such that A < p */ + rand_bytes(A_buf, bn_size(&A)); + for (i = A.dmax - 1; i >= 0; i--) { + while (BN_DIGIT(&A, i) > BN_DIGIT(p, i)) + BN_DIGIT(&A, i) = rand(); + if (BN_DIGIT(&A, i) < BN_DIGIT(p, i)) + break; + } + + /* y = a ^ r mod p */ + bn_modexp(&y, &A, &r, p); + if (bn_eq(&y, &ONE)) + continue; + bn_add(&y, &ONE); + if (bn_eq(&y, p)) + continue; + bn_sub(&y, &ONE); + + /* y = y ^ 2 mod p */ + for (i = 0; i < s - 1; i++) { + bn_copy(&A, &y); + bn_modexp_word(&y, &A, 2, p); + + if (bn_eq(&y, &ONE)) + return 0; + + bn_add(&y, &ONE); + if (bn_eq(&y, p)) { + bn_sub(&y, &ONE); + break; + } + bn_sub(&y, &ONE); + } + bn_add(&y, &ONE); + if (!bn_eq(&y, p)) + return 0; + } + + return 1; +} + +/* #define PRINT_PRIMES to enable printing predefined prime numbers' set. */ +static void print_primes(uint16_t prime) +{ +#ifdef PRINT_PRIMES + static uint16_t num_per_line; + static uint16_t max_printed; + + if (prime <= max_printed) + return; + + if (!(num_per_line++ % 8)) { + if (num_per_line == 1) + ccprintf("Prime numbers:"); + ccprintf("\n"); + cflush(); + } + max_printed = prime; + ccprintf(" %6d", prime); +#endif +} + +int DCRYPTO_bn_generate_prime(struct LITE_BIGNUM *p) +{ + int i; + int j; + /* Using a sieve size of 2048-bits results in a failure rate + * of ~0.5% @ 1024-bit candidates. The failure rate rises to ~6% + * if the sieve size is halved. */ + uint8_t composites_buf[256]; + struct LITE_BIGNUM composites; + uint16_t prime = PRIME1; + + /* Set top two bits, as well as LSB. */ + bn_set_bit(p, 0); + bn_set_bit(p, bn_bits(p) - 1); + bn_set_bit(p, bn_bits(p) - 2); + + /* Save on trial division by marking known composites. */ + bn_init(&composites, composites_buf, sizeof(composites_buf)); + for (i = 0; i < ARRAY_SIZE(PRIME_DELTAS); i++) { + uint16_t rem; + uint8_t unpacked_deltas[2]; + uint8_t packed_deltas = PRIME_DELTAS[i]; + int k; + int m; + + if (packed_deltas) { + unpacked_deltas[0] = (packed_deltas >> 4) << 1; + unpacked_deltas[1] = (packed_deltas & 0xf) << 1; + m = 2; + } else { + i += 1; + unpacked_deltas[0] = PRIME_DELTAS[i]; + m = 1; + } + + for (k = 0; k < m; k++) { + prime += unpacked_deltas[k]; + print_primes(prime); + rem = bn_mod_word16(p, prime); + /* Skip marking odd offsets (i.e. even candidates). */ + for (j = (rem == 0) ? 0 : prime - rem; + j < bn_bits(&composites) << 1; + j += prime) { + if ((j & 1) == 0) + bn_set_bit(&composites, j >> 1); + } + } + } + + /* composites now marked, apply Miller-Rabin to prime candidates. */ + j = 0; + for (i = 0; i < bn_bits(&composites); i++) { + uint32_t diff_buf; + struct LITE_BIGNUM diff; + + if (bn_is_bit_set(&composites, i)) + continue; + + /* Recover increment from the composites sieve. */ + diff_buf = (i << 1) - j; + j = (i << 1); + DCRYPTO_bn_wrap(&diff, &diff_buf, sizeof(diff_buf)); + bn_add(p, &diff); + /* Make sure prime will work with F4 public exponent. */ + if (bn_mod_f4(p) >= 2) { + if (bn_probable_prime(p)) + return 1; + } + } + + always_memset(composites_buf, 0, sizeof(composites_buf)); + return 0; +} |