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

   1 /* -*-c-*-
   2  *
   3  * $Id: pgen.c,v 1.7 2001/02/03 16:05:32 mdw Exp $
   4  *
   5  * Prime generation glue
   6  *
   7  * (c) 1999 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: pgen.c,v $
  33  * Revision 1.7  2001/02/03 16:05:32  mdw
  34  * Now @mp_drop@ checks its argument is non-NULL before attempting to free
  35  * it.  Note that the macro version @MP_DROP@ doesn't do this.
  36  *
  37  * Revision 1.6  2000/10/08 12:11:22  mdw
  38  * Use @MP_EQ@ instead of @MP_CMP@.
  39  *
  40  * Revision 1.5  2000/06/17 11:52:36  mdw
  41  * Signal a pgen abort if the jump and base share a common factor.
  42  *
  43  * Revision 1.4  1999/12/22 16:01:11  mdw
  44  * Same file, completely different code.  Main interface for new prime-
  45  * search system.
  46  *
  47  */
  48
  49 /*----- Header files ------------------------------------------------------*/
  50
  51 #include <assert.h>
  52 #include <stdio.h>
  53 #include <stdlib.h>
  54 #include <string.h>
  55
  56 #include "fibrand.h"
  57 #include "grand.h"
  58 #include "mp.h"
  59 #include "mprand.h"
  60 #include "pgen.h"
  61 #include "pfilt.h"
  62 #include "rabin.h"
  63
  64 /*----- Standard prime filter ---------------------------------------------*/
  65
  66 /* --- @pgen_filter@ --- */
  67
  68 int pgen_filter(int rq, pgen_event *ev, void *p)
  69 {
  70   pgen_filterctx *f = p;
  71   int rc = PGEN_ABORT;
  72
  73   switch (rq) {
  74     case PGEN_BEGIN:
  75       rc = pfilt_create(&f->f, ev->m);
  76       mp_drop(ev->m);
  77       break;
  78     case PGEN_TRY:
  79       mp_drop(ev->m);
  80       if (!((f->step | f->f.m->v[0]) & 1))
  81         rc = pfilt_step(&f->f, 1);
  82       else
  83         rc = pfilt_step(&f->f, f->step);
  84       break;
  85     case PGEN_DONE:
  86       pfilt_destroy(&f->f);
  87       return (PGEN_DONE);
  88   }
  89
  90   while (rc == PGEN_FAIL)
  91     rc = pfilt_step(&f->f, f->step);
  92   ev->m = MP_COPY(f->f.m);
  93   return (rc);
  94 }
  95
  96 /* --- @pgen_jump@ --- *
  97  *
  98  * Similar to the standard @pgen_filter@, but jumps in large steps rather
  99  * than small ones.
 100  */
 101
 102 int pgen_jump(int rq, pgen_event *ev, void *p)
 103 {
 104   pgen_jumpctx *f = p;
 105   int rc = PGEN_ABORT;
 106
 107   switch (rq) {
 108     case PGEN_BEGIN: {
 109       mp *g = MP_NEW;
 110       mp_gcd(&g, 0, 0, ev->m, f->j->m);
 111       if (MP_CMP(g, >, MP_ONE)) {
 112         mp_drop(g);
 113         return (PGEN_ABORT);
 114       }
 115       mp_drop(g);
 116       rc = pfilt_create(&f->f, ev->m);
 117       mp_drop(ev->m);
 118     } break;
 119     case PGEN_TRY:
 120       mp_drop(ev->m);
 121       rc = pfilt_jump(&f->f, f->j);
 122       break;
 123     case PGEN_DONE:
 124       pfilt_destroy(&f->f);
 125       return (PGEN_DONE);
 126   }
 127
 128   while (rc == PGEN_FAIL)
 129     rc = pfilt_jump(&f->f, f->j);
 130   ev->m = MP_COPY(f->f.m);
 131   return (rc);
 132 }
 133
 134 /*----- Standard prime test -----------------------------------------------*/
 135
 136 /* --- @pgen_test@ --- */
 137
 138 int pgen_test(int rq, pgen_event *ev, void *p)
 139 {
 140   rabin *r = p;
 141   int rc = PGEN_ABORT;
 142
 143   switch (rq) {
 144     case PGEN_BEGIN:
 145       rabin_create(r, ev->m);
 146       rc = PGEN_TRY;
 147       break;
 148     case PGEN_TRY: {
 149       mp *a = mprand_range(MP_NEW, ev->m, ev->r, 0);
 150       rc = rabin_test(r, a);
 151       mp_drop(a);
 152     } break;
 153     case PGEN_DONE:
 154       rabin_destroy(r);
 155       rc = PGEN_DONE;
 156       break;
 157   }
 158
 159   return (rc);
 160 }
 161
 162 /*----- The main driver ---------------------------------------------------*/
 163
 164 /* --- @pgen@ --- *
 165  *
 166  * Arguments:   @const char *name@ = name of the value being searched for
 167  *              @mp *d@ = destination for the result integer
 168  *              @mp *m@ = start value to pass to stepper
 169  *              @pgen_proc *event@ = event handler function
 170  *              @void *ectx@ = context argument for event andler
 171  *              @unsigned steps@ = number of steps to take in search
 172  *              @pgen_proc *step@ = stepper function to use
 173  *              @void *sctx@ = context argument for stepper
 174  *              @unsigned tests@ = number of tests to make
 175  *              @pgen_proc *test@ = tester function to use
 176  *              @void *tctx@ = context argument for tester
 177  *
 178  * Returns:     Pointer to final result, or null.
 179  *
 180  * Use:         A generalized prime-number search skeleton.  Yes, that's a
 181  *              scary number of arguments.
 182  */
 183
 184 mp *pgen(const char *name, mp *d, mp *m, pgen_proc *event, void *ectx,
 185          unsigned steps, pgen_proc *step, void *sctx,
 186          unsigned tests, pgen_proc *test, void *tctx)
 187 {
 188   pgen_event ev;
 189   int rq, rc;
 190   pgen_proc *proc;
 191   void *ctx;
 192
 193   /* --- Set up the initial event block --- */
 194
 195   ev.name = name;
 196   if (m)
 197     ev.m = MP_COPY(m);
 198   else
 199     ev.m = 0;
 200   ev.steps = steps;
 201   ev.tests = tests;
 202   ev.r = fibrand_create(0);
 203
 204   /* --- Tell the event handler we're under way --- */
 205
 206   if (event && event(PGEN_BEGIN, &ev, ectx) == PGEN_ABORT)
 207     return (0);
 208
 209   /* --- Set up for the initial call --- */
 210
 211   proc = step; ctx = sctx; rq = PGEN_BEGIN;
 212
 213   /* --- Enter the great maelstrom of state transitions --- */
 214
 215   for (;;) {
 216     unsigned act = 0;
 217
 218     enum {
 219       A_STEP = 1u,
 220       A_TEST = 2u,
 221       A_EVENT = 4u,
 222       A_ENDTEST = 8u,
 223       A_ENDSTEP = 16u,
 224       A_DONE = 32u
 225     };
 226
 227     /* --- Call the procedure and decide what to do next --- */
 228
 229     rc = proc(rq, &ev, ctx);
 230     switch (rc) {
 231       case PGEN_TRY:
 232         if (proc == test)
 233           rq = PGEN_TRY;
 234         else {
 235           act |= A_EVENT;
 236           proc = test; ctx = tctx;
 237           rq = PGEN_BEGIN;
 238         }
 239         break;
 240       case PGEN_PASS:
 241         act |= A_TEST | A_EVENT;
 242         if (proc == test)
 243           rq = PGEN_TRY;
 244         else {
 245           proc = test; ctx = tctx;
 246           rq = PGEN_BEGIN;
 247         }
 248         break;
 249       case PGEN_FAIL:
 250         act |= A_STEP;
 251         if (proc == test) {
 252           act |= A_ENDTEST | A_EVENT;
 253           proc = step; ctx = sctx;
 254         }
 255         rq = PGEN_TRY;
 256         break;
 257       case PGEN_DONE:
 258         act |= A_EVENT | A_DONE | A_ENDSTEP;
 259         if (proc == test)
 260           act |= A_ENDTEST;
 261         break;
 262       case PGEN_ABORT:
 263         act |= A_EVENT | A_DONE;
 264         if (proc == test || rq == PGEN_TRY)
 265           act |= A_ENDSTEP;
 266         if (proc == test && rq == PGEN_BEGIN)
 267           act |= A_ENDTEST;
 268         break;
 269       default:
 270         assert(((void)"Invalid response from function", 0));
 271         break;
 272     }
 273
 274     /* --- If decrementing counters is requested, do that --- */
 275
 276     if ((act & A_STEP) && steps) {
 277       ev.steps--;
 278       if (!ev.steps) {
 279         act |= A_EVENT | A_ENDSTEP | A_DONE;
 280         rc = PGEN_ABORT;
 281       }
 282       ev.tests = tests;
 283     }
 284
 285     if ((act & A_TEST) && tests) {
 286       ev.tests--;
 287       if (!ev.tests) {
 288         act |= A_ENDTEST | A_ENDSTEP | A_DONE;
 289         rc = PGEN_DONE;
 290       }
 291     }
 292
 293     /* --- Report an event if so directed --- */
 294
 295     if ((act & A_EVENT) && event && event(rc, &ev, ectx) == PGEN_ABORT) {
 296       rc = PGEN_ABORT;
 297       if (!(act & A_DONE)) {
 298         act |= A_ENDSTEP | A_DONE;
 299         if (proc == test)
 300           act |= A_ENDTEST;
 301       }
 302     }
 303
 304     /* --- Close down tester and stepper functions --- */
 305
 306     if (act & A_ENDTEST)
 307       test(PGEN_DONE, &ev, tctx);
 308     if (act & A_ENDSTEP)
 309       step(PGEN_DONE, &ev, sctx);
 310
 311     /* --- Stop the entire test if necessary --- */
 312
 313     if (act & A_DONE)
 314       break;
 315   }
 316
 317   /* --- Tidy up and return --- */
 318
 319   if (rc == PGEN_ABORT) {
 320     mp_drop(ev.m);
 321     ev.m = 0;
 322   }
 323   ev.r->ops->destroy(ev.r);
 324   mp_drop(d);
 325
 326   return (ev.m);
 327 }
 328
 329 /*----- Test rig ----------------------------------------------------------*/
 330
 331 #ifdef TEST_RIG
 332
 333 #include <mLib/testrig.h>
 334
 335 static int verify(dstr *v)
 336 {
 337   mp *m = *(mp **)v[0].buf;
 338   mp *q = *(mp **)v[1].buf;
 339   mp *p;
 340   int ok = 1;
 341
 342   pgen_filterctx pf;
 343   rabin r;
 344
 345   pf.step = 2;
 346   p = pgen("p", MP_NEW, m, pgen_evspin, 0, 0, pgen_filter, &pf,
 347            rabin_iters(mp_bits(m)), pgen_test, &r);
 348   if (!p || !MP_EQ(p, q)) {
 349     fputs("\n*** pgen failed", stderr);
 350     fputs("\nm = ", stderr); mp_writefile(m, stderr, 10);
 351     fputs("\np = ", stderr); mp_writefile(p, stderr, 10);
 352     fputs("\nq = ", stderr); mp_writefile(q, stderr, 10);
 353     fputc('\n', stderr);
 354     ok = 0;
 355   }
 356
 357   mp_drop(m);
 358   mp_drop(q);
 359   mp_drop(p);
 360   assert(mparena_count(MPARENA_GLOBAL) == 0);
 361   return (ok);
 362 }
 363
 364 static test_chunk tests[] = {
 365   { "pgen", verify, { &type_mp, &type_mp, 0 } },
 366   { 0, 0, { 0 } }
 367 };
 368
 369 int main(int argc, char *argv[])
 370 {
 371   sub_init();
 372   test_run(argc, argv, tests, SRCDIR "/tests/pgen");
 373   return (0);
 374 }
 375 #endif
 376
 377 /*----- That's all, folks -------------------------------------------------*/