Skip to content

Instantly share code, notes, and snippets.

@Hermann-SW
Last active Aug 12, 2021
Embed
What would you like to do?
gmp-6.2.1/demos/factorize.c modified to alway use Pollard Rho with tortoise and hare (hangs when called for primes)
/* gcc -O3 -Wall -Wextra -pedantic -fomit-frame-pointer -m64 -mtune=skylake -march=broadwell -o factorize factorize.c -lgmp
$ ./factorize
factorize [n [a] [floyd|brent]]
Determines 1st factor for composite n with timing;
hangs indefinitely for prime n.
$
*/
/* Factoring with Pollard's rho method.
Copyright 1995, 1997-2003, 2005, 2009, 2012, 2015 Free Software
Foundation, Inc.
This program is free software; you can redistribute it and/or modify it under
the terms of the GNU General Public License as published by the Free Software
Foundation; either version 3 of the License, or (at your option) any later
version.
This program 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 General Public License for more details.
You should have received a copy of the GNU General Public License along with
this program. If not, see https://www.gnu.org/licenses/. */
#include <assert.h>
#include <inttypes.h>
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <sys/time.h>
#include <time.h>
#include <gmp.h>
long c;
void
factor_using_pollard_rho (mpz_t n, unsigned long a, mpz_t s)
{
mpz_t x, z;
mpz_t t, t2;
// gmp_printf ("rho(%Zd,%lu)\n", n, a);
mpz_inits (t, t2, NULL);
mpz_init_set (x, s);
mpz_init_set (z, s);
while (mpz_cmp_ui (n, 1) != 0)
{
for (;;)
{
++c;
mpz_mul (t, x, x);
mpz_mod (x, t, n);
mpz_add_ui (x, x, a);
mpz_mul (t, z, z);
mpz_mod (z, t, n);
mpz_add_ui (z, z, a);
mpz_mul (t, z, z);
mpz_mod (z, t, n);
mpz_add_ui (z, z, a);
mpz_sub (t2, z, x);
mpz_gcd (t, t2, n);
// gmp_printf (":%Zd %Zd %Zd\n", x, z, t);
if (mpz_cmp_ui (t, 1) != 0) { goto factor_found; }
}
}
factor_found:
if (mpz_cmp (t, n) == 0)
{
factor_using_pollard_rho(n, a+1, s);
}
else
{
mpz_divexact(t2, n, t);
gmp_printf ("%Zd %Zd (%Zd %Zd)\n", t, t2, x, z);
}
}
void
factor_using_pollard_rho_brent (mpz_t n, unsigned long a, mpz_t s)
{
mpz_t x, z, y, P;
mpz_t t, t2, ta;
unsigned long long k, l, i;
// gmp_printf ("rho(%Zd,%lu)\n", n, a);
mpz_inits (t, t2, NULL);
mpz_init_set (y, s);
mpz_init_set (x, s);
mpz_init_set (z, s);
mpz_init_set_ui (P, 1);
k = 1;
l = 1;
mpz_init_set_ui (ta, a);
while (mpz_cmp_ui (n, 1) != 0)
{
for (;;)
{
do
{
++c;
mpz_mul (t, x, x);
mpz_mod (x, t, n);
mpz_add (x, x, ta);
mpz_sub (t, z, x);
mpz_mul (t2, P, t);
mpz_mod (P, t2, n);
if (k % 32 == 1)
{
mpz_gcd (t, P, n);
if (mpz_cmp_ui (t, 1) != 0) { goto factor_found; }
mpz_set (y, x);
}
}
while (--k != 0);
mpz_set (z, x);
k = l;
l = 2 * l;
for (i = 0; i < k; i++)
{
mpz_mul (t, x, x);
mpz_mod (x, t, n);
mpz_add (x, x, ta);
}
mpz_set (y, x);
}
factor_found:
do
{
mpz_mul (t, y, y);
mpz_mod (y, t, n);
mpz_add_ui (y, y, a);
mpz_sub (t, z, y);
mpz_gcd (t, t, n);
}
while (mpz_cmp_ui (t, 1) == 0);
if (mpz_cmp (t, n) == 0)
{
factor_using_pollard_rho_brent(n, a+1, s);
}
else
{
mpz_divexact(t2, n, t);
gmp_printf ("%Zd %Zd\n", t, t2);
break;
}
}
}
int
main (int argc, char *argv[])
{
mpz_t t, s;
int brent = 0;
if (strcmp(argv[argc-1], "brent") == 0) { brent = 1; --argc; }
else if (strcmp(argv[argc-1], "floyd") == 0) { --argc; }
mpz_inits (t, s, NULL);
if (argc > 1)
{
struct timespec t0, t1;
unsigned long a = 1;
if (argc>2)
{
assert(strtol(argv[2], NULL, 0) >= 0);
a = (unsigned long) strtol(argv[2], NULL, 0);
}
mpz_set_str (t, argv[1], 0);
if (argc>3)
{ mpz_set_str (s, argv[3], 0); }
else
{ mpz_set_ui (s, 2); }
gmp_printf ("%Zd(%lu): ", t, a);
c=0;
clock_gettime(CLOCK_REALTIME, &t0);
if (brent)
{ factor_using_pollard_rho_brent (t, a, s); }
else
{ factor_using_pollard_rho (t, a, s); }
clock_gettime(CLOCK_REALTIME, &t1);
printf(" [%ld] %ldns (%s)\n", c, (t1.tv_sec-t0.tv_sec) * 1000000000 + t1.tv_nsec - t0.tv_nsec, brent ? "brent" : "floyd");
}
else
{
printf("factorize [n [a] [floyd|brent]]\n");
printf(" Determines 1st factor for composite n with timing;\n");
printf(" hangs indefinitely for prime n.\n");
}
exit (0);
}
@Hermann-SW

This comment has been minimized.

Copy link
Owner Author

@Hermann-SW Hermann-SW commented Aug 6, 2021

Discussion of enhanced output by latest diff (revision 5):
https://www.raspberrypi.org/forums/viewtopic.php?f=33&t=312344&p=1898717#p1898717

@Hermann-SW

This comment has been minimized.

Copy link
Owner Author

@Hermann-SW Hermann-SW commented Aug 7, 2021

Allow for optional 3rd argument start value (default is 2), revision 6.

@Hermann-SW

This comment has been minimized.

Copy link
Owner Author

@Hermann-SW Hermann-SW commented Aug 8, 2021

@Hermann-SW

This comment has been minimized.

Copy link
Owner Author

@Hermann-SW Hermann-SW commented Aug 12, 2021

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment