![[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: }