tests/gdsa: Test from P1363.
[u/mdw/catacomb] / strongprime.c
1 /* -*-c-*-
2 *
3 * $Id: strongprime.c,v 1.5 2004/04/08 01:36:15 mdw Exp $
4 *
5 * Generate `strong' prime numbers
6 *
7 * (c) 1999 Straylight/Edgeware
8 */
9
10 /*----- Licensing notice --------------------------------------------------*
11 *
12 * This file is part of Catacomb.
13 *
14 * Catacomb is free software; you can redistribute it and/or modify
15 * it under the terms of the GNU Library General Public License as
16 * published by the Free Software Foundation; either version 2 of the
17 * License, or (at your option) any later version.
18 *
19 * Catacomb is distributed in the hope that it will be useful,
20 * but WITHOUT ANY WARRANTY; without even the implied warranty of
21 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
22 * GNU Library General Public License for more details.
23 *
24 * You should have received a copy of the GNU Library General Public
25 * License along with Catacomb; if not, write to the Free
26 * Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
27 * MA 02111-1307, USA.
28 */
29
30 /*----- Header files ------------------------------------------------------*/
31
32 #include <mLib/dstr.h>
33
34 #include "grand.h"
35 #include "rand.h"
36 #include "mp.h"
37 #include "mpmont.h"
38 #include "mprand.h"
39 #include "pgen.h"
40 #include "pfilt.h"
41 #include "rabin.h"
42
43 /*----- Main code ---------------------------------------------------------*/
44
45 /* --- @strongprime_setup@ --- *
46 *
47 * Arguments: @const char *name@ = pointer to name root
48 * @mp *d@ = destination for search start point
49 * @pfilt *f@ = where to store filter jump context
50 * @unsigned nbits@ = number of bits wanted
51 * @grand *r@ = random number source
52 * @unsigned n@ = number of attempts to make
53 * @pgen_proc *event@ = event handler function
54 * @void *ectx@ = argument for the event handler
55 *
56 * Returns: A starting point for a `strong' prime search, or zero.
57 *
58 * Use: Sets up for a strong prime search, so that primes with
59 * particular properties can be found. It's probably important
60 * to note that the number left in the filter context @f@ is
61 * congruent to 2 (mod 4).
62 */
63
64 mp *strongprime_setup(const char *name, mp *d, pfilt *f, unsigned nbits,
65 grand *r, unsigned n, pgen_proc *event, void *ectx)
66 {
67 mp *s, *t, *q;
68 dstr dn = DSTR_INIT;
69
70 mp *rr = d;
71 pgen_filterctx c;
72 pgen_jumpctx j;
73 rabin rb;
74
75 /* --- The bitslop parameter --- *
76 *
77 * There's quite a lot of prime searching to be done. The constant
78 * @BITSLOP@ is a (low) approximation to the base-2 log of the expected
79 * number of steps to find a prime number. Experimentation shows that
80 * numbers around 10 seem to be good.
81 */
82
83 #define BITSLOP 12
84
85 /* --- Choose two primes %$s$% and %$t$% of half the required size --- */
86
87 assert(((void)"nbits too small in strongprime_setup", nbits/2 > BITSLOP));
88 nbits = nbits/2 - BITSLOP;
89 c.step = 1;
90
91 rr = mprand(rr, nbits, r, 1);
92 DRESET(&dn); dstr_putf(&dn, "%s [s]", name);
93 if ((s = pgen(dn.buf, MP_NEWSEC, rr, event, ectx, n, pgen_filter, &c,
94 rabin_iters(nbits), pgen_test, &rb)) == 0)
95 goto fail_s;
96
97 rr = mprand(rr, nbits, r, 1);
98 DRESET(&dn); dstr_putf(&dn, "%s [t]", name);
99 if ((t = pgen(dn.buf, MP_NEWSEC, rr, event, ectx, n, pgen_filter, &c,
100 rabin_iters(nbits), pgen_test, &rb)) == 0)
101 goto fail_t;
102
103 /* --- Choose a suitable value for %$r = 2it + 1$% for some %$i$% --- */
104
105 rr = mp_lsl(rr, t, 1);
106 pfilt_create(&c.f, rr);
107 rr = mp_lsl(rr, rr, BITSLOP - 1);
108 rr = mp_add(rr, rr, MP_ONE);
109 DRESET(&dn); dstr_putf(&dn, "%s [r]", name);
110 j.j = &c.f;
111 nbits += BITSLOP;
112 q = pgen(dn.buf, MP_NEW, rr, event, ectx, n, pgen_jump, &j,
113 rabin_iters(nbits), pgen_test, &rb);
114 pfilt_destroy(&c.f);
115 if (!q)
116 goto fail_r;
117
118 /* --- Select a suitable starting-point for finding %$p$% --- *
119 *
120 * This computes %$p_0 = 2(s^{r - 2} \bmod r)s - 1$%.
121 */
122
123 {
124 mpmont mm;
125
126 mpmont_create(&mm, q);
127 rr = mp_sub(rr, q, MP_TWO);
128 rr = mpmont_exp(&mm, rr, s, rr);
129 mpmont_destroy(&mm);
130 rr = mp_mul(rr, rr, s);
131 rr = mp_lsl(rr, rr, 1);
132 rr = mp_sub(rr, rr, MP_ONE);
133 }
134
135 /* --- Now find %$p = p_0 + 2jrs$% for some %$j$% --- */
136
137 {
138 mp *x;
139 x = mp_mul(MP_NEW, q, s);
140 x = mp_lsl(x, x, 1);
141 pfilt_create(f, x);
142 x = mp_lsl(x, x, BITSLOP - 1);
143 rr = mp_add(rr, rr, x);
144 mp_drop(x);
145 }
146
147 /* --- Return the result --- */
148
149 mp_drop(q);
150 mp_drop(t);
151 mp_drop(s);
152 dstr_destroy(&dn);
153 return (rr);
154
155 /* --- Tidy up if something failed --- */
156
157 fail_r:
158 mp_drop(t);
159 fail_t:
160 mp_drop(s);
161 fail_s:
162 mp_drop(rr);
163 dstr_destroy(&dn);
164 return (0);
165
166 #undef BITSLOP
167 }
168
169 /* --- @strongprime@ --- *
170 *
171 * Arguments: @const char *name@ = pointer to name root
172 * @mp *d@ = destination integer
173 * @unsigned nbits@ = number of bits wanted
174 * @grand *r@ = random number source
175 * @unsigned n@ = number of attempts to make
176 * @pgen_proc *event@ = event handler function
177 * @void *ectx@ = argument for the event handler
178 *
179 * Returns: A `strong' prime, or zero.
180 *
181 * Use: Finds `strong' primes. A strong prime %$p$% is such that
182 *
183 * * %$p - 1$% has a large prime factor %$r$%,
184 * * %$p + 1$% has a large prime factor %$s$%, and
185 * * %$r - 1$% has a large prime factor %$t$%.
186 *
187 * The numbers produced may be slightly larger than requested,
188 * by a few bits.
189 */
190
191 mp *strongprime(const char *name, mp *d, unsigned nbits, grand *r,
192 unsigned n, pgen_proc *event, void *ectx)
193 {
194 pfilt f;
195 pgen_jumpctx j;
196 rabin rb;
197
198 d = strongprime_setup(name, d, &f, nbits, r, n, event, ectx);
199 j.j = &f;
200 d = pgen(name, d, d, event, ectx, n, pgen_jump, &j,
201 rabin_iters(nbits), pgen_test, &rb);
202 pfilt_destroy(&f);
203 return (d);
204 }
205
206 /*----- That's all, folks -------------------------------------------------*/