X-Git-Url: https://git.distorted.org.uk/~mdw/catacomb/blobdiff_plain/2af1930eafd4b2747c3736dec408be25e6a8b737..0f5ec153dcb1c20a765d0420cc189ebac1b9682d:/pgen.c diff --git a/pgen.c b/pgen.c new file mode 100644 index 00000000..503b153d --- /dev/null +++ b/pgen.c @@ -0,0 +1,268 @@ +/* -*-c-*- + * + * $Id: pgen.c,v 1.1 1999/11/19 13:17:57 mdw Exp $ + * + * Finding and testing prime numbers + * + * (c) 1999 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: pgen.c,v $ + * Revision 1.1 1999/11/19 13:17:57 mdw + * Prime number generator and tester. + * + */ + +/*----- Header files ------------------------------------------------------*/ + +#include "mp.h" +#include "mpmont.h" +#include "pgen.h" +#include "ptab.h" +#include "rabin.h" + +/*----- Main code ---------------------------------------------------------*/ + +/* --- @pgen_create@ --- * + * + * Arguments: @pgen *p@ = pointer to prime generation context + * @mp *m@ = pointer to initial number to test + * + * Returns: One of the @PGEN@ constants above. + * + * Use: Tests an initial number for primality by computing its + * residue modulo various small prime numbers. This is fairly + * quick, but not particularly certain. If a @PGEN_MAYBE@ + * result is returned, perform Rabin-Miller tests to confirm. + */ + +int pgen_create(pgen *p, mp *m) +{ + int rc = PGEN_MAYBE; + int i; + mp *r = MP_NEW; + mpw qw; + mp q; + + /* --- Take a copy of the number --- */ + + mp_shrink(m); + p->m = MP_COPY(m); + + /* --- Fill in the residues --- */ + + mp_build(&q, &qw, &qw + 1); + for (i = 0; i < NPRIME; i++) { + qw = ptab[i]; + mp_div(0, &r, m, &q); + p->r[i] = r->v[0]; + if (!p->r[i] && rc == PGEN_MAYBE) { + if (MP_LEN(m) == 1 && m->v[0] == ptab[i]) + rc = PGEN_PRIME; + else + rc = PGEN_COMPOSITE; + } + } + + /* --- Done --- */ + + mp_drop(r); + return (rc); +} + +/* --- @pgen_destroy@ --- * + * + * Arguments: @pgen *p@ = pointer to prime generation context + * + * Returns: --- + * + * Use: Discards a context and all the resources it holds. + */ + +void pgen_destroy(pgen *p) +{ + mp_drop(p->m); +} + +/* --- @pgen_step@ --- * + * + * Arguments: @pgen *p@ = pointer to prime generation context + * @mpw step@ = how much to step the number + * + * Returns: One of the @PGEN@ constants above. + * + * Use: Steps a number by a small amount. Stepping is much faster + * than initializing with a new number. The test performed is + * the same simple one used by @ptab_create@, so @PGEN_MAYBE@ + * results should be followed up by a Rabin-Miller test. + */ + +int pgen_step(pgen *p, mpw step) +{ + mp s; + int rc = PGEN_MAYBE; + int i; + + /* --- Add the step on to the number --- */ + + mp_build(&s, &step, &step + 1); + p->m = mp_add(p->m, p->m, &s); + + /* --- Update the residue table --- */ + + for (i = 0; i < NPRIME; i++) { + p->r[i] = (p->r[i] + step) % ptab[i]; + if (!p->r[i] && rc == PGEN_MAYBE) { + if (MP_LEN(p->m) == 1 && p->m->v[0] == ptab[i]) + rc = PGEN_PRIME; + else + rc = PGEN_COMPOSITE; + } + } + + /* --- Done --- */ + + return (rc); +} + +/* --- @pgen_jump@ --- * + * + * Arguments: @pgen *p@ = pointer to prime generation context + * @pgen *j@ = pointer to another generation context + * + * Returns: One of the @PGEN@ constants above. + * + * Use: Steps a number by a large amount. Even so, jumping is much + * faster than initializing a new number. The test peformed is + * the same simple one used by @ptab_create@, so @PGEN_MAYBE@ + * results should be followed up by a Rabin-Miller test. + * + * Note that the number stored in the @j@ context is probably + * better off being even than prime. The important thing is + * that all of the residues for the number have already been + * computed. + */ + +int pgen_jump(pgen *p, pgen *j) +{ + int rc = PGEN_MAYBE; + int i; + + /* --- Add the step on --- */ + + p->m = mp_add(p->m, p->m, j->m); + + /* --- Update the residue table --- */ + + for (i = 0; i < NPRIME; i++) { + p->r[i] = p->r[i] + j->r[i]; + if (p->r[i] > ptab[i]) + p->r[i] -= ptab[i]; + if (!p->r[i] && rc == PGEN_MAYBE) { + if (MP_LEN(p->m) == 1 && p->m->v[0] == ptab[i]) + rc = PGEN_PRIME; + else + rc = PGEN_COMPOSITE; + } + } + + /* --- Done --- */ + + return (rc); +} + +/*----- Test code ---------------------------------------------------------*/ + +#ifdef TEST_RIG + +#include + +static int verify(dstr *v) +{ + mp *m = *(mp **)v[0].buf; + mp *r = *(mp **)v[1].buf; + pgen p; + mpw gw; + mp g; + int rc; + static char baton[] = "/-\\|"; + char *b = baton; + int ok = 1; + + mp_build(&g, &gw, &gw + 1); + rc = pgen_create(&p, m); + for (;;) { + if (rc == PGEN_PRIME) + break; + if (rc == PGEN_MAYBE) { + int i; + rabin r; + rabin_create(&r, p.m); + for (i = 0; i < 5; i++) { + gw = rand(); + if (rabin_test(&r, &g) == PGEN_COMPOSITE) + break; + } + rabin_destroy(&r); + if (i == 5) + break; + putc(*b++, stderr); + putc('\b', stderr); + if (!*b) + b = baton; + } + rc = pgen_step(&p, 2); + } + + if (MP_CMP(p.m, !=, r)) { + fputs("\n*** failed.", stderr); + fputs("\nm = ", stderr); mp_writefile(m, stderr, 10); + fputs("\np = ", stderr); mp_writefile(p.m, stderr, 10); + fputs("\nr = ", stderr); mp_writefile(r, stderr, 10); + fputc('\n', stderr); + ok = 0; + } + + mp_drop(m); + mp_drop(r); + pgen_destroy(&p); + return (ok); +} + +static test_chunk tests[] = { + { "pgen", verify, { &type_mp, &type_mp, 0 } }, + { 0, 0, { 0 } } +}; + +int main(int argc, char *argv[]) +{ + sub_init(); + test_run(argc, argv, tests, SRCDIR "/tests/pgen"); + return (0); +} + +#endif + +/*----- That's all, folks -------------------------------------------------*/