math/mpreduce.h: Missing include files.
[u/mdw/catacomb] / math / pgen-gcd.c
1 /* -*-c-*-
2 *
3 * Prime search stepper ensuring a low GCD for %$(p - 1)/2$%
4 *
5 * (c) 2000 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 "pgen.h"
32
33 /*----- Main code ---------------------------------------------------------*/
34
35 int pgen_gcdstep(int rq, pgen_event *ev, void *p)
36 {
37 pgen_gcdstepctx *g = p;
38 int rc = PGEN_ABORT;
39
40 switch (rq) {
41
42 /* --- Set everything up --- *
43 *
44 * Call things off if @p@ and @jp@ have common factors, or if @q@, @r@
45 * and @jq@ have common factors greater than @max@.
46 */
47
48 case PGEN_BEGIN: {
49 mp *p = ev->m;
50 mp_gcd(&g->g, 0, 0, p, g->jp.m);
51 if (MP_CMP(g->g, >, MP_ONE))
52 return (PGEN_ABORT);
53 g->q = mp_lsr(MP_NEW, p, 1);
54 g->jq = mp_lsr(MP_NEW, g->jp.m, 1);
55 mp_gcd(&g->g, 0, 0, g->q, g->jq);
56 mp_gcd(&g->g, 0, 0, g->g, g->r);
57 if (MP_CMP(g->g, >, g->max)) {
58 mp_drop(g->q);
59 mp_drop(g->jq);
60 return (PGEN_ABORT);
61 }
62 rc = pfilt_create(&g->p, p);
63 mp_drop(p);
64 } break;
65
66 /* --- Grind through another iteration --- */
67
68 case PGEN_TRY:
69 mp_drop(ev->m);
70 rc = pfilt_jump(&g->p, &g->jp);
71 g->q = mp_add(g->q, g->q, g->jq);
72 break;
73
74 /* --- Finished --- */
75
76 case PGEN_DONE:
77 pfilt_destroy(&g->p);
78 mp_drop(g->q);
79 mp_drop(g->jq);
80 return (PGEN_DONE);
81 }
82
83 /* --- Step on until everything is OK --- */
84
85 for (;;) {
86 if (rc != PGEN_FAIL) {
87 mp_gcd(&g->g, 0, 0, g->r, g->q);
88 if (MP_CMP(g->g, >, g->max))
89 rc = PGEN_FAIL;
90 }
91 if (rc != PGEN_FAIL)
92 break;
93 rc = pfilt_jump(&g->p, &g->jp);
94 g->q = mp_add(g->q, g->q, g->jq);
95 }
96
97 ev->m = MP_COPY(g->p.m);
98 return (rc);
99 }
100
101 /*----- That's all, folks -------------------------------------------------*/