chiark - git - mdw - catacomb/blob - math/prim.c

   1 /* -*-c-*-
   2  *
   3  * Finding primitive elements
   4  *
   5  * (c) 1999 Straylight/Edgeware
   6  */
   7
   8 /*----- Licensing notice --------------------------------------------------*
   9  *
  10  * This file is part of Catacomb.
  11  *
  12  * Catacomb is free software; you can redistribute it and/or modify
  13  * it under the terms of the GNU Library General Public License as
  14  * published by the Free Software Foundation; either version 2 of the
  15  * License, or (at your option) any later version.
  16  *
  17  * Catacomb is distributed in the hope that it will be useful,
  18  * but WITHOUT ANY WARRANTY; without even the implied warranty of
  19  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  20  * GNU Library General Public License for more details.
  21  *
  22  * You should have received a copy of the GNU Library General Public
  23  * License along with Catacomb; if not, write to the Free
  24  * Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
  25  * MA 02111-1307, USA.
  26  */
  27
  28 /*----- Header files ------------------------------------------------------*/
  29
  30 #include "mp.h"
  31 #include "mpint.h"
  32 #include "mpmont.h"
  33 #include "mprand.h"
  34 #include "pgen.h"
  35 #include "prim.h"
  36
  37 /*----- Main code ---------------------------------------------------------*/
  38
  39 /* --- @prim_test@ --- */
  40
  41 int prim_test(int rq, pgen_event *ev, void *p)
  42 {
  43   prim_ctx *c = p;
  44   int rc = rq;
  45
  46   switch (rq) {
  47     case PGEN_BEGIN:
  48       return (PGEN_TRY);
  49     case PGEN_TRY: {
  50       mp *x;
  51       rc = PGEN_FAIL;
  52
  53       if (!c->exp)
  54         x = mp_copy(ev->m);
  55       else {
  56         x = mpmont_exp(&c->mm, MP_NEW, ev->m, c->exp);
  57         if (MP_EQ(x, MP_ONE))
  58           goto done;
  59       }
  60       if (c->n == 0)
  61         goto ok;
  62       else {
  63         size_t n = c->n;
  64         mp **f = c->f;
  65         mp *y = MP_NEW;
  66         while (n) {
  67           y = mpmont_exp(&c->mm, y, x, *f);
  68           if (MP_EQ(y, MP_ONE)) {
  69             mp_drop(y);
  70             goto done;
  71           }
  72           n--; f++;
  73         }
  74         mp_drop(y);
  75       }
  76     ok:
  77       rc = PGEN_DONE;
  78       mp_drop(ev->m);
  79       ev->m = x;
  80       break;
  81     done:
  82       mp_drop(x);
  83     } break;
  84   }
  85
  86   return (rc);
  87 }
  88
  89 /* --- Trivial stepping functions -----------------------------------------*/
  90
  91 /* --- @prim_step@ --- */
  92
  93 int prim_step(int rq, pgen_event *ev, void *p)
  94 {
  95   unsigned *i = p;
  96   switch (rq) {
  97     case PGEN_BEGIN:
  98     case PGEN_TRY:
  99       if (*i >= NPRIME)
 100         return PGEN_FAIL;
 101       ev->m = mp_fromint(ev->m, primetab[(*i)++]);
 102       return (PGEN_TRY);
 103   }
 104   return (0);
 105 }
 106
 107 /*----- That's all, folks -------------------------------------------------*/