git.distorted.org.uk Git - u/mdw/catacomb/blob - limlee.c

   1 /* -*-c-*-
   2  *
   3  * $Id: limlee.c,v 1.2 2000/07/26 18:00:00 mdw Exp $
   4  *
   5  * Generate Lim-Lee primes
   6  *
   7  * (c) 2000 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 /*----- Revision history --------------------------------------------------*
  31  *
  32  * $Log: limlee.c,v $
  33  * Revision 1.2  2000/07/26 18:00:00  mdw
  34  * No footer line!
  35  *
  36  * Revision 1.1  2000/07/09 21:30:58  mdw
  37  * Lim-Lee prime generation.
  38  *
  39  */
  40
  41 /*----- Header files ------------------------------------------------------*/
  42
  43 #include <mLib/alloc.h>
  44 #include <mLib/dstr.h>
  45
  46 #include "limlee.h"
  47 #include "mpmul.h"
  48 #include "mprand.h"
  49 #include "pgen.h"
  50 #include "primorial.h"
  51 #include "rabin.h"
  52
  53 /*----- Main code ---------------------------------------------------------*/
  54
  55 /* --- @limlee@ --- *
  56  *
  57  * Arguments:   @const char *name@ = pointer to name root
  58  *              @mp *d@ = pointer to destination integer
  59  *              @mp *newp@ = how to generate factor primes
  60  *              @unsigned ql@ = size of individual factors
  61  *              @unsigned pl@ = size of large prime
  62  *              @grand *r@ = a random number source
  63  *              @unsigned on@ = number of outer attempts to make
  64  *              @pgen_proc *oev@ = outer event handler function
  65  *              @void *oec@ = argument for the outer event handler
  66  *              @pgen_proc *iev@ = inner event handler function
  67  *              @void *iec@ = argument for the inner event handler
  68  *              @size_t *nf@, @mp ***f@ = output array for factors
  69  *
  70  * Returns:     A Lim-Lee prime, or null if generation failed.
  71  *
  72  * Use:         Generates Lim-Lee primes.  A Lim-Lee prime %$p$% is one which
  73  *              satisfies %$p = 2 \prod_i q_i + 1$%, where all of the %$q_i$%
  74  *              are large enough to resist square-root discrete log
  75  *              algorithms.
  76  *
  77  *              If we succeed, and @f@ is non-null, we write the array of
  78  *              factors chosen to @f@ for the benefit of the caller.
  79  */
  80
  81 static void comb_init(octet *c, unsigned n, unsigned r)
  82 {
  83   memset(c, 0, n - r);
  84   memset(c + (n - r), 1, r);
  85 }
  86
  87 static int comb_next(octet *c, unsigned n, unsigned r)
  88 {
  89   unsigned g = 0;
  90
  91   /* --- How the algorithm works --- *
  92    *
  93    * Set bits start at the end and work their way towards the start.
  94    * Excepting bits already at the start, we scan for the lowest set bit, and
  95    * move it one place nearer the start.  A group of bits at the start are
  96    * counted and reset just below the `moved' bit.  If there is no moved bit
  97    * then we're done.
  98    */
  99
 100   /* --- Count the group at the start --- */
 101
 102   for (; *c; c++) {
 103     g++;
 104     *c = 0;
 105   }
 106   if (g == r)
 107     return (0);
 108
 109   /* --- Move the next bit down one --- *
 110    *
 111    * There must be one, because otherwise we'd have counted %$r$% bits
 112    * earlier.
 113    */
 114
 115   for (; !*c; c++)
 116     ;
 117   *c = 0;
 118   g++;
 119   for (; g; g--)
 120     *--c = 1;
 121   return (1);
 122 }
 123
 124 mp *limlee(const char *name, mp *d, mp *newp,
 125            unsigned ql, unsigned pl, grand *r,
 126            unsigned on, pgen_proc *oev, void *oec,
 127            pgen_proc *iev, void *iec,
 128            size_t *nf, mp ***f)
 129 {
 130   dstr dn = DSTR_INIT;
 131   unsigned qql;
 132   mp *qq = 0;
 133   unsigned nn;
 134   unsigned mm;
 135   mp **v;
 136   octet *c;
 137   unsigned i;
 138   unsigned long seq = 0;
 139   pgen_event ev;
 140   unsigned ntest;
 141   rabin rb;
 142   pgen_filterctx pf;
 143
 144   /* --- First of all, decide on a number of factors to make --- */
 145
 146   nn = pl/ql;
 147   qql = pl%ql;
 148   if (!nn)
 149     return (0);
 150   else if (qql && nn > 1) {
 151     nn--;
 152     qql += ql;
 153   }
 154
 155   /* --- Now decide on how many primes I'll actually generate --- *
 156    *
 157    * The formula %$m = \max(3 n + 5, 25)$% comes from GPG's prime generation
 158    * library.
 159    */
 160
 161   mm = nn * 3 + 5;
 162   if (mm < 25)
 163     mm = 25;
 164
 165   /* --- Now allocate the working memory --- */
 166
 167   primorial_setup();
 168   v = xmalloc(mm * sizeof(mp *));
 169   c = xmalloc(mm);
 170
 171   /* --- Initialize everything and try to find a prime --- */
 172
 173   ev.name = name;
 174   ev.m = 0;
 175   ev.steps = on;
 176   ev.tests = ntest = rabin_iters(pl);
 177   ev.r = r;
 178
 179   if (oev && oev(PGEN_BEGIN, &ev, oec) == PGEN_ABORT)
 180     goto fail;
 181
 182   if (qql) {
 183     dstr_putf(&dn, "%s [+]", name);
 184     qq = mprand(d, qql, r, 1);
 185     pf.step = 2;
 186     qq = pgen(dn.buf, qq, qq, iev, iec,
 187               0, pgen_filter, &pf, rabin_iters(qql), pgen_test, &rb);
 188   }
 189
 190 again:
 191   comb_init(c, mm, nn);
 192   for (i = 0; i < mm; i++)
 193     v[i] = 0;
 194
 195   /* --- The main combinations loop --- */
 196
 197   do {
 198     mpmul mmul = MPMUL_INIT;
 199
 200     /* --- Multiply a bunch of primes together --- */
 201
 202     if (qq) {
 203       mpmul_add(&mmul, qq);
 204     }
 205     for (i = 0; i < mm; i++) {
 206       if (!c[i])
 207         continue;
 208       if (!v[i]) {
 209         mp *z;
 210
 211         DRESET(&dn);
 212         dstr_putf(&dn, "%s [%lu] = ", name, seq++);
 213         z = mprand(newp, ql, ev.r, 1);
 214         z = pgen(dn.buf, z, z, iev, iec,
 215                  0, pgen_filter, &pf, rabin_iters(ql), pgen_test, &rb);
 216         v[i] = z;
 217       }
 218       mpmul_add(&mmul, v[i]);
 219     }
 220
 221     /* --- Now do some testing --- */
 222
 223     {
 224       mp *p = mpmul_done(&mmul);
 225       mp *g = newp;
 226       int rc;
 227
 228       /* --- Check for small factors --- */
 229
 230       p = mp_lsl(p, p, 1);
 231       p = mp_add(p, p, MP_ONE);
 232       mp_gcd(&g, 0, 0, p, primorial);
 233       if (MP_CMP(g, !=, MP_ONE)) {
 234         mp_drop(g);
 235         mp_drop(p);
 236         continue;
 237       }
 238       mp_drop(g);
 239
 240       /* --- Send an event out --- */
 241
 242       ev.m = p;
 243       if (oev && oev(PGEN_TRY, &ev, oec) == PGEN_ABORT) {
 244         mp_drop(p);
 245         goto fail;
 246       }
 247
 248       /* --- Do the Rabin testing --- */
 249
 250       rabin_create(&rb, p);
 251       g = MP_NEW;
 252       do {
 253         g = mprand_range(g, p, ev.r, 1);
 254         rc = rabin_test(&rb, g);
 255         if (rc == PGEN_PASS) {
 256           ev.tests--;
 257           if (!ev.tests)
 258             rc = PGEN_DONE;
 259         }
 260         if (oev &&oev(rc, &ev, oec) == PGEN_ABORT)
 261           rc = PGEN_ABORT;
 262       } while (rc == PGEN_PASS);
 263
 264       rabin_destroy(&rb);
 265       mp_drop(g);
 266       if (rc == PGEN_DONE)
 267         d = p;
 268       else
 269         mp_drop(p);
 270       if (rc == PGEN_ABORT)
 271         goto fail;
 272       if (rc == PGEN_DONE)
 273         goto done;
 274       ev.tests = ntest;
 275       ev.m = 0;
 276     }
 277   } while (comb_next(c, mm, nn));
 278
 279   /* --- That failed --- */
 280
 281   if (ev.steps) {
 282     ev.steps--;
 283     if (!ev.steps) {
 284       if (oev)
 285         oev(PGEN_ABORT, &ev, &oec);
 286       goto fail;
 287     }
 288   }
 289
 290   for (i = 0; i < mm; i++)
 291     mp_drop(v[i]);
 292   goto again;
 293
 294   /* --- We did it! --- */
 295
 296 done: {
 297     mp **vv = 0;
 298     if (f) {
 299       if (qq)
 300         nn++;
 301       *nf = nn;
 302       *f = vv = xmalloc(nn * sizeof(mp *));
 303     }
 304
 305     for (i = 0; i < mm; i++) {
 306       if (c[i] && vv)
 307         *vv++ = v[i];
 308       else if (v[i])
 309         mp_drop(v[i]);
 310     }
 311     if (qq) {
 312       if (vv)
 313         *vv++ = qq;
 314       else
 315         mp_drop(qq);
 316     }
 317     xfree(v);
 318     xfree(c);
 319     dstr_destroy(&dn);
 320     return (d);
 321   }
 322
 323   /* --- We blew it --- */
 324
 325 fail:
 326   for (i = 0; i < mm; i++)
 327     mp_drop(v[i]);
 328   if (qq)
 329     mp_drop(qq);
 330   xfree(v);
 331   xfree(c);
 332   dstr_destroy(&dn);
 333   return (0);
 334 }
 335
 336 /*----- That's all, folks -------------------------------------------------*/