Rho method

Today I am going to be talking a little bit about Pollard's Rho algorithm, which, for such an apparently simple algorithm, is surprisingly good at finding moderate-size factors of general numbers.

Some reference links:
http://en.wikipedia.org/wiki/Pollard%27s_rho_algorithm
[as I mentioned briefly once before]
http://mersennewiki.org/index.php/Rho_Factorization_Method

Also, the pseudocode immediately below is taken from a PD document by Connelly Barnes of Oregon State University
http://oregonstate.edu/~barnesc/documents/factoring.pdf

function pollardRho(N)
# Initial values x(i) and x(2*i) for i = 0.
xi := 2
x2i := 2
do
# Find x(i+1) and x(2*(i+1))
xiPrime := xi ^ 2 + 1
x2iPrime := (x2i ^ 2 + 1) ^ 2 + 1
# Increment i: change our running values for x(i), x(2*i).
xi := xiPrime % N
x2i := x2iPrime % N
s := gcd(xi - x2i, N)
if s <> 1 and s <> N then
return s, N/s
end if
end do
end function


The Rho method is also the main method used by GMP's demo factorization sample program (or was last time I looked) - immediately below:

void
factor_using_pollard_rho (mpz_t n, int a_int, unsigned long p)
{
mpz_t x, x1, y, P;
mpz_t a;
mpz_t g;
mpz_t t1, t2;
int k, l, c, i;

if (flag_verbose)
{
printf ("[pollard-rho (%d)] ", a_int);
fflush (stdout);
}

mpz_init (g);
mpz_init (t1);
mpz_init (t2);

mpz_init_set_si (a, a_int);
mpz_init_set_si (y, 2);
mpz_init_set_si (x, 2);
mpz_init_set_si (x1, 2);
k = 1;
l = 1;
mpz_init_set_ui (P, 1);
c = 0;

while (mpz_cmp_ui (n, 1) != 0)
{
S2:
if (p != 0)
{
mpz_powm_ui (x, x, p, n); mpz_add (x, x, a);
}
else
{
mpz_mul (x, x, x); mpz_add (x, x, a); mpz_mod (x, x, n);
}
mpz_sub (t1, x1, x); mpz_mul (t2, P, t1); mpz_mod (P, t2, n);
c++;
if (c == 20)
{
c = 0;
mpz_gcd (g, P, n);
if (mpz_cmp_ui (g, 1) != 0)
goto S4;
mpz_set (y, x);
}
S3:
k--;
if (k > 0)
goto S2;

mpz_gcd (g, P, n);
if (mpz_cmp_ui (g, 1) != 0)
goto S4;

mpz_set (x1, x);
k = l;
l = 2 * l;
for (i = 0; i < k; i++)
{
if (p != 0)
{
mpz_powm_ui (x, x, p, n); mpz_add (x, x, a);
}
else
{
mpz_mul (x, x, x); mpz_add (x, x, a); mpz_mod (x, x, n);
}
}
mpz_set (y, x);
c = 0;
goto S2;
S4:
do
{
if (p != 0)
{
mpz_powm_ui (y, y, p, n); mpz_add (y, y, a);
}
else
{
mpz_mul (y, y, y); mpz_add (y, y, a); mpz_mod (y, y, n);
}
mpz_sub (t1, x1, y); mpz_gcd (g, t1, n);
}
while (mpz_cmp_ui (g, 1) == 0);

if (!mpz_probab_prime_p (g, 3))
{
do
{
mp_limb_t a_limb;
mpn_random (&a_limb, (mp_size_t) 1);
a_int = (int) a_limb;
}
while (a_int == -2 || a_int == 0);

if (flag_verbose)
{
printf ("[composite factor--restarting pollard-rho] ");
fflush (stdout);
}
factor_using_pollard_rho (g, a_int, p);
break;
}
else
{
mpz_out_str (stdout, 10, g);
fflush (stdout);
fputc (' ', stdout);
}
mpz_div (n, n, g);
mpz_mod (x, x, n);
mpz_mod (x1, x1, n);
mpz_mod (y, y, n);
if (mpz_probab_prime_p (n, 3))
{
mpz_out_str (stdout, 10, n);
fflush (stdout);
fputc (' ', stdout);
break;
}
}

mpz_clear (g);
mpz_clear (P);
mpz_clear (t2);
mpz_clear (t1);
mpz_clear (a);
mpz_clear (x1);
mpz_clear (x);
mpz_clear (y);
}


This Java code is taken directly from superfac9:

BigInteger factorizerho(BigInteger n) {
BigInteger loop = new BigInteger("1");
BigInteger x = new BigInteger("5");
BigInteger y = new BigInteger("2");

while (n.gcd(x.subtract(y)).compareTo(ONE) == 0 && loop.compareTo(TENTHOUSAND) < 0) {
x = x.multiply(x).add(ONE).mod(n);
x = x.multiply(x).add(ONE).mod(n);
y = y.multiply(y).add(ONE).mod(n);
loop = loop.add(ONE);
}

return n.gcd(x.subtract(y));

}


Note that there are several variants (given different start-point values for x and y for example), but the basic "Monte-Carlo" principle remains the same.

Notable successes with the Rho method or its variants include Brent and Pollard's factorization of the eighth Fermat number in 1980. The 16-digit factor of F8 took about 2 hours to find on a UNIVAC 1100/42, apparently.