1*8ccd4a63SDavid du Colombier #include "os.h"
2*8ccd4a63SDavid du Colombier #include <mp.h>
3*8ccd4a63SDavid du Colombier #include <libsec.h>
4*8ccd4a63SDavid du Colombier
5*8ccd4a63SDavid du Colombier // Gordon's algorithm for generating a strong prime
6*8ccd4a63SDavid du Colombier // Menezes et al () Handbook, p.150
7*8ccd4a63SDavid du Colombier void
genstrongprime(mpint * p,int n,int accuracy)8*8ccd4a63SDavid du Colombier genstrongprime(mpint *p, int n, int accuracy)
9*8ccd4a63SDavid du Colombier {
10*8ccd4a63SDavid du Colombier mpint *s, *t, *r, *i;
11*8ccd4a63SDavid du Colombier
12*8ccd4a63SDavid du Colombier if(n < 64)
13*8ccd4a63SDavid du Colombier n = 64;
14*8ccd4a63SDavid du Colombier
15*8ccd4a63SDavid du Colombier s = mpnew(n/2);
16*8ccd4a63SDavid du Colombier genprime(s, (n/2)-16, accuracy);
17*8ccd4a63SDavid du Colombier t = mpnew(n/2);
18*8ccd4a63SDavid du Colombier genprime(t, n-mpsignif(s)-32, accuracy);
19*8ccd4a63SDavid du Colombier
20*8ccd4a63SDavid du Colombier // first r = 2it + 1 that's prime
21*8ccd4a63SDavid du Colombier i = mpnew(16);
22*8ccd4a63SDavid du Colombier r = mpnew(0);
23*8ccd4a63SDavid du Colombier itomp(0x8000, i);
24*8ccd4a63SDavid du Colombier mpleft(t, 1, t); // 2t
25*8ccd4a63SDavid du Colombier mpmul(i, t, r); // 2it
26*8ccd4a63SDavid du Colombier mpadd(r, mpone, r); // 2it + 1
27*8ccd4a63SDavid du Colombier for(;;){
28*8ccd4a63SDavid du Colombier if(probably_prime(r, 18))
29*8ccd4a63SDavid du Colombier break;
30*8ccd4a63SDavid du Colombier mpadd(r, t, r); // r += 2t
31*8ccd4a63SDavid du Colombier }
32*8ccd4a63SDavid du Colombier
33*8ccd4a63SDavid du Colombier // p0 = 2(s**(r-2) mod r)s - 1
34*8ccd4a63SDavid du Colombier itomp(2, p);
35*8ccd4a63SDavid du Colombier mpsub(r, p, p);
36*8ccd4a63SDavid du Colombier mpexp(s, p, r, p);
37*8ccd4a63SDavid du Colombier mpmul(s, p, p);
38*8ccd4a63SDavid du Colombier mpleft(p, 1, p);
39*8ccd4a63SDavid du Colombier mpsub(p, mpone, p);
40*8ccd4a63SDavid du Colombier
41*8ccd4a63SDavid du Colombier // first p = p0 + 2irs that's prime
42*8ccd4a63SDavid du Colombier itomp(0x8000, i);
43*8ccd4a63SDavid du Colombier mpleft(r, 1, r); // 2r
44*8ccd4a63SDavid du Colombier mpmul(r, s, r); // 2rs
45*8ccd4a63SDavid du Colombier mpmul(r, i, i); // 2irs
46*8ccd4a63SDavid du Colombier mpadd(p, i, p); // p0 + 2irs
47*8ccd4a63SDavid du Colombier for(;;){
48*8ccd4a63SDavid du Colombier if(probably_prime(p, accuracy))
49*8ccd4a63SDavid du Colombier break;
50*8ccd4a63SDavid du Colombier mpadd(p, r, p); // p += 2rs
51*8ccd4a63SDavid du Colombier }
52*8ccd4a63SDavid du Colombier
53*8ccd4a63SDavid du Colombier mpfree(i);
54*8ccd4a63SDavid du Colombier mpfree(s);
55*8ccd4a63SDavid du Colombier mpfree(r);
56*8ccd4a63SDavid du Colombier mpfree(t);
57*8ccd4a63SDavid du Colombier }
58