Made gcry_prime_check more robust (and slower).
[libgcrypt.git] / cipher / primegen.c
1 /* primegen.c - prime number generator
2  * Copyright (C) 1998, 2000, 2001, 2002, 2003
3  *               2004 Free Software Foundation, Inc.
4  *
5  * This file is part of Libgcrypt.
6  *
7  * Libgcrypt is free software; you can redistribute it and/or modify
8  * it under the terms of the GNU Lesser general Public License as
9  * published by the Free Software Foundation; either version 2.1 of
10  * the License, or (at your option) any later version.
11  *
12  * Libgcrypt is distributed in the hope that it will be useful,
13  * but WITHOUT ANY WARRANTY; without even the implied warranty of
14  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
15  * GNU Lesser General Public License for more details.
16  *
17  * You should have received a copy of the GNU Lesser General Public
18  * License along with this program; if not, write to the Free Software
19  * Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA
20  *
21  * ***********************************************************************
22  * The algorithm used to generate practically save primes is due to
23  * Lim and Lee as described in the CRYPTO '97 proceedings (ISBN3540633847)
24  * page 260.
25  */
26
27 #include <config.h>
28
29 #include <stdio.h>
30 #include <stdlib.h>
31 #include <string.h>
32 #include <assert.h>
33 #include <errno.h>
34
35 #include "g10lib.h"
36 #include "mpi.h"
37 #include "cipher.h"
38
39 static gcry_mpi_t gen_prime (unsigned int nbits, int secret, int randomlevel, 
40                              int (*extra_check)(void *, gcry_mpi_t),
41                              void *extra_check_arg);
42 static int check_prime( gcry_mpi_t prime, gcry_mpi_t val_2, int rm_rounds,
43                         gcry_prime_check_func_t cb_func, void *cb_arg );
44 static int is_prime( gcry_mpi_t n, int steps, int *count );
45 static void m_out_of_n( char *array, int m, int n );
46
47 static void (*progress_cb) (void *,const char*,int,int, int );
48 static void *progress_cb_data;
49
50 /* Note: 2 is not included because it can be tested more easily by
51    looking at bit 0. The last entry in this list is marked by a zero */
52 static ushort small_prime_numbers[] = {
53     3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43,
54     47, 53, 59, 61, 67, 71, 73, 79, 83, 89, 97, 101,
55     103, 107, 109, 113, 127, 131, 137, 139, 149, 151,
56     157, 163, 167, 173, 179, 181, 191, 193, 197, 199,
57     211, 223, 227, 229, 233, 239, 241, 251, 257, 263,
58     269, 271, 277, 281, 283, 293, 307, 311, 313, 317,
59     331, 337, 347, 349, 353, 359, 367, 373, 379, 383,
60     389, 397, 401, 409, 419, 421, 431, 433, 439, 443,
61     449, 457, 461, 463, 467, 479, 487, 491, 499, 503,
62     509, 521, 523, 541, 547, 557, 563, 569, 571, 577,
63     587, 593, 599, 601, 607, 613, 617, 619, 631, 641,
64     643, 647, 653, 659, 661, 673, 677, 683, 691, 701,
65     709, 719, 727, 733, 739, 743, 751, 757, 761, 769,
66     773, 787, 797, 809, 811, 821, 823, 827, 829, 839,
67     853, 857, 859, 863, 877, 881, 883, 887, 907, 911,
68     919, 929, 937, 941, 947, 953, 967, 971, 977, 983,
69     991, 997, 1009, 1013, 1019, 1021, 1031, 1033,
70     1039, 1049, 1051, 1061, 1063, 1069, 1087, 1091,
71     1093, 1097, 1103, 1109, 1117, 1123, 1129, 1151,
72     1153, 1163, 1171, 1181, 1187, 1193, 1201, 1213,
73     1217, 1223, 1229, 1231, 1237, 1249, 1259, 1277,
74     1279, 1283, 1289, 1291, 1297, 1301, 1303, 1307,
75     1319, 1321, 1327, 1361, 1367, 1373, 1381, 1399,
76     1409, 1423, 1427, 1429, 1433, 1439, 1447, 1451,
77     1453, 1459, 1471, 1481, 1483, 1487, 1489, 1493,
78     1499, 1511, 1523, 1531, 1543, 1549, 1553, 1559,
79     1567, 1571, 1579, 1583, 1597, 1601, 1607, 1609,
80     1613, 1619, 1621, 1627, 1637, 1657, 1663, 1667,
81     1669, 1693, 1697, 1699, 1709, 1721, 1723, 1733,
82     1741, 1747, 1753, 1759, 1777, 1783, 1787, 1789,
83     1801, 1811, 1823, 1831, 1847, 1861, 1867, 1871,
84     1873, 1877, 1879, 1889, 1901, 1907, 1913, 1931,
85     1933, 1949, 1951, 1973, 1979, 1987, 1993, 1997,
86     1999, 2003, 2011, 2017, 2027, 2029, 2039, 2053,
87     2063, 2069, 2081, 2083, 2087, 2089, 2099, 2111,
88     2113, 2129, 2131, 2137, 2141, 2143, 2153, 2161,
89     2179, 2203, 2207, 2213, 2221, 2237, 2239, 2243,
90     2251, 2267, 2269, 2273, 2281, 2287, 2293, 2297,
91     2309, 2311, 2333, 2339, 2341, 2347, 2351, 2357,
92     2371, 2377, 2381, 2383, 2389, 2393, 2399, 2411,
93     2417, 2423, 2437, 2441, 2447, 2459, 2467, 2473,
94     2477, 2503, 2521, 2531, 2539, 2543, 2549, 2551,
95     2557, 2579, 2591, 2593, 2609, 2617, 2621, 2633,
96     2647, 2657, 2659, 2663, 2671, 2677, 2683, 2687,
97     2689, 2693, 2699, 2707, 2711, 2713, 2719, 2729,
98     2731, 2741, 2749, 2753, 2767, 2777, 2789, 2791,
99     2797, 2801, 2803, 2819, 2833, 2837, 2843, 2851,
100     2857, 2861, 2879, 2887, 2897, 2903, 2909, 2917,
101     2927, 2939, 2953, 2957, 2963, 2969, 2971, 2999,
102     3001, 3011, 3019, 3023, 3037, 3041, 3049, 3061,
103     3067, 3079, 3083, 3089, 3109, 3119, 3121, 3137,
104     3163, 3167, 3169, 3181, 3187, 3191, 3203, 3209,
105     3217, 3221, 3229, 3251, 3253, 3257, 3259, 3271,
106     3299, 3301, 3307, 3313, 3319, 3323, 3329, 3331,
107     3343, 3347, 3359, 3361, 3371, 3373, 3389, 3391,
108     3407, 3413, 3433, 3449, 3457, 3461, 3463, 3467,
109     3469, 3491, 3499, 3511, 3517, 3527, 3529, 3533,
110     3539, 3541, 3547, 3557, 3559, 3571, 3581, 3583,
111     3593, 3607, 3613, 3617, 3623, 3631, 3637, 3643,
112     3659, 3671, 3673, 3677, 3691, 3697, 3701, 3709,
113     3719, 3727, 3733, 3739, 3761, 3767, 3769, 3779,
114     3793, 3797, 3803, 3821, 3823, 3833, 3847, 3851,
115     3853, 3863, 3877, 3881, 3889, 3907, 3911, 3917,
116     3919, 3923, 3929, 3931, 3943, 3947, 3967, 3989,
117     4001, 4003, 4007, 4013, 4019, 4021, 4027, 4049,
118     4051, 4057, 4073, 4079, 4091, 4093, 4099, 4111,
119     4127, 4129, 4133, 4139, 4153, 4157, 4159, 4177,
120     4201, 4211, 4217, 4219, 4229, 4231, 4241, 4243,
121     4253, 4259, 4261, 4271, 4273, 4283, 4289, 4297,
122     4327, 4337, 4339, 4349, 4357, 4363, 4373, 4391,
123     4397, 4409, 4421, 4423, 4441, 4447, 4451, 4457,
124     4463, 4481, 4483, 4493, 4507, 4513, 4517, 4519,
125     4523, 4547, 4549, 4561, 4567, 4583, 4591, 4597,
126     4603, 4621, 4637, 4639, 4643, 4649, 4651, 4657,
127     4663, 4673, 4679, 4691, 4703, 4721, 4723, 4729,
128     4733, 4751, 4759, 4783, 4787, 4789, 4793, 4799,
129     4801, 4813, 4817, 4831, 4861, 4871, 4877, 4889,
130     4903, 4909, 4919, 4931, 4933, 4937, 4943, 4951,
131     4957, 4967, 4969, 4973, 4987, 4993, 4999,
132     0
133 };
134 static int no_of_small_prime_numbers = DIM (small_prime_numbers) - 1;
135
136 void
137 _gcry_register_primegen_progress ( void (*cb)(void *,const char*,int,int,int),
138                                    void *cb_data )
139 {
140   progress_cb = cb;
141   progress_cb_data = cb_data;
142 }
143
144
145 static void
146 progress( int c )
147 {
148   if ( progress_cb )
149     progress_cb ( progress_cb_data, "primegen", c, 0, 0 );
150 }
151
152
153 /****************
154  * Generate a prime number (stored in secure memory)
155  */
156 gcry_mpi_t
157 _gcry_generate_secret_prime (unsigned int nbits,
158                              int (*extra_check)(void*, gcry_mpi_t),
159                              void *extra_check_arg)
160 {
161   gcry_mpi_t prime;
162
163   prime = gen_prime( nbits, 1, 2, extra_check, extra_check_arg);
164   progress('\n');
165   return prime;
166 }
167
168 gcry_mpi_t
169 _gcry_generate_public_prime( unsigned int nbits,
170                              int (*extra_check)(void*, gcry_mpi_t),
171                              void *extra_check_arg)
172 {
173   gcry_mpi_t prime;
174
175   prime = gen_prime( nbits, 0, 2, extra_check, extra_check_arg );
176   progress('\n');
177   return prime;
178 }
179
180
181 /****************
182  * We do not need to use the strongest RNG because we gain no extra
183  * security from it - The prime number is public and we could also
184  * offer the factors for those who are willing to check that it is
185  * indeed a strong prime.  With ALL_FACTORS set to true all afcors of
186  * prime-1 are returned in FACTORS.
187  *
188  * mode 0: Standard
189  *      1: Make sure that at least one factor is of size qbits.
190  */
191 static gcry_err_code_t
192 prime_generate_internal (int mode,
193                          gcry_mpi_t *prime_generated, unsigned int pbits,
194                          unsigned int qbits, gcry_mpi_t g,
195                          gcry_mpi_t **ret_factors,
196                          gcry_random_level_t randomlevel, unsigned int flags,
197                          int all_factors,
198                          gcry_prime_check_func_t cb_func, void *cb_arg)
199 {
200   gcry_err_code_t err = 0;
201   gcry_mpi_t *factors_new = NULL; /* Factors to return to the
202                                      caller.  */
203   gcry_mpi_t *factors = NULL;   /* Current factors.  */
204   gcry_mpi_t *pool = NULL;      /* Pool of primes.  */
205   unsigned char *perms = NULL;  /* Permutations of POOL.  */
206   gcry_mpi_t q_factor = NULL;   /* Used if QBITS is non-zero.  */
207   unsigned int fbits = 0;       /* Length of prime factors.  */
208   unsigned int n = 0;           /* Number of factors.  */
209   unsigned int m = 0;           /* Number of primes in pool.  */
210   gcry_mpi_t q = NULL;          /* First prime factor.  */
211   gcry_mpi_t prime = NULL;      /* Prime candidate.  */
212   unsigned int nprime = 0;      /* Bits of PRIME.  */
213   unsigned int req_qbits;       /* The original QBITS value.  */
214   gcry_mpi_t val_2;             /* For check_prime().  */
215   unsigned int is_secret = (flags & GCRY_PRIME_FLAG_SECRET);
216   unsigned int count1 = 0, count2 = 0;
217   unsigned int i = 0, j = 0;
218
219   if (pbits < 48)
220     return GPG_ERR_INV_ARG;
221
222   /* If QBITS is not given, assume a reasonable value. */
223   if (!qbits)
224     qbits = pbits / 3;
225
226   req_qbits = qbits;
227
228   /* Find number of needed prime factors.  */
229   for (n = 1; (pbits - qbits - 1) / n  >= qbits; n++)
230     ;
231   n--;
232
233   val_2 = mpi_alloc_set_ui (2);
234
235   if ((! n) || ((mode == 1) && (n < 2)))
236     {
237       err = GPG_ERR_INV_ARG;
238       goto leave;
239     }
240
241   if (mode == 1)
242     {
243       n--;
244       fbits = (pbits - 2 * req_qbits -1) / n;
245       qbits =  pbits - req_qbits - n * fbits;
246     }
247   else
248     {
249       fbits = (pbits - req_qbits -1) / n;
250       qbits = pbits - n * fbits;
251     }
252   
253   if (DBG_CIPHER)
254     log_debug ("gen prime: pbits=%u qbits=%u fbits=%u/%u n=%d\n",
255                pbits, req_qbits, qbits, fbits, n);
256
257   prime = gcry_mpi_new (pbits);
258
259   /* Generate first prime factor.  */
260   q = gen_prime (qbits, is_secret, randomlevel, NULL, NULL);
261   
262   if (mode == 1)
263     q_factor = gen_prime (req_qbits, is_secret, randomlevel, NULL, NULL);
264   
265   /* Allocate an array to hold the factors + 2 for later usage.  */
266   factors = gcry_calloc (n + 2, sizeof (*factors));
267   if (!factors)
268     {
269       err = gpg_err_code_from_errno (errno);
270       goto leave;
271     }
272       
273   /* Make a pool of 3n+5 primes (this is an arbitrary value).  */
274   m = n * 3 + 5;
275   if (mode == 1) /* Need some more (for e.g. DSA).  */
276     m += 5;
277   if (m < 25)
278     m = 25;
279   pool = gcry_calloc (m , sizeof (*pool));
280   if (! pool)
281     {
282       err = gpg_err_code_from_errno (errno);
283       goto leave;
284     }
285
286   /* Permutate over the pool of primes.  */
287   do
288     {
289     next_try:
290       if (! perms)
291         {
292           /* Allocate new primes.  */
293           for(i = 0; i < m; i++)
294             {
295               mpi_free (pool[i]);
296               pool[i] = NULL;
297             }
298
299           /* Init m_out_of_n().  */
300           perms = gcry_calloc (1, m);
301           if (! perms)
302             {
303               err = gpg_err_code_from_errno (errno);
304               goto leave;
305             }
306           for(i = 0; i < n; i++)
307             {
308               perms[i] = 1;
309               pool[i] = gen_prime (fbits, is_secret,
310                                    randomlevel, NULL, NULL);
311               factors[i] = pool[i];
312             }
313         }
314       else
315         {
316           m_out_of_n (perms, n, m);
317           for (i = j = 0; (i < m) && (j < n); i++)
318             if (perms[i])
319               {
320                 if(! pool[i])
321                   pool[i] = gen_prime (fbits, 0, 1, NULL, NULL);
322                 factors[j++] = pool[i];
323               }
324           if (i == n)
325             {
326               gcry_free (perms);
327               perms = NULL;
328               progress ('!');
329               goto next_try;    /* Allocate new primes.  */
330             }
331         }
332
333         /* Generate next prime candidate:
334            p = 2 * q [ * q_factor] * factor_0 * factor_1 * ... * factor_n + 1. 
335         */
336         mpi_set (prime, q);
337         mpi_mul_ui (prime, prime, 2);
338         if (mode == 1)
339           mpi_mul (prime, prime, q_factor);
340         for(i = 0; i < n; i++)
341           mpi_mul (prime, prime, factors[i]);
342         mpi_add_ui (prime, prime, 1);
343         nprime = mpi_get_nbits (prime);
344
345         if (nprime < pbits)
346           {
347             if (++count1 > 20)
348               {
349                 count1 = 0;
350                 qbits++;
351                 progress('>');
352                 mpi_free (q);
353                 q = gen_prime (qbits, 0, 0, NULL, NULL);
354                 goto next_try;
355               }
356           }
357         else
358           count1 = 0;
359         
360         if (nprime > pbits)
361           {
362             if (++count2 > 20)
363               {
364                 count2 = 0;
365                 qbits--;
366                 progress('<');
367                 mpi_free (q);
368                 q = gen_prime (qbits, 0, 0, NULL, NULL);
369                 goto next_try;
370               }
371           }
372         else
373           count2 = 0;
374     }
375   while (! ((nprime == pbits) && check_prime (prime, val_2, 5,
376                                               cb_func, cb_arg)));
377
378   if (DBG_CIPHER)
379     {
380       progress ('\n');
381       log_mpidump ("prime    : ", prime);
382       log_mpidump ("factor  q: ", q);
383       if (mode == 1)
384         log_mpidump ("factor q0: ", q_factor);
385       for (i = 0; i < n; i++)
386         log_mpidump ("factor pi: ", factors[i]);
387       log_debug ("bit sizes: prime=%u, q=%u",
388                  mpi_get_nbits (prime), mpi_get_nbits (q));
389       if (mode == 1)
390         log_debug (", q0=%u", mpi_get_nbits (q_factor));
391       for (i = 0; i < n; i++)
392         log_debug (", p%d=%u", i, mpi_get_nbits (factors[i]));
393       progress('\n');
394     }
395
396   if (ret_factors)
397     {
398       /* Caller wants the factors.  */
399       factors_new = gcry_calloc (n + 4, sizeof (*factors_new));
400       if (! factors_new)
401         {
402           err = gpg_err_code_from_errno (errno);
403           goto leave;
404         }
405
406       if (all_factors)
407         {
408           i = 0;
409           factors_new[i++] = gcry_mpi_set_ui (NULL, 2);
410           factors_new[i++] = mpi_copy (q);
411           if (mode == 1)
412             factors_new[i++] = mpi_copy (q_factor);
413           for(j=0; j < n; j++)
414             factors_new[i++] = mpi_copy (factors[j]);
415         }
416       else
417         {
418           i = 0;
419           if (mode == 1)
420             {
421               factors_new[i++] = mpi_copy (q_factor);
422               for (; i <= n; i++)
423                 factors_new[i] = mpi_copy (factors[i]);
424             }
425           else
426             for (; i < n; i++ )
427               factors_new[i] = mpi_copy (factors[i]);
428         }
429     }
430   
431   if (g)
432     {
433       /* Create a generator (start with 3).  */
434       gcry_mpi_t tmp = mpi_alloc (mpi_get_nlimbs (prime));
435       gcry_mpi_t b = mpi_alloc (mpi_get_nlimbs (prime));
436       gcry_mpi_t pmin1 = mpi_alloc (mpi_get_nlimbs (prime));
437       
438       if (mode == 1)
439         err = GPG_ERR_NOT_IMPLEMENTED;
440       else
441         {
442           factors[n] = q;
443           factors[n + 1] = mpi_alloc_set_ui (2);
444           mpi_sub_ui (pmin1, prime, 1);
445           mpi_set_ui (g, 2);
446           do
447             {
448               mpi_add_ui (g, g, 1);
449               if (DBG_CIPHER)
450                 {
451                   log_debug ("checking g:");
452                   gcry_mpi_dump (g);
453                   log_printf ("\n");
454                 }
455               else
456                 progress('^');
457               for (i = 0; i < n + 2; i++)
458                 {
459                   mpi_fdiv_q (tmp, pmin1, factors[i]);
460                   /* No mpi_pow(), but it is okay to use this with mod
461                      prime.  */
462                   gcry_mpi_powm (b, g, tmp, prime);
463                   if (! mpi_cmp_ui (b, 1))
464                     break;
465                 }
466               if (DBG_CIPHER)
467                 progress('\n');
468             } 
469           while (i < n + 2);
470
471           mpi_free (factors[n+1]);
472           mpi_free (tmp);
473           mpi_free (b);
474           mpi_free (pmin1);
475         }
476     }
477   
478   if (! DBG_CIPHER)
479     progress ('\n');
480
481
482  leave:
483   if (pool)
484     {
485       for(i = 0; i < m; i++)
486         mpi_free (pool[i]);
487       gcry_free (pool);
488     }
489   if (factors)
490     gcry_free (factors);  /* Factors are shallow copies.  */
491   if (perms)
492     gcry_free (perms);
493
494   mpi_free (val_2);
495   mpi_free (q);
496   mpi_free (q_factor);
497
498   if (! err)
499     {
500       *prime_generated = prime;
501       if (ret_factors)
502         *ret_factors = factors_new;
503     }
504   else
505     {
506       if (factors_new)
507         {
508           for (i = 0; factors_new[i]; i++)
509             mpi_free (factors_new[i]);
510           gcry_free (factors_new);
511         }
512       mpi_free (prime);
513     }
514
515   return err;
516 }
517
518 gcry_mpi_t
519 _gcry_generate_elg_prime (int mode, unsigned pbits, unsigned qbits,
520                           gcry_mpi_t g, gcry_mpi_t **ret_factors)
521 {
522   gcry_err_code_t err = GPG_ERR_NO_ERROR;
523   gcry_mpi_t prime = NULL;
524   
525   err = prime_generate_internal (mode, &prime, pbits, qbits, g,
526                                  ret_factors, GCRY_WEAK_RANDOM, 0, 0,
527                                  NULL, NULL);
528
529   return prime;
530 }
531
532 static gcry_mpi_t
533 gen_prime (unsigned int nbits, int secret, int randomlevel, 
534            int (*extra_check)(void *, gcry_mpi_t), void *extra_check_arg)
535 {
536   gcry_mpi_t prime, ptest, pminus1, val_2, val_3, result;
537   int i;
538   unsigned x, step;
539   unsigned count1, count2;
540   int *mods;
541   
542 /*   if (  DBG_CIPHER ) */
543 /*     log_debug ("generate a prime of %u bits ", nbits ); */
544
545   if (nbits < 16)
546     log_fatal ("can't generate a prime with less than %d bits\n", 16);
547
548   mods = gcry_xmalloc( no_of_small_prime_numbers * sizeof *mods );
549   /* Make nbits fit into gcry_mpi_t implementation. */
550   val_2  = mpi_alloc_set_ui( 2 );
551   val_3 = mpi_alloc_set_ui( 3);
552   prime  = secret? gcry_mpi_snew ( nbits ): gcry_mpi_new ( nbits );
553   result = mpi_alloc_like( prime );
554   pminus1= mpi_alloc_like( prime );
555   ptest  = mpi_alloc_like( prime );
556   count1 = count2 = 0;
557   for (;;)
558     {  /* try forvever */
559       int dotcount=0;
560       
561       /* generate a random number */
562       gcry_mpi_randomize( prime, nbits, randomlevel );
563       
564       /* Set high order bit to 1, set low order bit to 1.  If we are
565          generating a secret prime we are most probably doing that
566          for RSA, to make sure that the modulus does have the
567          requested key size we set the 2 high order bits. */
568       mpi_set_highbit (prime, nbits-1);
569       if (secret)
570         mpi_set_bit (prime, nbits-2);
571       mpi_set_bit(prime, 0);
572       
573       /* Calculate all remainders. */
574       for (i=0; (x = small_prime_numbers[i]); i++ )
575         mods[i] = mpi_fdiv_r_ui(NULL, prime, x);
576       
577       /* Now try some primes starting with prime. */
578       for(step=0; step < 20000; step += 2 ) 
579         {
580           /* Check against all the small primes we have in mods. */
581           count1++;
582           for (i=0; (x = small_prime_numbers[i]); i++ ) 
583             {
584               while ( mods[i] + step >= x )
585                 mods[i] -= x;
586               if ( !(mods[i] + step) )
587                 break;
588             }
589           if ( x )
590             continue;   /* Found a multiple of an already known prime. */
591           
592           mpi_add_ui( ptest, prime, step );
593
594           /* Do a fast Fermat test now. */
595           count2++;
596           mpi_sub_ui( pminus1, ptest, 1);
597           gcry_mpi_powm( result, val_2, pminus1, ptest );
598           if ( !mpi_cmp_ui( result, 1 ) )
599             { 
600               /* Not composite, perform stronger tests */
601               if (is_prime(ptest, 5, &count2 ))
602                 {
603                   if (!mpi_test_bit( ptest, nbits-1-secret ))
604                     {
605                       progress('\n');
606                       log_debug ("overflow in prime generation\n");
607                       break; /* Stop loop, continue with a new prime. */
608                     }
609
610                   if (extra_check && extra_check (extra_check_arg, ptest))
611                     { 
612                       /* The extra check told us that this prime is
613                          not of the caller's taste. */
614                       progress ('/');
615                     }
616                   else
617                     { 
618                       /* Got it. */
619                       mpi_free(val_2);
620                       mpi_free(val_3);
621                       mpi_free(result);
622                       mpi_free(pminus1);
623                       mpi_free(prime);
624                       gcry_free(mods);
625                       return ptest; 
626                     }
627                 }
628             }
629           if (++dotcount == 10 )
630             {
631               progress('.');
632               dotcount = 0;
633             }
634         }
635       progress(':'); /* restart with a new random value */
636     }
637 }
638
639 /****************
640  * Returns: true if this may be a prime
641  * RM_ROUNDS gives the number of Rabin-Miller tests to run.
642  */
643 static int
644 check_prime( gcry_mpi_t prime, gcry_mpi_t val_2, int rm_rounds,
645              gcry_prime_check_func_t cb_func, void *cb_arg)
646 {
647   int i;
648   unsigned int x;
649   int count=0;
650
651   /* Check against small primes. */
652   for (i=0; (x = small_prime_numbers[i]); i++ )
653     {
654       if ( mpi_divisible_ui( prime, x ) )
655         return 0;
656     }
657
658   /* A quick Fermat test. */
659   {
660     gcry_mpi_t result = mpi_alloc_like( prime );
661     gcry_mpi_t pminus1 = mpi_alloc_like( prime );
662     mpi_sub_ui( pminus1, prime, 1);
663     gcry_mpi_powm( result, val_2, pminus1, prime );
664     mpi_free( pminus1 );
665     if ( mpi_cmp_ui( result, 1 ) )
666       { 
667         /* Is composite. */
668         mpi_free( result );
669         progress('.');
670         return 0;
671       }
672     mpi_free( result );
673   }
674
675   if (!cb_func || cb_func (cb_arg, GCRY_PRIME_CHECK_AT_MAYBE_PRIME, prime))
676     {
677       /* Perform stronger tests. */
678       if ( is_prime( prime, rm_rounds, &count ) )
679         {
680           if (!cb_func
681               || cb_func (cb_arg, GCRY_PRIME_CHECK_AT_GOT_PRIME, prime))
682             return 1; /* Probably a prime. */
683         }
684     }
685   progress('.');
686   return 0;
687 }
688
689
690 /*
691  * Return true if n is probably a prime
692  */
693 static int
694 is_prime (gcry_mpi_t n, int steps, int *count)
695 {
696   gcry_mpi_t x = mpi_alloc( mpi_get_nlimbs( n ) );
697   gcry_mpi_t y = mpi_alloc( mpi_get_nlimbs( n ) );
698   gcry_mpi_t z = mpi_alloc( mpi_get_nlimbs( n ) );
699   gcry_mpi_t nminus1 = mpi_alloc( mpi_get_nlimbs( n ) );
700   gcry_mpi_t a2 = mpi_alloc_set_ui( 2 );
701   gcry_mpi_t q;
702   unsigned i, j, k;
703   int rc = 0;
704   unsigned nbits = mpi_get_nbits( n );
705
706   if (steps < 5) /* Make sure that we do at least 5 rounds. */
707     steps = 5; 
708
709   mpi_sub_ui( nminus1, n, 1 );
710
711   /* Find q and k, so that n = 1 + 2^k * q . */
712   q = mpi_copy ( nminus1 );
713   k = mpi_trailing_zeros ( q );
714   mpi_tdiv_q_2exp (q, q, k);
715
716   for (i=0 ; i < steps; i++ )
717     {
718       ++*count;
719       if( !i )
720         {
721           mpi_set_ui( x, 2 );
722         }
723       else
724         {
725           gcry_mpi_randomize( x, nbits, GCRY_WEAK_RANDOM );
726
727           /* Make sure that the number is smaller than the prime and
728              keep the randomness of the high bit. */
729           if ( mpi_test_bit ( x, nbits-2) )
730             {
731               mpi_set_highbit ( x, nbits-2); /* Clear all higher bits. */
732             }
733           else
734             {
735               mpi_set_highbit( x, nbits-2 );
736               mpi_clear_bit( x, nbits-2 );
737             }
738           assert ( mpi_cmp( x, nminus1 ) < 0 && mpi_cmp_ui( x, 1 ) > 0 );
739         }
740       gcry_mpi_powm ( y, x, q, n);
741       if ( mpi_cmp_ui(y, 1) && mpi_cmp( y, nminus1 ) )
742         {
743           for ( j=1; j < k && mpi_cmp( y, nminus1 ); j++ )
744             {
745               gcry_mpi_powm(y, y, a2, n);
746               if( !mpi_cmp_ui( y, 1 ) )
747                 goto leave; /* Not a prime. */
748             }
749           if (mpi_cmp( y, nminus1 ) )
750             goto leave; /* Not a prime. */
751         }
752       progress('+');
753     }
754   rc = 1; /* May be a prime. */
755
756  leave:
757   mpi_free( x );
758   mpi_free( y );
759   mpi_free( z );
760   mpi_free( nminus1 );
761   mpi_free( q );
762   mpi_free( a2 );
763
764   return rc;
765 }
766
767
768 static void
769 m_out_of_n ( char *array, int m, int n )
770 {
771   int i=0, i1=0, j=0, jp=0,  j1=0, k1=0, k2=0;
772
773   if( !m || m >= n )
774     return;
775
776   if( m == 1 )
777     { 
778       /* Special case. */
779       for (i=0; i < n; i++ )
780         {
781           if( array[i] )
782             {
783               array[i++] = 0;
784               if( i >= n )
785                 i = 0;
786               array[i] = 1;
787               return;
788             }
789         }
790       BUG();
791     }
792
793   for (j=1; j < n; j++ )
794     {
795       if ( array[n-1] == array[n-j-1])
796         continue;
797       j1 = j;
798       break;
799     }
800
801   if ( (m & 1) )
802     {
803       /* M is odd. */
804       if( array[n-1] )
805         {
806           if( j1 & 1 )
807             {
808               k1 = n - j1;
809               k2 = k1+2;
810               if( k2 > n )
811                 k2 = n;
812               goto leave;
813             }
814           goto scan;
815         }
816       k2 = n - j1 - 1;
817       if( k2 == 0 )
818         {
819           k1 = i;
820           k2 = n - j1;
821         }
822       else if( array[k2] && array[k2-1] )
823         k1 = n;
824       else
825         k1 = k2 + 1;
826     }
827   else 
828     {
829       /* M is even. */
830       if( !array[n-1] )
831         {
832           k1 = n - j1;
833           k2 = k1 + 1;
834           goto leave;
835         }
836         
837       if( !(j1 & 1) )
838         {
839           k1 = n - j1;
840           k2 = k1+2;
841           if( k2 > n )
842             k2 = n;
843           goto leave;
844         }
845     scan:
846       jp = n - j1 - 1;
847       for (i=1; i <= jp; i++ ) 
848         {
849           i1 = jp + 2 - i;
850           if( array[i1-1]  )
851             {
852               if( array[i1-2] )
853                 {
854                   k1 = i1 - 1;
855                   k2 = n - j1;
856                 }
857               else
858                 {
859                   k1 = i1 - 1;
860                   k2 = n + 1 - j1;
861                 }
862               goto leave;
863             }
864         }
865       k1 = 1;
866       k2 = n + 1 - m;
867     }
868  leave:
869   array[k1-1] = !array[k1-1];
870   array[k2-1] = !array[k2-1];
871 }
872
873
874 /* Generate a new prime number of PRIME_BITS bits and store it in
875    PRIME.  If FACTOR_BITS is non-zero, one of the prime factors of
876    (prime - 1) / 2 must be FACTOR_BITS bits long.  If FACTORS is
877    non-zero, allocate a new, NULL-terminated array holding the prime
878    factors and store it in FACTORS.  FLAGS might be used to influence
879    the prime number generation process.  */
880 gcry_error_t
881 gcry_prime_generate (gcry_mpi_t *prime, unsigned int prime_bits,
882                      unsigned int factor_bits, gcry_mpi_t **factors,
883                      gcry_prime_check_func_t cb_func, void *cb_arg,
884                      gcry_random_level_t random_level,
885                      unsigned int flags)
886 {
887   gcry_err_code_t err = GPG_ERR_NO_ERROR;
888   gcry_mpi_t *factors_generated = NULL;
889   gcry_mpi_t prime_generated = NULL;
890   unsigned int mode = 0;
891
892   if (!prime)
893     return gpg_error (GPG_ERR_INV_ARG);
894   *prime = NULL; 
895
896   if (flags & GCRY_PRIME_FLAG_SPECIAL_FACTOR)
897     mode = 1;
898
899   /* Generate.  */
900   err = prime_generate_internal (mode, &prime_generated, prime_bits,
901                                  factor_bits, NULL,
902                                  factors? &factors_generated : NULL,
903                                  random_level, flags, 1,
904                                  cb_func, cb_arg);
905
906   if (! err)
907     if (cb_func)
908       {
909         /* Additional check. */
910         if ( !cb_func (cb_arg, GCRY_PRIME_CHECK_AT_FINISH, prime_generated))
911           {
912             /* Failed, deallocate resources.  */
913             unsigned int i;
914
915             mpi_free (prime_generated);
916             if (factors)
917               {
918                 for (i = 0; factors_generated[i]; i++)
919                   mpi_free (factors_generated[i]);
920                 gcry_free (factors_generated);
921               }
922             err = GPG_ERR_GENERAL; 
923           }
924       }
925
926   if (! err)
927     {
928       if (factors)
929         *factors = factors_generated;
930       *prime = prime_generated;
931     }
932
933   return gcry_error (err);
934 }
935
936 /* Check wether the number X is prime.  */
937 gcry_error_t
938 gcry_prime_check (gcry_mpi_t x, unsigned int flags)
939 {
940   gcry_err_code_t err = GPG_ERR_NO_ERROR;
941   gcry_mpi_t val_2 = mpi_alloc_set_ui (2); /* Used by the Fermat test. */
942
943   /* We use 64 rounds because the prime we are going to test is not
944      guaranteed to be a random one. */
945   if (! check_prime (x, val_2, 64, NULL, NULL))
946     err = GPG_ERR_NO_PRIME;
947
948   mpi_free (val_2);
949
950   return gcry_error (err);
951 }
952
953 /* Find a generator for PRIME where the factorization of (prime-1) is
954    in the NULL terminated array FACTORS. Return the generator as a
955    newly allocated MPI in R_G.  If START_G is not NULL, use this as s
956    atart for the search. Returns 0 on success.*/
957 gcry_error_t
958 gcry_prime_group_generator (gcry_mpi_t *r_g,
959                             gcry_mpi_t prime, gcry_mpi_t *factors,
960                             gcry_mpi_t start_g)
961 {
962   gcry_mpi_t tmp = gcry_mpi_new (0);
963   gcry_mpi_t b = gcry_mpi_new (0);
964   gcry_mpi_t pmin1 = gcry_mpi_new (0);
965   gcry_mpi_t g = start_g? gcry_mpi_copy (start_g) : gcry_mpi_set_ui (NULL, 3);
966   int first = 1;
967   int i, n;
968
969   if (!factors || !r_g || !prime)
970     return gpg_error (GPG_ERR_INV_ARG);
971   *r_g = NULL; 
972
973   for (n=0; factors[n]; n++)
974     ;
975   if (n < 2)
976     return gpg_error (GPG_ERR_INV_ARG);
977
978   /* Extra sanity check - usually disabled. */  
979 /*   mpi_set (tmp, factors[0]); */
980 /*   for(i = 1; i < n; i++) */
981 /*     mpi_mul (tmp, tmp, factors[i]); */
982 /*   mpi_add_ui (tmp, tmp, 1); */
983 /*   if (mpi_cmp (prime, tmp)) */
984 /*     return gpg_error (GPG_ERR_INV_ARG); */
985   
986   gcry_mpi_sub_ui (pmin1, prime, 1);      
987   do         
988     {
989       if (first)
990         first = 0;
991       else
992         gcry_mpi_add_ui (g, g, 1);
993       
994       if (DBG_CIPHER)
995         {
996           log_debug ("checking g:");
997           gcry_mpi_dump (g);
998           log_debug ("\n");
999         }
1000       else
1001         progress('^');
1002       
1003       for (i = 0; i < n; i++)
1004         {
1005           mpi_fdiv_q (tmp, pmin1, factors[i]);
1006           gcry_mpi_powm (b, g, tmp, prime);
1007           if (! mpi_cmp_ui (b, 1))
1008             break;
1009         }
1010       if (DBG_CIPHER)
1011         progress('\n');
1012     }
1013   while (i < n);
1014   
1015   gcry_mpi_release (tmp);
1016   gcry_mpi_release (b); 
1017   gcry_mpi_release (pmin1); 
1018   *r_g = g; 
1019
1020   return 0; 
1021 }
1022
1023 /* Convenience function to release the factors array. */
1024 void
1025 gcry_prime_release_factors (gcry_mpi_t *factors)
1026 {
1027   if (factors)
1028     {
1029       int i;
1030       
1031       for (i=0; factors[i]; i++)
1032         mpi_free (factors[i]);
1033       gcry_free (factors);
1034     }
1035 }