Rho assignments

Knuth's version of Pollard's rho algorithm

In section 4.5.4 of Knuth's TAOCP, he gives this version of Pollard's famous "rho" algorithm in Algorithm B (page 385 of the 3rd edition); the variable "N" indicates the number to be factored (the input argument, if you like.) (Note that "%" below represents the modulo operator).

B1: Initialization:

n ← N
x ← 5
x' ← 2
k ← 1
l ← 1

B2: Primality test

If n is prime, output n and terminate.

B3: Factor found?

g ← gcd(x' - x, n)
if g == 1, go to B4
if g == n, output "failed" and stop
output g (note that if g is not prime, its factors will not be found by this method)
n ← n / g
x ← x % n
x' ← x' % n
Go to B2

B4: Advance

k ← k - 1
if k == 0

x' ← x
l ← 2*l
k ← l

x ← (x*x + 1) % n
Go to B3

The details

Your assignment is to write programs implementing this function in C, Racket, and Gauche. Please don't use the factorization routines included in Racket's number theory module, such as factorize or prime-divisors, or those included in Gauche's math.prime module, such as naive-factorize or mc-factorize; instead, you must implement Knuth's version of the Rho algorithm given above.

Your programs should be able to factor the following numbers in a reasonable time (a few seconds, not a few minutes):

35
367756293167964701
2688307493077732358348134727
12902144921579660477018619711827

The programming languages that you are to use for the assignment are:

C, using GMP. Please name this source file "kr.c".
Gauche (/usr/bin/gosh), using math.prime. Please name this source file "kr.gosh"
Racket (/usr/bin/racket), using math/number-theory. Please name this source file "kr.racket"

Note that Gauche and Racket are both dialects of Scheme; once you write your code in one, you can largely re-use the same code for the other.

Your output should look like the following (order of factors reported is not important, though this algorithm tends to be better at finding smaller factors):

COP4020_test ~% ./kr  35
Attempting to factor N = 35 with Pollard rho algorithm from Knuth's Algorithm B in section 4.5.4 of TAOCP
Factors are:
                 7
                 5
COP4020_test ~% ./kr.gosh  35
Attempting to factor N = 35 with Pollard rho algorithm from Knuth's Algorithm B in section 4.5.4 of TAOCP
Factors are: 
 		 7
 		 5
COP4020_test ~% ./kr.racket  35
Attempting to factor N = 35 with Pollard rho algorithm from Knuth's Algorithm B in section 4.5.4 of TAOCP
Factors are: 
                 7
                 5
COP4020_test ~% ./kr  367756293167964701
Attempting to factor N = 367756293167964701 with Pollard rho algorithm from Knuth's Algorithm B in section 4.5.4 of TAOCP
Factors are:
                 7
                 11
                 13
                 107
                 101
                 103
                 109
                 3001
                 1009
COP4020_test ~% ./kr.gosh  367756293167964701
Attempting to factor N = 367756293167964701 with Pollard rho algorithm from Knuth's Algorithm B in section 4.5.4 of TAOCP
Factors are: 
 		 7
 		 11
 		 13
 		 107
 		 101
 		 103
 		 109
 		 3001
 		 1009
COP4020_test ~% ./kr.racket  367756293167964701
Attempting to factor N = 367756293167964701 with Pollard rho algorithm from Knuth's Algorithm B in section 4.5.4 of TAOCP
Factors are: 
                 7
                 11
                 13
                 107
                 101
                 103
                 109
                 3001
                 1009
COP4020_test ~% ./kr  2688307493077732358348134727
Attempting to factor N = 2688307493077732358348134727 with Pollard rho algorithm from Knuth's Algorithm B in section 4.5.4 of TAOCP
Factors are:
                 689816117
                 7688151521
                 506901611
COP4020_test ~% ./kr.gosh  2688307493077732358348134727
Attempting to factor N = 2688307493077732358348134727 with Pollard rho algorithm from Knuth's Algorithm B in section 4.5.4 of TAOCP
Factors are: 
 		 689816117
 		 7688151521
 		 506901611
COP4020_test ~% ./kr.racket  2688307493077732358348134727
Attempting to factor N = 2688307493077732358348134727 with Pollard rho algorithm from Knuth's Algorithm B in section 4.5.4 of TAOCP
Factors are: 
                 689816117
                 7688151521
                 506901611
COP4020_test ~% ./kr  12902144921579660477018619711827
Attempting to factor N = 12902144921579660477018619711827 with Pollard rho algorithm from Knuth's Algorithm B in section 4.5.4 of TAOCP
Factors are:
                 359195201
                 35919591591591616161427
COP4020_test ~% ./kr.gosh  12902144921579660477018619711827
Attempting to factor N = 12902144921579660477018619711827 with Pollard rho algorithm from Knuth's Algorithm B in section 4.5.4 of TAOCP
Factors are: 
 		 359195201
 		 35919591591591616161427
COP4020_test ~% ./kr.racket  12902144921579660477018619711827
Attempting to factor N = 12902144921579660477018619711827 with Pollard rho algorithm from Knuth's Algorithm B in section 4.5.4 of TAOCP
Factors are: 
                 359195201
                 35919591591591616161427

What makes these assignments interesting from a programming languages perspective?

C language implementations do not provide a native type that will allow you to do this with larger numbers; the question here is "How can we reasonably extend C to do multiprecision arithmetic?" The best answers so far are, rather than hacking something into the syntax of the language, to instead add the functionality via a library. This is a far more general answer, and leads to interesting questions about language design and the desirability of "yet more syntax" (a la, say, Perl.)
Scheme does not provide a convenient way to implement this particular algorithm, since Knuth gives it with several "go to" statements. Reworking the algorithm into a recursive formulation is not obvious, and takes some careful thought.

Submission

There are separate submission links for each of these assignments, and these each will be graded as a separate assignment. Please just turn in one file for each; for the C program, turn in your source file called "kr.c" and which should require the GMP library. For the Gauche version, turn in your source file called "kr.gosh". For the Racket version, turn in your source file called "kr.racket".

The assignment is due by 11:59pm on Sunday, October 27.