Lim-Lee prime generation.

diff --git a/limlee.c b/limlee.c
new file mode 100644 (file)
index 0000000..a25c034
--- /dev/null
+++ b/limlee.c
@@ -0,0 +1,332 @@
+/* -*-c-*-
+ *
+ * $Id: limlee.c,v 1.1 2000/07/09 21:30:58 mdw Exp $
+ *
+ * Generate Lim-Lee primes
+ *
+ * (c) 2000 Straylight/Edgeware
+ */
+
+/*----- Licensing notice --------------------------------------------------* 
+ *
+ * This file is part of Catacomb.
+ *
+ * Catacomb is free software; you can redistribute it and/or modify
+ * it under the terms of the GNU Library General Public License as
+ * published by the Free Software Foundation; either version 2 of the
+ * License, or (at your option) any later version.
+ * 
+ * Catacomb is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+ * GNU Library General Public License for more details.
+ * 
+ * You should have received a copy of the GNU Library General Public
+ * License along with Catacomb; if not, write to the Free
+ * Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
+ * MA 02111-1307, USA.
+ */
+
+/*----- Revision history --------------------------------------------------* 
+ *
+ * $Log: limlee.c,v $
+ * Revision 1.1  2000/07/09 21:30:58  mdw
+ * Lim-Lee prime generation.
+ *
+ */
+
+/*----- Header files ------------------------------------------------------*/
+
+#include <mLib/alloc.h>
+#include <mLib/dstr.h>
+
+#include "limlee.h"
+#include "mpmul.h"
+#include "mprand.h"
+#include "pgen.h"
+#include "primorial.h"
+#include "rabin.h"
+
+/*----- Main code ---------------------------------------------------------*/
+
+/* --- @limlee@ --- *
+ *
+ * Arguments:  @const char *name@ = pointer to name root
+ *             @mp *d@ = pointer to destination integer
+ *             @mp *newp@ = how to generate factor primes
+ *             @unsigned ql@ = size of individual factors
+ *             @unsigned pl@ = size of large prime
+ *             @grand *r@ = a random number source
+ *             @unsigned on@ = number of outer attempts to make
+ *             @pgen_proc *oev@ = outer event handler function
+ *             @void *oec@ = argument for the outer event handler
+ *             @pgen_proc *iev@ = inner event handler function
+ *             @void *iec@ = argument for the inner event handler
+ *             @size_t *nf@, @mp ***f@ = output array for factors
+ *
+ * Returns:    A Lim-Lee prime, or null if generation failed.
+ *
+ * Use:                Generates Lim-Lee primes.  A Lim-Lee prime %$p$% is one which
+ *             satisfies %$p = 2 \prod_i q_i + 1$%, where all of the %$q_i$%
+ *             are large enough to resist square-root discrete log
+ *             algorithms.
+ *
+ *             If we succeed, and @f@ is non-null, we write the array of
+ *             factors chosen to @f@ for the benefit of the caller.
+ */
+
+static void comb_init(octet *c, unsigned n, unsigned r)
+{
+  memset(c, 0, n - r);
+  memset(c + (n - r), 1, r);
+}
+
+static int comb_next(octet *c, unsigned n, unsigned r)
+{
+  unsigned g = 0;
+
+  /* --- How the algorithm works --- *
+   *
+   * Set bits start at the end and work their way towards the start.
+   * Excepting bits already at the start, we scan for the lowest set bit, and
+   * move it one place nearer the start.  A group of bits at the start are
+   * counted and reset just below the `moved' bit.  If there is no moved bit
+   * then we're done.
+   */
+
+  /* --- Count the group at the start --- */
+
+  for (; *c; c++) {
+    g++;
+    *c = 0;
+  }
+  if (g == r)
+    return (0);
+
+  /* --- Move the next bit down one --- *
+   *
+   * There must be one, because otherwise we'd have counted %$r$% bits
+   * earlier.
+   */
+
+  for (; !*c; c++)
+    ;
+  *c = 0;
+  g++;
+  for (; g; g--)
+    *--c = 1;
+  return (1);
+}
+
+mp *limlee(const char *name, mp *d, mp *newp,
+          unsigned ql, unsigned pl, grand *r,
+          unsigned on, pgen_proc *oev, void *oec,
+          pgen_proc *iev, void *iec,
+          size_t *nf, mp ***f)
+{
+  dstr dn = DSTR_INIT;
+  unsigned qql;
+  mp *qq = 0;
+  unsigned nn;
+  unsigned mm;
+  mp **v;
+  octet *c;
+  unsigned i;
+  unsigned long seq = 0;
+  pgen_event ev;
+  unsigned ntest;
+  rabin rb;
+  pgen_filterctx pf;
+
+  /* --- First of all, decide on a number of factors to make --- */
+
+  nn = pl/ql;
+  qql = pl%ql;
+  if (!nn)
+    return (0);
+  else if (qql && nn > 1) {
+    nn--;
+    qql += ql;
+  }
+
+  /* --- Now decide on how many primes I'll actually generate --- *
+   *
+   * The formula %$m = \max(3 n + 5, 25)$% comes from GPG's prime generation
+   * library.
+   */
+
+  mm = nn * 3 + 5;
+  if (mm < 25)
+    mm = 25;
+
+  /* --- Now allocate the working memory --- */
+
+  primorial_setup();
+  v = xmalloc(mm * sizeof(mp *));
+  c = xmalloc(mm);
+
+  /* --- Initialize everything and try to find a prime --- */
+
+  ev.name = name;
+  ev.m = 0;
+  ev.steps = on;
+  ev.tests = ntest = rabin_iters(pl);
+  ev.r = r;
+
+  if (oev && oev(PGEN_BEGIN, &ev, oec) == PGEN_ABORT)
+    goto fail;
+
+  if (qql) {
+    dstr_putf(&dn, "%s [+]", name);
+    qq = mprand(d, qql, r, 1);
+    pf.step = 2;
+    qq = pgen(dn.buf, qq, qq, iev, iec,
+             0, pgen_filter, &pf, rabin_iters(qql), pgen_test, &rb);
+  }
+
+again:
+  comb_init(c, mm, nn);
+  for (i = 0; i < mm; i++)
+    v[i] = 0;
+
+  /* --- The main combinations loop --- */
+
+  do {
+    mpmul mmul = MPMUL_INIT;
+
+    /* --- Multiply a bunch of primes together --- */
+
+    if (qq) {
+      mpmul_add(&mmul, qq);
+    }
+    for (i = 0; i < mm; i++) {
+      if (!c[i])
+       continue;
+      if (!v[i]) {
+       mp *z;
+
+       DRESET(&dn);
+       dstr_putf(&dn, "%s [%lu] = ", name, seq++);
+       z = mprand(newp, ql, ev.r, 1);
+       z = pgen(dn.buf, z, z, iev, iec,
+                0, pgen_filter, &pf, rabin_iters(ql), pgen_test, &rb);
+       v[i] = z;
+      }
+      mpmul_add(&mmul, v[i]);
+    }
+
+    /* --- Now do some testing --- */
+
+    {
+      mp *p = mpmul_done(&mmul);
+      mp *g = newp;
+      int rc;
+
+      /* --- Check for small factors --- */
+
+      p = mp_lsl(p, p, 1);
+      p = mp_add(p, p, MP_ONE);
+      mp_gcd(&g, 0, 0, p, primorial);
+      if (MP_CMP(g, !=, MP_ONE)) {
+       mp_drop(g);
+       mp_drop(p);
+       continue;
+      }
+      mp_drop(g);
+
+      /* --- Send an event out --- */
+
+      ev.m = p;
+      if (oev && oev(PGEN_TRY, &ev, oec) == PGEN_ABORT) {
+       mp_drop(p);
+       goto fail;
+      }
+
+      /* --- Do the Rabin testing --- */
+
+      rabin_create(&rb, p);
+      g = MP_NEW;
+      do {
+       g = mprand_range(g, p, ev.r, 1);
+       rc = rabin_test(&rb, g);
+       if (rc == PGEN_PASS) {
+         ev.tests--;
+         if (!ev.tests)
+           rc = PGEN_DONE;
+       }
+       if (oev &&oev(rc, &ev, oec) == PGEN_ABORT)
+         rc = PGEN_ABORT;
+      } while (rc == PGEN_PASS);
+
+      rabin_destroy(&rb);
+      mp_drop(g);
+      if (rc == PGEN_DONE)
+       d = p;
+      else
+       mp_drop(p);
+      if (rc == PGEN_ABORT)
+       goto fail;
+      if (rc == PGEN_DONE)
+       goto done;
+      ev.tests = ntest;
+      ev.m = 0;
+    }
+  } while (comb_next(c, mm, nn));
+
+  /* --- That failed --- */
+
+  if (ev.steps) {
+    ev.steps--;
+    if (!ev.steps) {
+      if (oev)
+       oev(PGEN_ABORT, &ev, &oec);
+      goto fail;
+    }
+  }
+
+  for (i = 0; i < mm; i++)
+    mp_drop(v[i]);
+  goto again;
+
+  /* --- We did it! --- */
+
+done: {
+    mp **vv = 0;
+    if (f) {
+      if (qq)
+       nn++;
+      *nf = nn;
+      *f = vv = xmalloc(nn * sizeof(mp *));
+    }
+  
+    for (i = 0; i < mm; i++) {
+      if (c[i] && vv)
+       *vv++ = v[i];
+      else if (v[i])
+       mp_drop(v[i]);
+    }
+    if (qq) {
+      if (vv)
+       *vv++ = qq;
+      else
+       mp_drop(qq);
+    }
+    xfree(v);
+    xfree(c);
+    dstr_destroy(&dn);
+    return (d);
+  }
+
+  /* --- We blew it --- */
+
+fail:
+  for (i = 0; i < mm; i++)
+    mp_drop(v[i]);
+  if (qq)
+    mp_drop(qq);
+  xfree(v);
+  xfree(c);
+  dstr_destroy(&dn);
+  return (0);
+}
+

diff --git a/limlee.h b/limlee.h
new file mode 100644 (file)
index 0000000..d6fcc4b
--- /dev/null
+++ b/limlee.h
@@ -0,0 +1,99 @@
+/* -*-c-*-
+ *
+ * $Id: limlee.h,v 1.1 2000/07/09 21:30:58 mdw Exp $
+ *
+ * Generate Lim-Lee primes
+ *
+ * (c) 2000 Straylight/Edgeware
+ */
+
+/*----- Licensing notice --------------------------------------------------* 
+ *
+ * This file is part of Catacomb.
+ *
+ * Catacomb is free software; you can redistribute it and/or modify
+ * it under the terms of the GNU Library General Public License as
+ * published by the Free Software Foundation; either version 2 of the
+ * License, or (at your option) any later version.
+ * 
+ * Catacomb is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+ * GNU Library General Public License for more details.
+ * 
+ * You should have received a copy of the GNU Library General Public
+ * License along with Catacomb; if not, write to the Free
+ * Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
+ * MA 02111-1307, USA.
+ */
+
+/*----- Revision history --------------------------------------------------* 
+ *
+ * $Log: limlee.h,v $
+ * Revision 1.1  2000/07/09 21:30:58  mdw
+ * Lim-Lee prime generation.
+ *
+ */
+
+#ifndef CATACOMB_LIMLEE_H
+#define CATACOMB_LIMLEE_H
+
+#ifdef __cplusplus
+  extern "C" {
+#endif
+
+/*----- Header files ------------------------------------------------------*/
+
+#ifndef CATACOMB_GRAND_H
+#  include "grand.h"
+#endif
+
+#ifndef CATACOMB_MP_H
+#  include "mp.h"
+#endif
+
+#ifndef CATACOMB_PGEN_H
+#  include "pgen.h"
+#endif
+
+/*----- Functions provided ------------------------------------------------*/
+
+/* --- @limlee@ --- *
+ *
+ * Arguments:  @const char *name@ = pointer to name root
+ *             @mp *d@ = pointer to destination integer
+ *             @mp *newp@ = how to generate factor primes
+ *             @unsigned ql@ = size of individual factors
+ *             @unsigned pl@ = size of large prime
+ *             @grand *r@ = a random number source
+ *             @unsigned on@ = number of outer attempts to make
+ *             @pgen_proc *oev@ = outer event handler function
+ *             @void *oec@ = argument for the outer event handler
+ *             @pgen_proc *iev@ = inner event handler function
+ *             @void *iec@ = argument for the inner event handler
+ *             @size_t *nf@, @mp ***f@ = output array for factors
+ *
+ * Returns:    A Lim-Lee prime, or null if generation failed.
+ *
+ * Use:                Generates Lim-Lee primes.  A Lim-Lee prime %$p$% is one which
+ *             satisfies %$p = 2 \prod_i q_i + 1$%, where all of the %$q_i$%
+ *             are large enough to resist square-root discrete log
+ *             algorithms.
+ *
+ *             If we succeed, and @f@ is non-null, we write the array of
+ *             factors chosen to @f@ for the benefit of the caller.
+ */
+
+extern mp *limlee(const char */*name*/, mp */*d*/, mp */*newp*/,
+                 unsigned /*ql*/, unsigned /*pl*/, grand */*r*/,
+                 unsigned /*on*/, pgen_proc */*oev*/, void */*oec*/,
+                 pgen_proc */*iev*/, void */*iec*/,
+                 size_t */*nf*/, mp ***/*f*/);
+
+/*----- That's all, folks -------------------------------------------------*/
+
+#ifdef __cplusplus
+  }
+#endif
+
+#endif