COS 226 Programming Assignment 8


Write a program to factor huge integers. Read an integer N from the command line, then print it as the product of two integers (or note that it is prime) as in the following brute-force solution:
     #include <stdio.h>
     #include <limits.h>
     main(int argc, char *argv[])
       { int d;
         int N = atoi(argv[1]);
         for (d = 2; d < N/d; d++)
             if (N % d == 0) break;
         if (N % d == 0) 
              printf("%d = %d*%d\n", N, d, N/d);
         else printf("%d is prime\n", N);
This code simply checks each integer greater than 1 to see whether it divides evenly into N. The header limits.h gives access to constants like INT_MAX that are not used by this program, but which you may need to avoid overflow in your own code.

We stop after finding any factor since we know that we could find all factors, if desired, by using the same method for the factors.

The solution just given is not quite a brute-force solution because it already contains a test that improves performance by a substantial amount: If we do not find a factor by the time that we reach the square root of N, then we certainly will not find a factor larger than that value, and we can declare N to be prime. A truly naive brute-force implementation might use time proportional to N (try every integer less than N; this improvement improves the running time to be proportional to the square root of N, so we can use it to factor any 32-bit integer in a fraction of a second. In this assignment, we explore the problem of factoring larger integers.

There are two reasons that the above code is not useful for factoring huge integers. First, the standard int data type supports 32-bit integers, so we would get overflow and/or incorrect results if we were to try to use it for larger numbers. Second, it is too slow. Even if we could deal with (say) 128-bit integers, we could never afford to try 2^64 divide operations to check for factors.

For this assignment, you will develop an implementation of a faster algorithm for factoring that is known as Pollard's rho method. The method is brought to us by the magic of number theory that is beyond the scope of this course, but it is an easy computation, based on iterating the function f(x) = x*x + c (mod N). Specifically, we start at some value x = x0 and maintain, for i = 1, 2, 3, ... , one variable a that tracks the result of i iterations of f and another variable b that tracks the result of 2i iterations of f, computing the greatest common divisor d of |a - b| and N. It turns out that whenever we find that d is not 1, we can infer that it is a factor of N, so we are done. (If you know number theory and are interested in details, see Pollard's paper in BIT 15, 1975, pp. 331-334.) For example, the following table traces the computation when we use it to factor 299.
x0 = 1
 c = 1
             a    b  |a-b|  d
             2    5    3    1
             5   79   74    1
            26  174  148    1
            79  105   26   13

          299 = 13*23
The key to a simple implementation is to note that we can maintain the values of a and b as required by simply setting a = f(a) and b = f(f(b)) in the inner loop. The loop always terminates, since a and b must eventually be equal (think about why this must be true!) and the greatest common divisor of 0 and N is N. Pollard's method is a randomized algorithm that factors with high probability if the starting value of a and the additive constant c are chosen at random.

Your first task for this assignment is to write a program small.c to factor int N that uses Pollard's method to compute the factorization. Instrument your code to also print the chosen values of x0 and c, the number of times the loop iterates, and, for N less than 1000, a table like the one shown above. Use your program to factor 713, 6319, 10609, and 42037. Completing this part of the assignment successfully is worth 6 points. Note that the numbers you have to factor are so small that you do not have to worry about overflow (see the first Extra Credit below).

You will need to implement a function to compute the greatest common divisor of two integers. For this purpose, use a recursive gcd function based on the fact that the greatest common divisor of a and b is the same as the greatest common divisor of b and a % b. (This classic method is known as Euclid's algorithm.) Incidentally, use a simple if test, not a library function, to compute the absolute value.

We will hardly notice the difference in performance between Pollard's algorithm and the brute-force method unless we want to factor huge integers. To develop an implementation for this purpose, we need an abstract data type (ADT) that supports all the usual arithmetic operations on such objects. This is an exercise in ADT design and implementation, a topic which is covered briefly in COS 126; in some detail in Chapter 4 of Algorithms in C; and in great detail in Hanson's book C Interfaces and Implementations and in COS 217. For this assignment, you can use the interface bignum.h. At Princeton, there are hundreds of implementations of this interface lying around, because it is a COS 217 assignment. A trustworthy implementation of this ADT is provided in ~cs226/prog8_files/ by including the objects bignum.o and libgmp.a during your compilation.

Your second task for this assignment is to write a program big.c that uses the BigNum ADT to factor huge integers, using Pollard's algorithm. Use your program to factor 5352499, 670726081, 9449868410449, 1082154235955237, 121932633334857493, and 13565005454706599869.

Experiment with other huge integers, and include in your readme file a discussion of the number of iterations that you think are required, in general.

Extra credit A. Implement small.c such that it is functionally equivalent to the sample program given above (but much faster). That is, your program needs to factor any int, avoiding overflow in the intermediate calculations.

Extra credit B. Find the largest prime pair (N and N+2 both prime) that you can in about ten seconds of computer time. An impressively large value of N will earn you an extra point.

Due: 11:59PM Sunday, April 18.