SICP - Solution: Exercise 1.28
October 21, 2018
Exercise 1.28 #
One variant of the Fermat test that cannot be fooled is called the Miller-Rabin test (Miller 1976; Rabin 1980). This starts from an alternate form of Fermat’s Little Theorem, which states that if $n$ is a prime number and $a$ is any positive integer less than $n$, then $a$ raised to the $(n−1)$-st power is congruent to 1 modulo n. To test the primality of a number $n$ by the Miller-Rabin test, we pick a random number ${a<n}$ and raise a to the ${(n-1)}$-st power modulo $n$ using the
expmod
procedure. However, whenever we perform the squaring step in expmod, we check to see if we have discovered a “nontrivial square root of 1 modulo n,” that is, a number not equal to 1 or n−1 whose square is equal to 1 modulo n. It is possible to prove that if such a nontrivial square root of 1 exists, then n is not prime. It is also possible to prove that if n is an odd number that is not prime, then, for at least half the numbers $a<n$, computing $a^{n-1}$ in this way will reveal a nontrivial square root of 1 modulo n. (This is why the Miller-Rabin test cannot be fooled.) Modify the expmod procedure to signal if it discovers a nontrivial square root of 1, and use this to implement the Miller-Rabin test with a procedure analogous to fermat-test. Check your procedure by testing various known primes and non-primes. Hint: One convenient way to make expmod signal is to have it return 0.
Solution #
Let’s break it down into pieces.
This starts from an alternate form of Fermat’s Little Theorem, which states that if $n$ is a prime number and $a$ is any positive integer less than $n$, then $a$ raised to the $(n−1)$-st power is congruent to 1 modulo n. To test the primality of a number $n$ by the Miller-Rabin test, we pick a random number ${a<n}$ and raise a to the ${(n-1)}$-st power modulo $n$ using the
expmod
procedure.
If we express this in terms of modular arithmetic, we have a congruence relation:
$$a^{n-1}\equiv1 (mod n)$$
Which will be implemented by updating the function fermat-test
into:
(define (miller-rabin-test n)
(define (try-it a)
(= (expmod-checked a (- n 1) n) 1))
(try-it (+ 1 (random (- n 1)))))
However, whenever we perform the squaring step in
expmod
, we check to see if we have discovered a “nontrivial square root of 1 modulo n,” that is, a number not equal to 1 or n−1 whose square is equal to 1 modulo n.
This can be translated into:
(define (remainder-square-checked x m)
(if (and (not (or (= x 1)
(= x (- m 1))))
(= (remainder (* x x) m) 1))
0
(remainder (* x x) m)))
by returning 0, when checking that (expmod-checked a (- n 1) n)
is equal to 1 in miller-rabin-test
will fail and indicate that the number is not prime.
Modify the expmod procedure to signal if it discovers a nontrivial square root of 1, and use this to implement the Miller-Rabin test with a procedure analogous to fermat-test. Check your procedure by testing various known primes and non-primes. Hint: One convenient way to make expmod signal is to have it return 0.
Putting all that together so that expmod-checked
works with our new functions:
(define (remainder-square-checked x m)
(if (and (not (or (= x 1)
(= x (- m 1))))
(= (remainder (* x x) m) 1))
0
(remainder (* x x) m)))
(define (expmod-checked base exp m)
(cond ((= exp 0) 1)
((even? exp)
(remainder-square-checked (expmod-checked base (/ exp 2) m) m))
(else
(remainder (* base (expmod-checked base (- exp 1) m))
m))))
(define (miller-rabin-test n)
(define (try-it a)
(= (expmod-checked a (- n 1) n) 1))
(try-it (+ 1 (random (- n 1)))))
(define (miller-rabin-prime? n times)
(cond ((= times 0) #t)
((miller-rabin-test n)
(miller-rabin-prime? n (- times 1)))
(else #f)))
(define (assert-result test-name actual expected)
(display (if (eq? actual expected) "pass: " "fail: "))
(display test-name)
(newline))
(assert-result " 2 is prime" (miller-rabin-prime? 2 10) #t)
(assert-result "1009 is prime" (miller-rabin-prime? 1009 10) #t)
(assert-result " 4 is not prime" (miller-rabin-prime? 4 10) #f)
(assert-result " 99 is not prime" (miller-rabin-prime? 99 10) #f)
(assert-result " 561 is not prime " (miller-rabin-prime? 561 10) #f) ; (Carmichael number)
That will output:
pass: 2 is prime
pass: 1009 is prime
pass: 4 is not prime
pass: 99 is not prime
pass: 561 is not prime
Open questions #
- Running it without
square-checked
(using normalsquare
) didn’t change the passing of these specific tests. How can we check that we have a correct implementation?