{
mp *s, *t, *q;
dstr dn = DSTR_INIT;
+ unsigned slop, nb, u, i;
mp *rr = d;
pgen_filterctx c;
pgen_jumpctx j;
rabin rb;
- /* --- The bitslop parameter --- *
+ /* --- Figure out how large the smaller primes should be --- *
*
- * There's quite a lot of prime searching to be done. The constant
- * @BITSLOP@ is a (low) approximation to the base-2 log of the expected
- * number of steps to find a prime number. Experimentation shows that
- * numbers around 10 seem to be good.
+ * We want them to be `as large as possible', subject to the constraint
+ * that we produce a number of the requested size at the end. This is
+ * tricky, because the final prime search is going to involve quite large
+ * jumps from its starting point; the size of the jumps are basically
+ * determined by our choice here, and if they're too big then we won't find
+ * a prime in time.
+ *
+ * Let's suppose we're trying to make an %$N$%-bit prime. The expected
+ * number of steps tends to increase linearly with size, i.e., we need to
+ * take about %2^k N$% steps for some %$k$%. If we're jumping by a
+ * %$J$%-bit quantity each time, from an %$N$%-bit starting point, then we
+ * will only be able to find a match if %$2^k N 2^{J-1} \le 2^{N-1}$%,
+ * i.e., if %$J \le N - (k + \log_2 N)$%.
+ *
+ * Experimentation shows that taking %$k + \log_2 N = 12$% works well for
+ * %$N = 1024$%, so %$k = 2$%.
*/
-#define BITSLOP 12
+ for (i = 1; i && nbits >> i; i <<= 1); assert(i);
+ for (slop = 2, nb = nbits; nb > 1; i >>= 1) {
+ u = nb >> i;
+ if (u) { slop += i; nb = u; }
+ }
+ if (nbits/2 <= slop) return (0);
/* --- Choose two primes %$s$% and %$t$% of half the required size --- */
- assert(((void)"nbits too small in strongprime_setup", nbits/2 > BITSLOP));
- nbits = nbits/2 - BITSLOP;
+ nb = nbits/2 - slop;
c.step = 1;
- rr = mprand(rr, nbits, r, 1);
+ rr = mprand(rr, nb, r, 1);
DRESET(&dn); dstr_putf(&dn, "%s [s]", name);
if ((s = pgen(dn.buf, MP_NEWSEC, rr, event, ectx, n, pgen_filter, &c,
- rabin_iters(nbits), pgen_test, &rb)) == 0)
+ rabin_iters(nb), pgen_test, &rb)) == 0)
goto fail_s;
- rr = mprand(rr, nbits, r, 1);
+ rr = mprand(rr, nb, r, 1);
DRESET(&dn); dstr_putf(&dn, "%s [t]", name);
if ((t = pgen(dn.buf, MP_NEWSEC, rr, event, ectx, n, pgen_filter, &c,
- rabin_iters(nbits), pgen_test, &rb)) == 0)
+ rabin_iters(nb), pgen_test, &rb)) == 0)
goto fail_t;
/* --- Choose a suitable value for %$r = 2it + 1$% for some %$i$% --- */
rr = mp_lsl(rr, t, 1);
pfilt_create(&c.f, rr);
- rr = mp_lsl(rr, rr, BITSLOP - 1);
+ rr = mp_lsl(rr, rr, slop - 1);
rr = mp_add(rr, rr, MP_ONE);
DRESET(&dn); dstr_putf(&dn, "%s [r]", name);
j.j = &c.f;
- nbits += BITSLOP;
q = pgen(dn.buf, MP_NEW, rr, event, ectx, n, pgen_jump, &j,
- rabin_iters(nbits), pgen_test, &rb);
+ rabin_iters(nb + slop), pgen_test, &rb);
pfilt_destroy(&c.f);
if (!q)
goto fail_r;
/* --- Now find %$p = p_0 + 2jrs$% for some %$j$% --- */
{
- mp *x;
+ mp *x, *y;
x = mp_mul(MP_NEW, q, s);
x = mp_lsl(x, x, 1);
pfilt_create(f, x);
- x = mp_lsl(x, x, BITSLOP - 1);
- rr = mp_add(rr, rr, x);
- mp_drop(x);
+ y = mp_lsl(MP_NEW, MP_ONE, nbits - 1);
+ rr = mp_leastcongruent(rr, y, rr, x);
+ mp_drop(x); mp_drop(y);
}
/* --- Return the result --- */
mp_drop(rr);
dstr_destroy(&dn);
return (0);
-
-#undef BITSLOP
}
/* --- @strongprime@ --- *
* * %$p - 1$% has a large prime factor %$r$%,
* * %$p + 1$% has a large prime factor %$s$%, and
* * %$r - 1$% has a large prime factor %$t$%.
- *
- * The numbers produced may be slightly larger than requested,
- * by a few bits.
*/
mp *strongprime(const char *name, mp *d, unsigned nbits, grand *r,
unsigned n, pgen_proc *event, void *ectx)
{
+ mp *p;
pfilt f;
pgen_jumpctx j;
rabin rb;
- d = strongprime_setup(name, d, &f, nbits, r, n, event, ectx);
+ if (d) mp_copy(d);
+ p = strongprime_setup(name, d, &f, nbits, r, n, event, ectx);
+ if (!p) { mp_drop(d); return (0); }
j.j = &f;
- d = pgen(name, d, d, event, ectx, n, pgen_jump, &j,
+ p = pgen(name, p, p, event, ectx, n, pgen_jump, &j,
rabin_iters(nbits), pgen_test, &rb);
+ if (mp_bits(p) != nbits) { mp_drop(p); return (0); }
pfilt_destroy(&f);
- return (d);
+ mp_drop(d);
+ return (p);
}
/*----- That's all, folks -------------------------------------------------*/