[u/mdw/catacomb] / limlee.c

/* -*-c-*-
 *
 * $Id: limlee.c,v 1.6 2001/01/25 21:16:20 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.6  2001/01/25 21:16:20  mdw
 * Boring cosmetic stuff.
 *
 * Revision 1.5  2000/08/18 19:16:51  mdw
 * New stepper interface for constructing Lim-Lee primes.
 *
 * Revision 1.4  2000/08/15 21:45:05  mdw
 * Use the new trial division equipment in pfilt.  This gives a 10%
 * performance improvement in dsa-gen.t.
 *
 * Revision 1.3  2000/07/29 09:58:32  mdw
 * (limlee): Bug fix.  Old versions didn't set the filter step if @ql@ was
 * an exact divisor of @pl@.
 *
 * Revision 1.2  2000/07/26 18:00:00  mdw
 * No footer line!
 *
 * 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 "rabin.h"

/*----- Stepping through combinations -------------------------------------*/

/* --- @comb_init@ --- *
 *
 * Arguments:	@octet *c@ = pointer to byte-flag array
 *		@unsigned n@ = number of items in the array
 *		@unsigned r@ = number of desired items
 *
 * Returns:	---
 *
 * Use:		Initializes a byte-flag array which, under the control of
 *		@comb_next@, will step through all combinations of @r@ chosen
 *		elements.
 */

static void comb_init(octet *c, unsigned n, unsigned r)
{
  memset(c, 0, n - r);
  memset(c + (n - r), 1, r);
}

/* --- @comb_next@ --- *
 *
 * Arguments:	@octet *c@ = pointer to byte-flag array
 *		@unsigned n@ = number of items in the array
 *		@unsigned r@ = number of desired items
 *
 * Returns:	Nonzero if another combination was returned, zero if we've
 *		reached the end.
 *
 * Use:		Steps on to the next combination in sequence.
 */

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);
}

/*----- Default prime generator -------------------------------------------*/

static void llgen(limlee_factor *f, unsigned pl, limlee_stepctx *l)
{
  pgen_filterctx pf;
  rabin r;
  mp *p;

again:
  p = mprand(l->newp, pl, l->r, 1);
  pf.step = 2;
  p = pgen(l->d.buf, p, p, l->iev, l->iec, 0, pgen_filter, &pf,
	   rabin_iters(pl), pgen_test, &r);
  if (!p)
    goto again;
  f->p = p;
}

static void llfree(limlee_factor *f, limlee_stepctx *l)
{
  if (f->p)
    mp_drop(f->p);
}

static const limlee_primeops primeops_simple = { llgen, llfree };

/*----- Lim-Lee stepper ---------------------------------------------------*/

/* --- @init@ --- *
 *
 * Arguments:	@pgen_event *ev@ = pointer to event block
 *		@limlee_stepctx *l@ = pointer to Lim-Lee context
 *
 * Returns:	A @PGEN@ result code.
 *
 * Use:		Initializes the stepper.
 */

static int init(pgen_event *ev, limlee_stepctx *l)
{
  size_t i;
  unsigned qql;

  /* --- First of all, decide on a number of factors to make --- */

  l->nf = l->pl / l->ql;
  qql = l->pl % l->ql;
  if (!l->nf)
    return (PGEN_ABORT);
  else if (qql && l->nf > 1) {
    l->nf--;
    qql += l->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.
   */

  l->poolsz = l->nf * 3 + 5;
  if (l->poolsz < 25)
    l->poolsz = 25;

  /* --- Allocate and initialize the various tables --- */

  l->c = xmalloc(l->poolsz);
  l->v = xmalloc(l->poolsz * sizeof(limlee_factor));
  comb_init(l->c, l->poolsz, l->nf);
  for (i = 0; i < l->poolsz; i++)
    l->v[i].p = 0;

  /* --- Other bits of initialization --- */

  l->seq = 0;
  l->r = ev->r;
  dstr_create(&l->d);
  if (!l->pops) {
    l->pops = &primeops_simple;
    l->pc = 0;
  }

  /* --- Find a big prime --- */

  if (!qql)
    l->qq.p = 0;
  else {
    dstr_putf(&l->d, "%s*", ev->name);
    l->pops->pgen(&l->qq, qql, l);
  }       

  return (PGEN_TRY);
}

/* --- @next@ --- *
 *
 * Arguments:	@int rq@ = request which triggered this call
 *		@pgen_event *ev@ = pointer to event block
 *		@limlee_stepctx *l@ = pointer to Lim-Lee context
 *
 * Returns:	A @PGEN@ result code.
 *
 * Use:		Initializes the stepper.
 */

static int next(int rq, pgen_event *ev, limlee_stepctx *l)
{
  mp *p;
  int rc;

  if (ev->m)
    mp_drop(ev->m);
  l->r = ev->r;

  for (;;) {
    size_t i;
    mpmul mm = MPMUL_INIT;

    /* --- Step on to next combination --- */

    if (rq == PGEN_TRY && !comb_next(l->c, l->poolsz, l->nf)) {
      for (i = 0; i < l->poolsz; i++) {
	l->pops->pfree(&l->v[i], l);
	l->v[i].p = 0;
      }
    }
    rq = PGEN_TRY; /* For next time through */

    /* --- Gather up some factors --- */

    if (l->qq.p)
      mpmul_add(&mm, l->qq.p);
    for (i = 0; i < l->poolsz; i++) {
      if (!l->c[i])
	continue;
      if (!l->v[i].p) {
	DRESET(&l->d);
	dstr_putf(&l->d, "%s_%lu", ev->name, l->seq++);
	l->pops->pgen(&l->v[i], l->ql, l);
      }
      mpmul_add(&mm, l->v[i].p);
    }

    /* --- Check it for small factors --- */

    p = mpmul_done(&mm);
    p = mp_lsl(p, p, 1);
    p->v[0] |= 1;
    if ((rc = pfilt_smallfactor(p)) != PGEN_FAIL)
      break;
    mp_drop(p);
  }

  ev->m = p;
  return (rc);
}

/* --- @done@ --- *
 *
 * Arguments:	@pgen_event *ev@ = pointer to event block
 *		@limlee_stepctx *l@ = pointer to Lim-Lee context
 *
 * Returns:	A @PGEN@ result code.
 *
 * Use:		Finalizes the stepper.  The output values in the context
 *		take on their final results; other resources are discarded.
 */

static int done(pgen_event *ev, limlee_stepctx *l)
{
  size_t i, j;
  limlee_factor *v;

  /* --- If an output vector of factors is wanted, produce one --- */

  if (!(l->f & LIMLEE_KEEPFACTORS))
    v = 0;
  else {
    if (l->qq.p)
      l->nf++;
    v = xmalloc(l->nf * sizeof(limlee_factor));
  }

  for (i = 0, j = 0; i < l->poolsz; i++) {
    if (v && l->c[i])
      v[j++] = l->v[i];
    else if (l->v[i].p)
      l->pops->pfree(&l->v[i], l);
  }

  if (l->qq.p) {
    if (v)
      v[j++] = l->qq;
    else
      l->pops->pfree(&l->qq, l);
  }

  xfree(l->v);
  l->v = v;

  /* --- Free other resources --- */

  xfree(l->c);
  dstr_destroy(&l->d);

  /* --- Done --- */

  return (PGEN_DONE);
}

/* --- @limlee_step@ --- */

int limlee_step(int rq, pgen_event *ev, void *p)
{
  limlee_stepctx *l = p;
  int rc;

  switch (rq) {
    case PGEN_BEGIN:
      if ((rc = init(ev, l)) != PGEN_TRY)
	return (rc);
    case PGEN_TRY:
      return (next(rq, ev, l));
    case PGEN_DONE:
      return (done(ev, l));
  }
  return (PGEN_ABORT);
}      

/*----- 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.
 */

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)
{
#ifdef notdef
  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 --- */

  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;

  pf.step = 2;
  if (qql) {
    dstr_putf(&dn, "%s*", name);
    qq = mprand(d, qql, r, 1);
    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;
      int rc;

      /* --- Check for small factors --- */

      p = mp_lsl(p, p, 1);
      p = mp_add(p, p, MP_ONE);
      rc = pfilt_smallfactor(p);
      if (rc == PGEN_FAIL) {
	mp_drop(p);
	continue;
      }

      /* --- 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);
#else
  limlee_stepctx l;
  rabin rr;

  l.f = 0; if (f) l.f |= LIMLEE_KEEPFACTORS;
  l.newp = newp;
  l.pl = pl; l.ql = ql;
  l.pops = 0;
  l.iev = iev;
  l.iec = iec;

  d = pgen(name, d, 0, oev, oec, on, limlee_step, &l,
	   rabin_iters(pl), pgen_test, &rr);

  if (f) {
    mp **v;
    size_t i;
    v = xmalloc(l.nf * sizeof(mp *));
    for (i = 0; i < l.nf; i++)
      v[i] = l.v[i].p;
    xfree(l.v);
    *f = v;
    *nf = l.nf;
  }

  return (d);
#endif
}

/*----- That's all, folks -------------------------------------------------*/
Commit	Line	Data
04361334	1	/* --c--
04361334	2	*
3c0e5a02	3	* $Id: limlee.c,v 1.6 2001/01/25 21:16:20 mdw Exp $
04361334	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 $
3c0e5a02	33	* Revision 1.6 2001/01/25 21:16:20 mdw
	34	* Boring cosmetic stuff.
	35	*
10217a5c	36	* Revision 1.5 2000/08/18 19:16:51 mdw
	37	* New stepper interface for constructing Lim-Lee primes.
	38	*
8a33545f	39	* Revision 1.4 2000/08/15 21:45:05 mdw
	40	* Use the new trial division equipment in pfilt. This gives a 10%
	41	* performance improvement in dsa-gen.t.
	42	*
2f627546	43	* Revision 1.3 2000/07/29 09:58:32 mdw
	44	* (limlee): Bug fix. Old versions didn't set the filter step if @ql@ was
	45	* an exact divisor of @pl@.
	46	*
d28d625a	47	* Revision 1.2 2000/07/26 18:00:00 mdw
	48	* No footer line!
	49	*
04361334	50	* Revision 1.1 2000/07/09 21:30:58 mdw
	51	* Lim-Lee prime generation.
	52	*
	53	*/
	54
	55	/----- Header files ------------------------------------------------------/
	56
	57	#include <mLib/alloc.h>
	58	#include <mLib/dstr.h>
	59
	60	#include "limlee.h"
	61	#include "mpmul.h"
	62	#include "mprand.h"
	63	#include "pgen.h"
04361334	64	#include "rabin.h"
04361334	65
10217a5c	66	/----- Stepping through combinations -------------------------------------/
04361334	67
10217a5c	68	/* --- @comb_init@ --- *
04361334	69	*
10217a5c	70	* Arguments: @octet *c@ = pointer to byte-flag array
	71	* @unsigned n@ = number of items in the array
	72	* @unsigned r@ = number of desired items
04361334	73	*
10217a5c	74	* Returns: ---
04361334	75	*
10217a5c	76	* Use: Initializes a byte-flag array which, under the control of
	77	* @comb_next@, will step through all combinations of @r@ chosen
	78	* elements.
04361334	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
10217a5c	87	/* --- @comb_next@ --- *
	88	*
	89	* Arguments: @octet *c@ = pointer to byte-flag array
	90	* @unsigned n@ = number of items in the array
	91	* @unsigned r@ = number of desired items
	92	*
	93	* Returns: Nonzero if another combination was returned, zero if we've
	94	* reached the end.
	95	*
	96	* Use: Steps on to the next combination in sequence.
	97	*/
	98
04361334	99	static int comb_next(octet *c, unsigned n, unsigned r)
	100	{
	101	unsigned g = 0;
	102
	103	/* --- How the algorithm works --- *
	104	*
	105	* Set bits start at the end and work their way towards the start.
	106	* Excepting bits already at the start, we scan for the lowest set bit, and
	107	* move it one place nearer the start. A group of bits at the start are
	108	* counted and reset just below the `moved' bit. If there is no moved bit
	109	* then we're done.
	110	*/
	111
	112	/* --- Count the group at the start --- */
	113
	114	for (; *c; c++) {
	115	g++;
	116	*c = 0;
	117	}
	118	if (g == r)
	119	return (0);
	120
	121	/* --- Move the next bit down one --- *
	122	*
	123	* There must be one, because otherwise we'd have counted %$r$% bits
	124	* earlier.
	125	*/
	126
	127	for (; !*c; c++)
	128	;
	129	*c = 0;
	130	g++;
	131	for (; g; g--)
	132	*--c = 1;
	133	return (1);
	134	}
	135
10217a5c	136	/----- Default prime generator -------------------------------------------/
	137
	138	static void llgen(limlee_factor f, unsigned pl, limlee_stepctx l)
	139	{
	140	pgen_filterctx pf;
	141	rabin r;
	142	mp *p;
	143
	144	again:
	145	p = mprand(l->newp, pl, l->r, 1);
	146	pf.step = 2;
	147	p = pgen(l->d.buf, p, p, l->iev, l->iec, 0, pgen_filter, &pf,
	148	rabin_iters(pl), pgen_test, &r);
	149	if (!p)
	150	goto again;
	151	f->p = p;
	152	}
	153
	154	static void llfree(limlee_factor f, limlee_stepctx l)
	155	{
	156	if (f->p)
	157	mp_drop(f->p);
	158	}
	159
	160	static const limlee_primeops primeops_simple = { llgen, llfree };
	161
	162	/----- Lim-Lee stepper ---------------------------------------------------/
	163
	164	/* --- @init@ --- *
	165	*
	166	* Arguments: @pgen_event *ev@ = pointer to event block
	167	* @limlee_stepctx *l@ = pointer to Lim-Lee context
	168	*
	169	* Returns: A @PGEN@ result code.
	170	*
	171	* Use: Initializes the stepper.
	172	*/
	173
	174	static int init(pgen_event ev, limlee_stepctx l)
	175	{
	176	size_t i;
	177	unsigned qql;
	178
	179	/* --- First of all, decide on a number of factors to make --- */
	180
	181	l->nf = l->pl / l->ql;
	182	qql = l->pl % l->ql;
	183	if (!l->nf)
	184	return (PGEN_ABORT);
	185	else if (qql && l->nf > 1) {
	186	l->nf--;
	187	qql += l->ql;
	188	}
	189
	190	/* --- Now decide on how many primes I'll actually generate --- *
	191	*
	192	* The formula %$m = \max(3 n + 5, 25)$% comes from GPG's prime generation
	193	* library.
	194	*/
	195
	196	l->poolsz = l->nf * 3 + 5;
	197	if (l->poolsz < 25)
	198	l->poolsz = 25;
	199
200	/* --- Allocate and initialize the various tables --- */
201
202	l->c = xmalloc(l->poolsz);
203	l->v = xmalloc(l->poolsz * sizeof(limlee_factor));
204	comb_init(l->c, l->poolsz, l->nf);
205	for (i = 0; i < l->poolsz; i++)
206	l->v[i].p = 0;
207
208	/* --- Other bits of initialization --- */
209
210	l->seq = 0;
211	l->r = ev->r;
212	dstr_create(&l->d);
213	if (!l->pops) {
214	l->pops = &primeops_simple;
215	l->pc = 0;
216	}
217
218	/* --- Find a big prime --- */
219
220	if (!qql)
221	l->qq.p = 0;
222	else {
223	dstr_putf(&l->d, "%s*", ev->name);
224	l->pops->pgen(&l->qq, qql, l);
225	}
226
227	return (PGEN_TRY);
228	}
229
230	/* --- @next@ --- *
231	*
232	* Arguments: @int rq@ = request which triggered this call
233	* @pgen_event *ev@ = pointer to event block
234	* @limlee_stepctx *l@ = pointer to Lim-Lee context
235	*
236	* Returns: A @PGEN@ result code.
237	*
238	* Use: Initializes the stepper.
239	*/
240
241	static int next(int rq, pgen_event ev, limlee_stepctx l)
242	{
243	mp *p;
244	int rc;
245
246	if (ev->m)
247	mp_drop(ev->m);
248	l->r = ev->r;
249
250	for (;;) {
251	size_t i;
252	mpmul mm = MPMUL_INIT;
253
254	/* --- Step on to next combination --- */
255
256	if (rq == PGEN_TRY && !comb_next(l->c, l->poolsz, l->nf)) {
257	for (i = 0; i < l->poolsz; i++) {
258	l->pops->pfree(&l->v[i], l);
259	l->v[i].p = 0;
260	}
261	}
262	rq = PGEN_TRY; /* For next time through */
263
264	/* --- Gather up some factors --- */
265
266	if (l->qq.p)
267	mpmul_add(&mm, l->qq.p);
268	for (i = 0; i < l->poolsz; i++) {
269	if (!l->c[i])
270	continue;
271	if (!l->v[i].p) {
272	DRESET(&l->d);
273	dstr_putf(&l->d, "%s_%lu", ev->name, l->seq++);
274	l->pops->pgen(&l->v[i], l->ql, l);
275	}
276	mpmul_add(&mm, l->v[i].p);
277	}
278
279	/* --- Check it for small factors --- */
280
281	p = mpmul_done(&mm);
282	p = mp_lsl(p, p, 1);
283	p->v[0] \|= 1;
284	if ((rc = pfilt_smallfactor(p)) != PGEN_FAIL)
285	break;
286	mp_drop(p);
287	}
288
289	ev->m = p;
290	return (rc);
291	}
292
293	/* --- @done@ --- *
294	*
295	* Arguments: @pgen_event *ev@ = pointer to event block
296	* @limlee_stepctx *l@ = pointer to Lim-Lee context
297	*
298	* Returns: A @PGEN@ result code.
299	*
300	* Use: Finalizes the stepper. The output values in the context
301	* take on their final results; other resources are discarded.
302	*/
303
304	static int done(pgen_event ev, limlee_stepctx l)
305	{
306	size_t i, j;
307	limlee_factor *v;
308
309	/* --- If an output vector of factors is wanted, produce one --- */
310
311	if (!(l->f & LIMLEE_KEEPFACTORS))
312	v = 0;
313	else {
314	if (l->qq.p)
315	l->nf++;
316	v = xmalloc(l->nf * sizeof(limlee_factor));
317	}
318
319	for (i = 0, j = 0; i < l->poolsz; i++) {
320	if (v && l->c[i])
321	v[j++] = l->v[i];
322	else if (l->v[i].p)
323	l->pops->pfree(&l->v[i], l);
324	}
325
326	if (l->qq.p) {
327	if (v)
328	v[j++] = l->qq;
329	else
330	l->pops->pfree(&l->qq, l);
331	}
332
333	xfree(l->v);
334	l->v = v;
335
336	/* --- Free other resources --- */
337
338	xfree(l->c);
339	dstr_destroy(&l->d);
340
341	/* --- Done --- */
342
343	return (PGEN_DONE);
344	}
345
346	/* --- @limlee_step@ --- */
347
348	int limlee_step(int rq, pgen_event ev, void p)
349	{
350	limlee_stepctx *l = p;
351	int rc;
352
353	switch (rq) {
354	case PGEN_BEGIN:
355	if ((rc = init(ev, l)) != PGEN_TRY)
356	return (rc);
357	case PGEN_TRY:
358	return (next(rq, ev, l));
359	case PGEN_DONE:
360	return (done(ev, l));
361	}
362	return (PGEN_ABORT);
363	}
364
365	/----- Main code ---------------------------------------------------------/
366
367	/* --- @limlee@ --- *
368	*
369	* Arguments: @const char *name@ = pointer to name root
370	* @mp *d@ = pointer to destination integer
371	* @mp *newp@ = how to generate factor primes
372	* @unsigned ql@ = size of individual factors
373	* @unsigned pl@ = size of large prime
374	* @grand *r@ = a random number source
375	* @unsigned on@ = number of outer attempts to make
376	* @pgen_proc *oev@ = outer event handler function
377	* @void *oec@ = argument for the outer event handler
378	* @pgen_proc *iev@ = inner event handler function
379	* @void *iec@ = argument for the inner event handler
380	* @size_t nf@, @mp **f@ = output array for factors
381	*
382	* Returns: A Lim-Lee prime, or null if generation failed.
383	*
384	* Use: Generates Lim-Lee primes. A Lim-Lee prime %$p$% is one which
385	* satisfies %$p = 2 \prod_i q_i + 1$%, where all of the %$q_i$%
386	* are large enough to resist square-root discrete log
387	* algorithms.
388	*
389	* If we succeed, and @f@ is non-null, we write the array of
390	* factors chosen to @f@ for the benefit of the caller.
391	*/
392
04361334	393	mp limlee(const char name, mp d, mp newp,
	394	unsigned ql, unsigned pl, grand *r,
	395	unsigned on, pgen_proc oev, void oec,
	396	pgen_proc iev, void iec,
	397	size_t nf, mp **f)
	398	{
10217a5c	399	#ifdef notdef
04361334	400	dstr dn = DSTR_INIT;
	401	unsigned qql;
	402	mp *qq = 0;
	403	unsigned nn;
	404	unsigned mm;
	405	mp **v;
	406	octet *c;
	407	unsigned i;
	408	unsigned long seq = 0;
	409	pgen_event ev;
	410	unsigned ntest;
	411	rabin rb;
	412	pgen_filterctx pf;
	413
	414	/* --- First of all, decide on a number of factors to make --- */
	415
	416	nn = pl/ql;
	417	qql = pl%ql;
	418	if (!nn)
	419	return (0);
	420	else if (qql && nn > 1) {
	421	nn--;
	422	qql += ql;
	423	}
	424
	425	/* --- Now decide on how many primes I'll actually generate --- *
	426	*
	427	* The formula %$m = \max(3 n + 5, 25)$% comes from GPG's prime generation
	428	* library.
	429	*/
	430
	431	mm = nn * 3 + 5;
	432	if (mm < 25)
	433	mm = 25;
	434
	435	/* --- Now allocate the working memory --- */
	436
04361334	437	v = xmalloc(mm * sizeof(mp *));
	438	c = xmalloc(mm);
	439
	440	/* --- Initialize everything and try to find a prime --- */
	441
	442	ev.name = name;
	443	ev.m = 0;
	444	ev.steps = on;
	445	ev.tests = ntest = rabin_iters(pl);
	446	ev.r = r;
	447
	448	if (oev && oev(PGEN_BEGIN, &ev, oec) == PGEN_ABORT)
	449	goto fail;
	450
2f627546	451	pf.step = 2;
04361334	452	if (qql) {
10217a5c	453	dstr_putf(&dn, "%s*", name);
04361334	454	qq = mprand(d, qql, r, 1);
04361334	455	qq = pgen(dn.buf, qq, qq, iev, iec,
	456	0, pgen_filter, &pf, rabin_iters(qql), pgen_test, &rb);
	457	}
	458
	459	again:
	460	comb_init(c, mm, nn);
	461	for (i = 0; i < mm; i++)
	462	v[i] = 0;
	463
	464	/* --- The main combinations loop --- */
	465
	466	do {
	467	mpmul mmul = MPMUL_INIT;
	468
	469	/* --- Multiply a bunch of primes together --- */
	470
3c0e5a02	471	if (qq)
04361334	472	mpmul_add(&mmul, qq);
04361334	473	for (i = 0; i < mm; i++) {
	474	if (!c[i])
	475	continue;
	476	if (!v[i]) {
	477	mp *z;
	478
	479	DRESET(&dn);
10217a5c	480	dstr_putf(&dn, "%s_%lu] = ", name, seq++);
04361334	481	z = mprand(newp, ql, ev.r, 1);
	482	z = pgen(dn.buf, z, z, iev, iec,
	483	0, pgen_filter, &pf, rabin_iters(ql), pgen_test, &rb);
	484	v[i] = z;
	485	}
	486	mpmul_add(&mmul, v[i]);
	487	}
	488
	489	/* --- Now do some testing --- */
	490
	491	{
	492	mp *p = mpmul_done(&mmul);
8a33545f	493	mp *g;
04361334	494	int rc;
	495
	496	/* --- Check for small factors --- */
	497
	498	p = mp_lsl(p, p, 1);
	499	p = mp_add(p, p, MP_ONE);
8a33545f	500	rc = pfilt_smallfactor(p);
8a33545f	501	if (rc == PGEN_FAIL) {
04361334	502	mp_drop(p);
	503	continue;
	504	}
04361334	505
	506	/* --- Send an event out --- */
	507
	508	ev.m = p;
	509	if (oev && oev(PGEN_TRY, &ev, oec) == PGEN_ABORT) {
	510	mp_drop(p);
	511	goto fail;
	512	}
	513
	514	/* --- Do the Rabin testing --- */
	515
	516	rabin_create(&rb, p);
	517	g = MP_NEW;
	518	do {
	519	g = mprand_range(g, p, ev.r, 1);
	520	rc = rabin_test(&rb, g);
	521	if (rc == PGEN_PASS) {
	522	ev.tests--;
	523	if (!ev.tests)
	524	rc = PGEN_DONE;
	525	}
3c0e5a02	526	if (oev && oev(rc, &ev, oec) == PGEN_ABORT)
04361334	527	rc = PGEN_ABORT;
	528	} while (rc == PGEN_PASS);
	529
	530	rabin_destroy(&rb);
	531	mp_drop(g);
	532	if (rc == PGEN_DONE)
	533	d = p;
	534	else
	535	mp_drop(p);
	536	if (rc == PGEN_ABORT)
	537	goto fail;
	538	if (rc == PGEN_DONE)
	539	goto done;
	540	ev.tests = ntest;
	541	ev.m = 0;
	542	}
	543	} while (comb_next(c, mm, nn));
	544
	545	/* --- That failed --- */
	546
	547	if (ev.steps) {
	548	ev.steps--;
	549	if (!ev.steps) {
	550	if (oev)
	551	oev(PGEN_ABORT, &ev, &oec);
	552	goto fail;
	553	}
	554	}
	555
	556	for (i = 0; i < mm; i++)
	557	mp_drop(v[i]);
	558	goto again;
	559
	560	/* --- We did it! --- */
	561
	562	done: {
	563	mp **vv = 0;
	564	if (f) {
	565	if (qq)
	566	nn++;
	567	*nf = nn;
	568	f = vv = xmalloc(nn sizeof(mp *));
	569	}
2f627546	570
04361334	571	for (i = 0; i < mm; i++) {
	572	if (c[i] && vv)
	573	*vv++ = v[i];
	574	else if (v[i])
	575	mp_drop(v[i]);
	576	}
	577	if (qq) {
	578	if (vv)
	579	*vv++ = qq;
	580	else
	581	mp_drop(qq);
	582	}
	583	xfree(v);
	584	xfree(c);
	585	dstr_destroy(&dn);
	586	return (d);
	587	}
	588
	589	/* --- We blew it --- */
	590
	591	fail:
	592	for (i = 0; i < mm; i++)
	593	mp_drop(v[i]);
	594	if (qq)
	595	mp_drop(qq);
	596	xfree(v);
	597	xfree(c);
	598	dstr_destroy(&dn);
	599	return (0);
10217a5c	600	#else
	601	limlee_stepctx l;
	602	rabin rr;
	603
	604	l.f = 0; if (f) l.f \|= LIMLEE_KEEPFACTORS;
	605	l.newp = newp;
	606	l.pl = pl; l.ql = ql;
	607	l.pops = 0;
	608	l.iev = iev;
	609	l.iec = iec;
	610
	611	d = pgen(name, d, 0, oev, oec, on, limlee_step, &l,
	612	rabin_iters(pl), pgen_test, &rr);
	613
	614	if (f) {
	615	mp **v;
	616	size_t i;
	617	v = xmalloc(l.nf * sizeof(mp *));
	618	for (i = 0; i < l.nf; i++)
	619	v[i] = l.v[i].p;
	620	xfree(l.v);
	621	*f = v;
	622	*nf = l.nf;
	623	}
	624
	625	return (d);
	626	#endif
04361334	627	}
04361334	628
d28d625a	629	/----- That's all, folks -------------------------------------------------/