X-Git-Url: https://git.distorted.org.uk/u/mdw/catacomb/blobdiff_plain/eab06f16ebf8a3c42b4a372e5fe5d14fefc05b24..5278d9afdf1aff9fd6f64073ea42395d756ee58c:/pgen.c diff --git a/pgen.c b/pgen.c index 74a10bc..1a2a8ac 100644 --- a/pgen.c +++ b/pgen.c @@ -1,13 +1,13 @@ /* -*-c-*- * - * $Id: pgen.c,v 1.7 2001/02/03 16:05:32 mdw Exp $ + * $Id$ * * Prime generation glue * * (c) 1999 Straylight/Edgeware */ -/*----- Licensing notice --------------------------------------------------* +/*----- Licensing notice --------------------------------------------------* * * This file is part of Catacomb. * @@ -15,37 +15,18 @@ * 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: pgen.c,v $ - * Revision 1.7 2001/02/03 16:05:32 mdw - * Now @mp_drop@ checks its argument is non-NULL before attempting to free - * it. Note that the macro version @MP_DROP@ doesn't do this. - * - * Revision 1.6 2000/10/08 12:11:22 mdw - * Use @MP_EQ@ instead of @MP_CMP@. - * - * 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 @@ -68,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: @@ -77,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); @@ -124,7 +106,7 @@ int pgen_jump(int rq, pgen_event *ev, void *p) pfilt_destroy(&f->f); return (PGEN_DONE); } - + while (rc == PGEN_FAIL) rc = pfilt_jump(&f->f, f->j); ev->m = MP_COPY(f->f.m); @@ -145,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; @@ -189,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 --- */ @@ -197,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: @@ -274,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; } @@ -296,10 +285,10 @@ 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; } - } + } /* --- Close down tester and stepper functions --- */ @@ -326,12 +315,67 @@ mp *pgen(const char *name, mp *d, mp *m, pgen_proc *event, void *ectx, 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 #include +static int t_primep(dstr *v) +{ + mp *m = *(mp **)v[0].buf; + int e = *(int *)v[1].buf; + int r; + grand *rng; + int ok = 1; + + rng = fibrand_create(0); + r = pgen_primep(m, rng); + GR_DESTROY(rng); + if (e != r) { + fputs("\n*** primep failed", stderr); + fputs("\nm = ", stderr); mp_writefile(m, stderr, 10); + fprintf(stderr, "\nexpected %d", e); + fprintf(stderr, "\nreported %d", r); + fputc('\n', stderr); + ok = 0; + } + + mp_drop(m); + assert(mparena_count(MPARENA_GLOBAL) == 0); + return (ok); +} + static int verify(dstr *v) { mp *m = *(mp **)v[0].buf; @@ -363,6 +407,7 @@ static int verify(dstr *v) static test_chunk tests[] = { { "pgen", verify, { &type_mp, &type_mp, 0 } }, + { "primep", t_primep, { &type_mp, &type_int, 0 } }, { 0, 0, { 0 } } };