![[BACK]](/icons/cvsweb/back.gif){kind=icon}
Return to factorize.c CVS log ![[TXT]](/icons/cvsweb/text.gif){kind=icon}
![[DIR]](/icons/cvsweb/dir.gif){kind=icon}
Up to [local] / OpenXM_contrib / gmp / demos|
Return to factorize.c CVS log | Up to [local] / OpenXM_contrib / gmp / demos |
1.1 maekawa 1: /* Factoring with Pollard's rho method.
2:
3: Copyright (C) 1995 Free Software Foundation, Inc.
4:
5: This program is free software; you can redistribute it and/or modify it
6: under the terms of the GNU General Public License as published by the
7: Free Software Foundation; either version 2, or (at your option) any
8: later version.
9:
10: This program is distributed in the hope that it will be useful, but
11: WITHOUT ANY WARRANTY; without even the implied warranty of
12: MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
13: General Public License for more details.
14:
15: You should have received a copy of the GNU General Public License along
16: with this program; see the file COPYING. If not, write to the Free Software
17: Foundation, 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA. */
18:
19: #include <stdio.h>
20: #include "gmp.h"
21:
22: int flag_mersenne = 0;
23:
24: static unsigned add[] = {4, 2, 4, 2, 4, 6, 2, 6};
25:
26: factor_using_division (t, limit)
27: mpz_t t;
28: unsigned int limit;
29: {
30: mpz_t q, r;
31: unsigned long int f;
32: int i, ai;
33: unsigned *addv = add;
34:
35: mpz_init (q);
36: mpz_init (r);
37:
38: if (mpz_probab_prime_p (t, 50))
39: goto ready;
40:
41: for (;;)
42: {
43: mpz_tdiv_qr_ui (q, r, t, 2);
44: if (mpz_cmp_ui (r, 0) != 0)
45: break;
46: mpz_set (t, q);
47: printf ("2 ");
48: fflush (stdout);
49: if (mpz_probab_prime_p (t, 50))
50: goto ready;
51: }
52:
53: for (;;)
54: {
55: mpz_tdiv_qr_ui (q, r, t, 3);
56: if (mpz_cmp_ui (r, 0) != 0)
57: break;
58: mpz_set (t, q);
59: printf ("3 ");
60: fflush (stdout);
61: if (mpz_probab_prime_p (t, 50))
62: goto ready;
63: }
64:
65: for (;;)
66: {
67: mpz_tdiv_qr_ui (q, r, t, 5);
68: if (mpz_cmp_ui (r, 0) != 0)
69: break;
70: mpz_set (t, q);
71: printf ("5 ");
72: fflush (stdout);
73: if (mpz_probab_prime_p (t, 50))
74: goto ready;
75: }
76:
77: f = 7;
78: ai = 0;
79: for (;;)
80: {
81: mpz_tdiv_qr_ui (q, r, t, f);
82: if (mpz_cmp_ui (r, 0) != 0)
83: {
84: f += addv[ai];
85: if (f > limit)
86: goto ret;
87: ai = (ai + 1) & 7;
88: }
89: else
90: {
91: mpz_set (t, q);
92: printf ("%lu ", f);
93: fflush (stdout);
94: if (mpz_probab_prime_p (t, 50))
95: goto ready;
96: }
97: }
98:
99: ready:
100: mpz_out_str (stdout, 10, t);
101: fflush (stdout);
102: mpz_set_ui (t, 1);
103: fputc (' ', stdout);
104: ret:
105: mpz_clear (q);
106: mpz_clear (r);
107: }
108:
109: void
110: factor_using_pollard_rho (m, a_int, x0, p)
111: mpz_t m;
112: long a_int;
113: long x0;
114: unsigned long p;
115: {
116: mpz_t x, y, q;
117: mpz_t a;
118: mpz_t d;
119: mpz_t tmp;
120: mpz_t n;
121: int i = 1;
122: int j = 1;
123:
124: mpz_init_set (n, m);
125:
126: mpz_init (d);
127: mpz_init_set_ui (q, 1);
128: mpz_init (tmp);
129:
130: mpz_init_set_si (a, a_int);
131: mpz_init_set_si (x, x0);
132: mpz_init_set_si (y, x0);
133:
134: while (mpz_cmp_ui (n, 1) != 0)
135: {
136: if (flag_mersenne)
137: {
138: mpz_powm_ui (x, x, p, n); mpz_add (x, x, a);
139: mpz_powm_ui (y, y, p, n); mpz_add (y, y, a);
140: mpz_powm_ui (y, y, p, n); mpz_add (y, y, a);
141: }
142: else
143: {
144: mpz_mul (x, x, x); mpz_add (x, x, a); mpz_mod (x, x, n);
145: mpz_mul (y, y, y); mpz_add (y, y, a); mpz_mod (y, y, n);
146: mpz_mul (y, y, y); mpz_add (y, y, a); mpz_mod (y, y, n);
147: }
148:
149: if (mpz_cmp (x, y) > 0)
150: mpz_sub (tmp, x, y);
151: else
152: mpz_sub (tmp, y, x);
153: mpz_mul (q, q, tmp);
154: mpz_mod (q, q, n);
155:
156: if (++i % j == 0)
157: {
158: j += 1;
159: mpz_gcd (d, q, n);
160: if (mpz_cmp_ui (d, 1) != 0)
161: {
162: if (!mpz_probab_prime_p (d, 50))
163: factor_using_pollard_rho (d, (random () & 31) - 16,
164: (random () & 31), p);
165: else
166: {
167: mpz_out_str (stdout, 10, d);
168: fflush (stdout);
169: fputc (' ', stdout);
170: }
171: mpz_div (n, n, d);
172: if (mpz_probab_prime_p (n, 50))
173: {
174: mpz_out_str (stdout, 10, n);
175: fflush (stdout);
176: fputc (' ', stdout);
177: break;
178: }
179: }
180: }
181: }
182:
183: mpz_clear (n);
184: mpz_clear (d);
185: mpz_clear (q);
186: mpz_clear (tmp);
187: mpz_clear (a);
188: mpz_clear (x);
189: mpz_clear (y);
190: }
191:
192: factor (t, a, x0, p)
193: mpz_t t;
194: long a;
195: long x0;
196: unsigned long p;
197: {
198: factor_using_division (t, 1000000);
199: factor_using_pollard_rho (t, a, x0, p);
200: }
201:
202: main (argc, argv)
203: int argc;
204: char *argv[];
205: {
206: mpz_t t;
207: long x0, a;
208: unsigned long p;
209: int i;
210:
211: for (i = 1; i < argc; i++)
212: {
213: if (!strncmp (argv[i], "-Mp", 3))
214: {
215: p = atoi (argv[i] + 3);
216: mpz_init_set_ui (t, 1);
217: mpz_mul_2exp (t, t, p);
218: mpz_sub_ui (t, t, 1);
219: flag_mersenne = 1;
220: }
221: else
222: {
223: p = 0;
224: mpz_init_set_str (t, argv[i], 0);
225: }
226:
227: a = -1;
228: x0 = 3;
229:
230: factor (t, a, x0, p);
231: puts ("");
232: }
233: }