chiark - git - mdw - catacomb/blob - limlee.c

   1 /* -*-c-*-
   2  *
   3  * $Id: limlee.c,v 1.9 2004/04/08 01:36:15 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 /*----- Header files ------------------------------------------------------*/
  31
  32 #include <mLib/alloc.h>
  33 #include <mLib/dstr.h>
  34
  35 #include "limlee.h"
  36 #include "mpmul.h"
  37 #include "mprand.h"
  38 #include "pgen.h"
  39 #include "rabin.h"
  40
  41 /*----- Stepping through combinations -------------------------------------*/
  42
  43 /* --- @comb_init@ --- *
  44  *
  45  * Arguments:   @octet *c@ = pointer to byte-flag array
  46  *              @unsigned n@ = number of items in the array
  47  *              @unsigned r@ = number of desired items
  48  *
  49  * Returns:     ---
  50  *
  51  * Use:         Initializes a byte-flag array which, under the control of
  52  *              @comb_next@, will step through all combinations of @r@ chosen
  53  *              elements.
  54  */
  55
  56 static void comb_init(octet *c, unsigned n, unsigned r)
  57 {
  58   memset(c, 0, n - r);
  59   memset(c + (n - r), 1, r);
  60 }
  61
  62 /* --- @comb_next@ --- *
  63  *
  64  * Arguments:   @octet *c@ = pointer to byte-flag array
  65  *              @unsigned n@ = number of items in the array
  66  *              @unsigned r@ = number of desired items
  67  *
  68  * Returns:     Nonzero if another combination was returned, zero if we've
  69  *              reached the end.
  70  *
  71  * Use:         Steps on to the next combination in sequence.
  72  */
  73
  74 static int comb_next(octet *c, unsigned n, unsigned r)
  75 {
  76   unsigned g = 0;
  77
  78   /* --- How the algorithm works --- *
  79    *
  80    * Set bits start at the end and work their way towards the start.
  81    * Excepting bits already at the start, we scan for the lowest set bit, and
  82    * move it one place nearer the start.  A group of bits at the start are
  83    * counted and reset just below the `moved' bit.  If there is no moved bit
  84    * then we're done.
  85    */
  86
  87   /* --- Count the group at the start --- */
  88
  89   for (; *c; c++) {
  90     g++;
  91     *c = 0;
  92   }
  93   if (g == r)
  94     return (0);
  95
  96   /* --- Move the next bit down one --- *
  97    *
  98    * There must be one, because otherwise we'd have counted %$r$% bits
  99    * earlier.
 100    */
 101
 102   for (; !*c; c++)
 103     ;
 104   *c = 0;
 105   g++;
 106   for (; g; g--)
 107     *--c = 1;
 108   return (1);
 109 }
 110
 111 /*----- Default prime generator -------------------------------------------*/
 112
 113 static void llgen(limlee_factor *f, unsigned pl, limlee_stepctx *l)
 114 {
 115   pgen_filterctx pf;
 116   rabin r;
 117   mp *p;
 118
 119 again:
 120   p = mprand(l->newp, pl, l->r, 1);
 121   pf.step = 2;
 122   p = pgen(l->d.buf, p, p, l->iev, l->iec, 0, pgen_filter, &pf,
 123            rabin_iters(pl), pgen_test, &r);
 124   if (!p)
 125     goto again;
 126   f->p = p;
 127 }
 128
 129 static void llfree(limlee_factor *f, limlee_stepctx *l)
 130 {
 131   mp_drop(f->p);
 132 }
 133
 134 static const limlee_primeops primeops_simple = { llgen, llfree };
 135
 136 /*----- Lim-Lee stepper ---------------------------------------------------*/
 137
 138 /* --- @init@ --- *
 139  *
 140  * Arguments:   @pgen_event *ev@ = pointer to event block
 141  *              @limlee_stepctx *l@ = pointer to Lim-Lee context
 142  *
 143  * Returns:     A @PGEN@ result code.
 144  *
 145  * Use:         Initializes the stepper.
 146  */
 147
 148 static int init(pgen_event *ev, limlee_stepctx *l)
 149 {
 150   size_t i;
 151   unsigned qql;
 152
 153   /* --- First of all, decide on a number of factors to make --- */
 154
 155   l->nf = l->pl / l->ql;
 156   qql = l->pl % l->ql;
 157   if (!l->nf)
 158     return (PGEN_ABORT);
 159   else if (qql && l->nf > 1) {
 160     l->nf--;
 161     qql += l->ql;
 162   }
 163
 164   /* --- Now decide on how many primes I'll actually generate --- *
 165    *
 166    * The formula %$m = \max(3 n + 5, 25)$% comes from GPG's prime generation
 167    * library.
 168    */
 169
 170   l->poolsz = l->nf * 3 + 5;
 171   if (l->poolsz < 25)
 172     l->poolsz = 25;
 173
 174   /* --- Allocate and initialize the various tables --- */
 175
 176   l->c = xmalloc(l->poolsz);
 177   l->v = xmalloc(l->poolsz * sizeof(limlee_factor));
 178   comb_init(l->c, l->poolsz, l->nf);
 179   for (i = 0; i < l->poolsz; i++)
 180     l->v[i].p = 0;
 181
 182   /* --- Other bits of initialization --- */
 183
 184   l->seq = 0;
 185   dstr_create(&l->d);
 186   if (!l->pops) {
 187     l->pops = &primeops_simple;
 188     l->pc = 0;
 189   }
 190
 191   /* --- Find a big prime --- */
 192
 193   if (!qql)
 194     l->qq.p = 0;
 195   else {
 196     dstr_putf(&l->d, "%s*", ev->name);
 197     l->pops->pgen(&l->qq, qql, l);
 198   }
 199
 200   return (PGEN_TRY);
 201 }
 202
 203 /* --- @next@ --- *
 204  *
 205  * Arguments:   @int rq@ = request which triggered this call
 206  *              @pgen_event *ev@ = pointer to event block
 207  *              @limlee_stepctx *l@ = pointer to Lim-Lee context
 208  *
 209  * Returns:     A @PGEN@ result code.
 210  *
 211  * Use:         Initializes the stepper.
 212  */
 213
 214 static int next(int rq, pgen_event *ev, limlee_stepctx *l)
 215 {
 216   mp *p;
 217   int rc;
 218
 219   mp_drop(ev->m);
 220
 221   for (;;) {
 222     size_t i;
 223     mpmul mm = MPMUL_INIT;
 224
 225     /* --- Step on to next combination --- */
 226
 227     if (rq == PGEN_TRY && !comb_next(l->c, l->poolsz, l->nf)) {
 228       for (i = 0; i < l->poolsz; i++) {
 229         l->pops->pfree(&l->v[i], l);
 230         l->v[i].p = 0;
 231       }
 232     }
 233     rq = PGEN_TRY; /* For next time through */
 234
 235     /* --- Gather up some factors --- */
 236
 237     if (l->qq.p)
 238       mpmul_add(&mm, l->qq.p);
 239     for (i = 0; i < l->poolsz; i++) {
 240       if (!l->c[i])
 241         continue;
 242       if (!l->v[i].p) {
 243         DRESET(&l->d);
 244         dstr_putf(&l->d, "%s_%lu", ev->name, l->seq++);
 245         l->pops->pgen(&l->v[i], l->ql, l);
 246       }
 247       mpmul_add(&mm, l->v[i].p);
 248     }
 249
 250     /* --- Check it for small factors --- */
 251
 252     p = mpmul_done(&mm);
 253     p = mp_lsl(p, p, 1);
 254     p->v[0] |= 1;
 255     if ((rc = pfilt_smallfactor(p)) != PGEN_FAIL)
 256       break;
 257     mp_drop(p);
 258   }
 259
 260   ev->m = p;
 261   return (rc);
 262 }
 263
 264 /* --- @done@ --- *
 265  *
 266  * Arguments:   @pgen_event *ev@ = pointer to event block
 267  *              @limlee_stepctx *l@ = pointer to Lim-Lee context
 268  *
 269  * Returns:     A @PGEN@ result code.
 270  *
 271  * Use:         Finalizes the stepper.  The output values in the context
 272  *              take on their final results; other resources are discarded.
 273  */
 274
 275 static int done(pgen_event *ev, limlee_stepctx *l)
 276 {
 277   size_t i, j;
 278   limlee_factor *v;
 279
 280   /* --- If an output vector of factors is wanted, produce one --- */
 281
 282   if (!(l->f & LIMLEE_KEEPFACTORS))
 283     v = 0;
 284   else {
 285     if (l->qq.p)
 286       l->nf++;
 287     v = xmalloc(l->nf * sizeof(limlee_factor));
 288   }
 289
 290   for (i = 0, j = 0; i < l->poolsz; i++) {
 291     if (v && l->c[i])
 292       v[j++] = l->v[i];
 293     else if (l->v[i].p)
 294       l->pops->pfree(&l->v[i], l);
 295   }
 296
 297   if (l->qq.p) {
 298     if (v)
 299       v[j++] = l->qq;
 300     else
 301       l->pops->pfree(&l->qq, l);
 302   }
 303
 304   xfree(l->v);
 305   l->v = v;
 306
 307   /* --- Free other resources --- */
 308
 309   xfree(l->c);
 310   dstr_destroy(&l->d);
 311
 312   /* --- Done --- */
 313
 314   return (PGEN_DONE);
 315 }
 316
 317 /* --- @limlee_step@ --- */
 318
 319 int limlee_step(int rq, pgen_event *ev, void *p)
 320 {
 321   limlee_stepctx *l = p;
 322   int rc;
 323
 324   switch (rq) {
 325     case PGEN_BEGIN:
 326       if ((rc = init(ev, l)) != PGEN_TRY)
 327         return (rc);
 328     case PGEN_TRY:
 329       return (next(rq, ev, l));
 330     case PGEN_DONE:
 331       return (done(ev, l));
 332   }
 333   return (PGEN_ABORT);
 334 }
 335
 336 /*----- Main code ---------------------------------------------------------*/
 337
 338 /* --- @limlee@ --- *
 339  *
 340  * Arguments:   @const char *name@ = pointer to name root
 341  *              @mp *d@ = pointer to destination integer
 342  *              @mp *newp@ = how to generate factor primes
 343  *              @unsigned ql@ = size of individual factors
 344  *              @unsigned pl@ = size of large prime
 345  *              @grand *r@ = a random number source
 346  *              @unsigned on@ = number of outer attempts to make
 347  *              @pgen_proc *oev@ = outer event handler function
 348  *              @void *oec@ = argument for the outer event handler
 349  *              @pgen_proc *iev@ = inner event handler function
 350  *              @void *iec@ = argument for the inner event handler
 351  *              @size_t *nf@, @mp ***f@ = output array for factors
 352  *
 353  * Returns:     A Lim-Lee prime, or null if generation failed.
 354  *
 355  * Use:         Generates Lim-Lee primes.  A Lim-Lee prime %$p$% is one which
 356  *              satisfies %$p = 2 \prod_i q_i + 1$%, where all of the %$q_i$%
 357  *              are large enough to resist square-root discrete log
 358  *              algorithms.
 359  *
 360  *              If we succeed, and @f@ is non-null, we write the array of
 361  *              factors chosen to @f@ for the benefit of the caller.
 362  */
 363
 364 mp *limlee(const char *name, mp *d, mp *newp,
 365            unsigned ql, unsigned pl, grand *r,
 366            unsigned on, pgen_proc *oev, void *oec,
 367            pgen_proc *iev, void *iec,
 368            size_t *nf, mp ***f)
 369 {
 370   limlee_stepctx l;
 371   rabin rr;
 372
 373   l.f = 0; if (f) l.f |= LIMLEE_KEEPFACTORS;
 374   l.newp = newp;
 375   l.pl = pl; l.ql = ql;
 376   l.pops = 0;
 377   l.iev = iev;
 378   l.iec = iec;
 379   l.r = r;
 380
 381   d = pgen(name, d, 0, oev, oec, on, limlee_step, &l,
 382            rabin_iters(pl), pgen_test, &rr);
 383
 384   if (d && f) {
 385     mp **v;
 386     size_t i;
 387     v = xmalloc(l.nf * sizeof(mp *));
 388     for (i = 0; i < l.nf; i++)
 389       v[i] = l.v[i].p;
 390     xfree(l.v);
 391     *f = v;
 392     *nf = l.nf;
 393   }
 394
 395   return (d);
 396 }
 397
 398 /*----- That's all, folks -------------------------------------------------*/