X-Git-Url: https://git.distorted.org.uk/u/mdw/catacomb/blobdiff_plain/8b021c3f89a78c3006ffc5d480feca6ef86d544e..025c5f4aa5ffbf8948482a4233318db81c2df5d2:/pgen.c diff --git a/pgen.c b/pgen.c index 0e67c76..86e31c1 100644 --- a/pgen.c +++ b/pgen.c @@ -1,6 +1,6 @@ /* -*-c-*- * - * $Id: pgen.c,v 1.5 2000/06/17 11:52:36 mdw Exp $ + * $Id$ * * Prime generation glue * @@ -27,18 +27,6 @@ * MA 02111-1307, USA. */ -/*----- Revision history --------------------------------------------------* - * - * $Log: pgen.c,v $ - * Revision 1.5 2000/06/17 11:52:36 mdw - * Signal a pgen abort if the jump and base share a common factor. - * - * Revision 1.4 1999/12/22 16:01:11 mdw - * Same file, completely different code. Main interface for new prime- - * search system. - * - */ - /*----- Header files ------------------------------------------------------*/ #include @@ -61,7 +49,7 @@ int pgen_filter(int rq, pgen_event *ev, void *p) { pgen_filterctx *f = p; - int rc = PGEN_ABORT; + int rc = PGEN_FAIL; switch (rq) { case PGEN_BEGIN: @@ -70,16 +58,17 @@ int pgen_filter(int rq, pgen_event *ev, void *p) break; case PGEN_TRY: mp_drop(ev->m); - if (!((f->step | f->f.m->v[0]) & 1)) - rc = pfilt_step(&f->f, 1); - else - rc = pfilt_step(&f->f, f->step); break; case PGEN_DONE: pfilt_destroy(&f->f); return (PGEN_DONE); + default: + rc = PGEN_ABORT; + break; } + if (rc == PGEN_FAIL && !((f->step | f->f.m->v[0]) & 1)) + rc = pfilt_step(&f->f, 1); while (rc == PGEN_FAIL) rc = pfilt_step(&f->f, f->step); ev->m = MP_COPY(f->f.m); @@ -138,11 +127,15 @@ int pgen_test(int rq, pgen_event *ev, void *p) rabin_create(r, ev->m); rc = PGEN_TRY; break; - case PGEN_TRY: { - mp *a = mprand_range(MP_NEW, ev->m, ev->r, 0); - rc = rabin_test(r, a); - mp_drop(a); - } break; + case PGEN_TRY: + if (!ev->tests) + rc = rabin_rtest(r, MP_TWO); + else { + mp *a = mprand_range(MP_NEW, ev->m, ev->r, 0); + rc = rabin_rtest(r, a); + mp_drop(a); + } + break; case PGEN_DONE: rabin_destroy(r); rc = PGEN_DONE; @@ -182,6 +175,9 @@ mp *pgen(const char *name, mp *d, mp *m, pgen_proc *event, void *ectx, int rq, rc; pgen_proc *proc; void *ctx; + int p; + + enum { P_STEP, P_TEST }; /* --- Set up the initial event block --- */ @@ -190,73 +186,73 @@ mp *pgen(const char *name, mp *d, mp *m, pgen_proc *event, void *ectx, ev.m = MP_COPY(m); else ev.m = 0; - ev.steps = steps; - ev.tests = tests; + ev.steps = 0; + ev.tests = 0; ev.r = fibrand_create(0); /* --- Tell the event handler we're under way --- */ - if (event && event(PGEN_BEGIN, &ev, ectx) == PGEN_ABORT) + if (event && event(PGEN_BEGIN, &ev, ectx) == PGEN_ABORT) { + ev.r->ops->destroy(ev.r); return (0); + } /* --- Set up for the initial call --- */ - proc = step; ctx = sctx; rq = PGEN_BEGIN; + proc = step; ctx = sctx; p = P_STEP; rq = PGEN_BEGIN; /* --- Enter the great maelstrom of state transitions --- */ for (;;) { unsigned act = 0; - enum { - A_STEP = 1u, - A_TEST = 2u, - A_EVENT = 4u, - A_ENDTEST = 8u, - A_ENDSTEP = 16u, - A_DONE = 32u - }; +#define A_STEP 1u +#define A_TEST 2u +#define A_EVENT 4u +#define A_ENDTEST 8u +#define A_ENDSTEP 16u +#define A_DONE 32u /* --- Call the procedure and decide what to do next --- */ rc = proc(rq, &ev, ctx); switch (rc) { case PGEN_TRY: - if (proc == test) + if (p == P_TEST) rq = PGEN_TRY; else { act |= A_EVENT; - proc = test; ctx = tctx; + proc = test; ctx = tctx; p = P_TEST; rq = PGEN_BEGIN; } break; case PGEN_PASS: act |= A_TEST | A_EVENT; - if (proc == test) + if (p == P_TEST) rq = PGEN_TRY; else { - proc = test; ctx = tctx; + proc = test; ctx = tctx; p = P_TEST; rq = PGEN_BEGIN; } break; case PGEN_FAIL: act |= A_STEP; - if (proc == test) { + if (p == P_TEST) { act |= A_ENDTEST | A_EVENT; - proc = step; ctx = sctx; + proc = step; ctx = sctx; p = P_STEP; } rq = PGEN_TRY; break; case PGEN_DONE: act |= A_EVENT | A_DONE | A_ENDSTEP; - if (proc == test) + if (p == P_TEST) act |= A_ENDTEST; break; case PGEN_ABORT: act |= A_EVENT | A_DONE; - if (proc == test || rq == PGEN_TRY) + if (p == P_TEST || rq == PGEN_TRY) act |= A_ENDSTEP; - if (proc == test && rq == PGEN_BEGIN) + if (p == P_TEST && rq != PGEN_BEGIN) act |= A_ENDTEST; break; default: @@ -267,17 +263,17 @@ mp *pgen(const char *name, mp *d, mp *m, pgen_proc *event, void *ectx, /* --- If decrementing counters is requested, do that --- */ if ((act & A_STEP) && steps) { - ev.steps--; - if (!ev.steps) { + ev.steps++; + if (ev.steps == steps) { act |= A_EVENT | A_ENDSTEP | A_DONE; rc = PGEN_ABORT; } - ev.tests = tests; + ev.tests = 0; } if ((act & A_TEST) && tests) { - ev.tests--; - if (!ev.tests) { + ev.tests++; + if (ev.tests == tests) { act |= A_ENDTEST | A_ENDSTEP | A_DONE; rc = PGEN_DONE; } @@ -289,7 +285,7 @@ mp *pgen(const char *name, mp *d, mp *m, pgen_proc *event, void *ectx, rc = PGEN_ABORT; if (!(act & A_DONE)) { act |= A_ENDSTEP | A_DONE; - if (proc == test) + if (p == P_TEST) act |= A_ENDTEST; } } @@ -314,12 +310,41 @@ mp *pgen(const char *name, mp *d, mp *m, pgen_proc *event, void *ectx, ev.m = 0; } ev.r->ops->destroy(ev.r); - if (d != MP_NEW) - mp_drop(d); + mp_drop(d); return (ev.m); } +/* --- @pgen_primep@ --- * + * + * Arguments: @mp *p@ = a number to check + * @grand *gr@ = a random number source + * + * Returns: Nonzero if @p@ is really prime. + */ + +int pgen_primep(mp *p, grand *gr) +{ + int i; + rabin r; + mp *x = MP_NEW; + + if (MP_NEGP(p)) return (0); + switch (pfilt_smallfactor(p)) { + case PGEN_DONE: return (1); + case PGEN_FAIL: return (0); + } + rabin_create(&r, p); + for (i = 32; i; i--) { + x = mprand_range(x, p, gr, 0); + if (rabin_rtest(&r, x) == PGEN_FAIL) + break; + } + MP_DROP(x); + rabin_destroy(&r); + return (!i); +} + /*----- Test rig ----------------------------------------------------------*/ #ifdef TEST_RIG @@ -339,7 +364,7 @@ static int verify(dstr *v) pf.step = 2; p = pgen("p", MP_NEW, m, pgen_evspin, 0, 0, pgen_filter, &pf, rabin_iters(mp_bits(m)), pgen_test, &r); - if (!p || MP_CMP(p, !=, q)) { + if (!p || !MP_EQ(p, q)) { fputs("\n*** pgen failed", stderr); fputs("\nm = ", stderr); mp_writefile(m, stderr, 10); fputs("\np = ", stderr); mp_writefile(p, stderr, 10); @@ -350,8 +375,7 @@ static int verify(dstr *v) mp_drop(m); mp_drop(q); - if (p) - mp_drop(p); + mp_drop(p); assert(mparena_count(MPARENA_GLOBAL) == 0); return (ok); }