Lecture 17 -  Pollard’s factoring algorithms

This lecture (actually 2 in class) concern two clever factoring algorithms introduced by J. Pollard.  The first, Pollard's ρ algorithm for integer factorization (as distinct from his ρ algorithm for discrete logarithms) factors an integer n in time O((n^(1/4))), as opposed to the naive  O((n^(1/2))) trial division algorithm.  The second is a special purpose algorithm, and only works well when p-1 is smooth for one of the factors p of the integer to be factored.

The idea of both algorithms is to hypothetically consider a factor p of n, and derive properties it satisfies from some experimentation.  Enough clues are found to obtain the value of p.

Pollard’ s ρ algorithm for integer factorization

Let p be the smallest factor of an integer n (p is unknown, of course).  We consider a random-looking function that maps integers mod n to other integers mod n.  The classic example is to take the polynomial f(x) = (x^2)+5 and iterate it. 
Let (Subscript[x, 0])=1, (Subscript[x, 1])= f((Subscript[x, 0])), (Subscript[x, 2])= f((Subscript[x, 1])), ..., (Subscript[x, k+1])=f((Subscript[x, k])), ....

 The Birthday Paradox tells us that it is unlikely two of the values will overlap until we have iterated this function about (n^(1/2))times.  However, the same Birthday Paradox tells us we do expect an overlap of values modulo the prime p within about (p^(1/2))≤ (n^(1/4)) steps.   Suppose that (Subscript[x, i])=(Subscript[x, j])is such a collision (mod p).  Of course we don't know which i and j these are, but if we took the GCD((Subscript[x, i])-(Subscript[x, j]),n) it would very likely be p.  After all, it has to be a divisor of n, and p divides both.  It's not likely to include a bigger factor, such as n itself, because it would mean a collision mod n (which we don't expect to happen until much later on).  

The algorithm thus tests for such collisions.  However, storing them all and comparing them all would be prohibitive.  Instead, it uses the fact that once there is a collision of values (Subscript[x, i])=(Subscript[x, j]), then there are collisions every successive step:
(Subscript[x, i+1])=(Subscript[x, j+1])
(Subscript[x, i+2])=(Subscript[x, j+2])
(Subscript[x, i+3])=(Subscript[x, j+3])
etc...

We don't need to see the first collision, just some collision.  In general, (Subscript[x, r ])= (Subscript[x, r+|i-j|])for any sufficiently large r.  In particular, for sufficiently large multiples k of |i-j| we have (Subscript[x, k])= (Subscript[x, 2k]).  We thus only look at differences of these (which requires little storage), and take GCD((Subscript[x, k])- (Subscript[x, 2k]),n).

Here is an example:

f[x_]=x^2+5

p=Prime[4];q=Prime[6];n= p q

(* This is the sequence *)
tab={1};Do[AppendTo[tab,Mod[f[tab[[-1]]],n]],{100}];tab

(* These are the GCDs *)
Table[GCD[%[[2i]]-%[[i]],n],{i,50}]

We see 7 is pulled out .
  
  Here is a much bigger number, but with a moderately - sized smallest prime factor.  This determines how long the algorithm takes.

p=Prime[100];q=Prime[10000000000];n=p q

tab={1};Do[AppendTo[tab,Mod[f[tab[[-1]]],n]],{100}];tab

Table[GCD[%[[2i]]-%[[i]],n],{i,50}]

Sure enough, 541 is a prime factor!

Pollard' s p - 1 algorithm

As we mentioned above, this algorithm only works well when p-1 is smooth, meaning that it is the product of small primes (exactly how small can be quantified in terms of its running time).   It uses the same idea of trying to find two numbers which are equal mod p, but very unlikely to be equal mod n.  (As before, we don't know p, so this method doesn't work very well when the numbers are random.)

Let us recall the generalization of Fermat's Little Theorem to composite numbers, (a^(ϕ(m)  ))is congruent to 1 modulo n, if a and n are relatively prime.  In practice, we can assume a and n are relative prime, for if they had a common factor, we'd've already factored n.

For simplicity, let's look at n=p q, where p and q are both primes.  It's easiest to motivate this with an example, with p =37 (here p-1=(2^2)(3^2)).  Let's think of a number a mod p as (g^x), where x is the exponent of a generator g (mod p).  Then (a^2)is (g^(2x))mod p, and repeated squaring gives us (a^4)=(g^((2^2)x)) (mod p).   Since x is determined mod p-1=36, this tells us that (a^4)=(g^y), where y is a multiple of (2^2).  Now by repeated cubing, we get (a^12)=(g^(3y))and (a^36)=(g^((3^2)y)).  So we'd get a multiple of 36 by some repeated squarings and cubings.

Of course, here we knew the value of p-1 to begin with.  However, even if we didn't, we could do some repeated squarings and cubings to get it quickly -- simply because p-1 is smooth, and hence a product of small powers of 2 and small powers of 3.  It doesn't take too many tries to find all of them: if p-1=(2^x)(3^y), you never need to try x bigger than the logarithm of p to the base 2, or try y bigger than the logarithm of p to the base 3. If you threw in more small primes (5,7,11,...) this wouldn't change much -- as long as you didn't throw in too many more, it still wouldn't take too many combinations to hit it.

Of course, in practice, we don't know p-1, or even its factorization.  But if it is indeed smooth, then we could take some random number a and keep taking small powers of it, hoping to get back to the identity.  Then by looking at the GCD of this power-1, we'd expect to get p.

In this example, I picked the 1076th prime because it, minus 1, doesn't have a factor bigger than 5.  Again, in practice you would not know this!

n=Prime[1076]Prime[351321]

a=2;Table[a=PowerMod[a,k,n];a,{k,1,10}]

a=2;Table[a=PowerMod[a,k,n];GCD[a-1,n],{k,1,10}]

FactorInteger[n]

n=Prime[30]Prime[10]

a=2;Table[a=PowerMod[a,k,n];a,{k,1,10}]

a=2;Table[a=PowerMod[a,k,n];GCD[a-1,n],{k,1,10}]

FactorInteger[n]