Structure and Interpretation of Computer Programs (11 page)

The greatest common divisor (GCD) of two integers
a
and
b
is
defined to be the largest integer that divides both
a
and
b
with no remainder. For example, the GCD of 16 and 28 is 4. In chapter 2,
when we investigate how to implement rational-number arithmetic, we
will need to be able to compute GCDs in order to reduce
rational numbers to lowest terms. (To reduce a rational number to
lowest terms, we must divide both the numerator and the denominator by their
GCD. For example, 16/28 reduces to 4/7.) One way to find the
GCD of two integers is to factor them and search for common
factors, but there is a famous algorithm that is much more efficient.

The idea of the algorithm is based on the observation that, if
r
is
the remainder when
a
is divided by
b
, then the common divisors of
a
and
b
are precisely the same as the common divisors of
b
and
r
. Thus, we can use the equation

GCD(
a
,
b
) = GCD(
b
,
r
)


GCD(206, 40) = GCD(40, 6)
= GCD(6, 4)
= GCD(4, 2)
= GCD(2, 0)
= 2

reduces GCD(206,40) to GCD(2,0), which is 2. It is
possible to show that starting with any two positive integers and
performing repeated reductions will always eventually produce a pair
where the second number is 0. Then the GCD is the other
number in the pair. This method for computing the GCD is
known as
Euclid's Algorithm
.
⁴²

It is easy to express Euclid's Algorithm as a procedure:

(define (gcd a b)
  (if (= b 0)
      a
      (gcd b (remainder a b))))

This generates an iterative process, whose number of steps grows as
the logarithm of the numbers involved.

The fact that the number of steps required by Euclid's Algorithm has
logarithmic growth bears an interesting relation to the Fibonacci
numbers:

Lamé's Theorem:
If Euclid's Algorithm requires
k
steps to
compute the GCD of some pair, then the smaller number in the pair
must be greater than or equal to the
k
th Fibonacci
number.
⁴³

We can use this theorem to get an order-of-growth estimate for Euclid's
Algorithm. Let
n
be the smaller of the two inputs to the
procedure. If the process takes
k
steps, then we must have
n
>
F
i
b
(
k
) ≈
φ
^k
/√5. Therefore
the number of steps
k
grows as the logarithm (to the base
φ
) of
n
. Hence, the order of growth is θ(
log
n
).

Exercise 1.20.
  
The process that a procedure generates is of course dependent on the
rules used by the interpreter. As an example, consider the iterative
gcd
procedure given above.
Suppose we were to interpret this procedure using normal-order
evaluation, as discussed in section 
1.1.5
.
(The normal-order-evaluation rule for
if
is described in
exercise 
1.5
.) Using the
substitution method (for normal order), illustrate the process
generated in evaluating
(gcd 206 40)
and indicate the
remainder
operations that are actually performed.
How many
remainder
operations are actually performed
in the normal-order evaluation of
(gcd 206 40)
?
In the applicative-order evaluation?

This section describes two methods for checking the primality of an
integer
n
, one with order of growth θ(√
n
), and a
“probabilistic” algorithm with order of growth θ(
log
n
). The
exercises at the end of this section suggest programming
projects based on these algorithms.

Since ancient times, mathematicians have been fascinated by problems
concerning prime numbers, and many people have worked on the problem
of determining ways to test if numbers are prime. One way
to test if a number is prime is to find the number's divisors. The
following program finds the smallest integral divisor (greater than 1)
of a given number 
n
. It does this in a straightforward way, by
testing
n
for divisibility by successive integers starting with 2.

(define (smallest-divisor n)
  (find-divisor n 2))
(define (find-divisor n test-divisor)
  (cond ((> (square test-divisor) n) n)
        ((divides? test-divisor n) test-divisor)
        (else (find-divisor n (+ test-divisor 1)))))
(define (divides? a b)
  (= (remainder b a) 0))

We can test whether a number is prime as follows:
n
is prime if
and only if
n
is its own smallest divisor.

(define (prime? n)
  (= n (smallest-divisor n)))

The end test for
find-divisor
is based on the fact that if
n
is not prime it must have a divisor less than or equal to
√
n
.
⁴⁴
This
means that the algorithm need only test divisors between 1 and
√
n
. Consequently, the number of steps required to identify
n
as prime will have order of growth θ(√
n
).

The θ(
log
n
) primality test is based on a result from number
theory known as Fermat's Little Theorem.
⁴⁵

Fermat's Little Theorem:
If
n
is a prime number and
a
is any positive integer less than
n
, then
a
raised to the
n
th power is congruent to
a
modulo
n
.

(Two numbers are said to be
congruent modulo
n
if
they both have the same remainder when divided by
n
. The
remainder of a number
a
when divided by
n
is also referred to as
the
remainder of
a
modulo
n
, or simply as
a
modulo
n
.)

If
n
is not prime, then, in general, most of the numbers
a
<
n
will not
satisfy the above relation. This leads to the following algorithm for
testing primality: Given a number
n
, pick a
random number
a
<
n
and
compute the remainder of
a
ⁿ
modulo
n
. If the result is not equal to
a
, then
n
is certainly not prime. If it is
a
, then chances are good
that
n
is prime. Now pick another random number
a
and test it with the
same method. If it also satisfies the equation, then we can be even more
confident that
n
is prime. By trying more and more values of
a
, we can
increase our confidence in the result. This algorithm is known as the
Fermat test.

To implement the Fermat test, we need a procedure that computes the
exponential of a number modulo another number:

(define (expmod base exp m)
  (cond ((= exp 0) 1)
        ((even? exp)
         (remainder (square (expmod base (/ exp 2) m))
                    m))
        (else
         (remainder (* base (expmod base (- exp 1) m))
                    m))))        

This is very similar to the
fast-expt
procedure of
section 
1.2.4
. It uses successive squaring, so
that the number of steps grows logarithmically with the
exponent.
⁴⁶

The Fermat test is performed by choosing at random a number
a
between 1 and
n
- 1 inclusive and checking whether the remainder
modulo
n
of the
n
th power of
a
is equal to
a
. The random
number
a
is chosen using the procedure
random
, which we assume is
included as a primitive in Scheme.
Random
returns a
nonnegative integer less than its integer input. Hence, to obtain a random
number between 1 and
n
- 1, we call
random
with an input of
n
- 1 and add 1 to the result:

(define (fermat-test n)
  (define (try-it a)
    (= (expmod a n n) a))
  (try-it (+ 1 (random (- n 1)))))

The following procedure runs the test a given number of times, as
specified by a parameter. Its value is true if the test succeeds
every time, and false otherwise.

(define (fast-prime? n times)
  (cond ((= times 0) true)
        ((fermat-test n) (fast-prime? n (- times 1)))
        (else false)))

The Fermat test differs in character from most familiar algorithms, in
which one computes an answer that is guaranteed to be correct. Here,
the answer obtained is only probably correct. More precisely, if
n
ever fails the Fermat test, we can be certain that
n
is not prime.
But the fact that
n
passes the test, while an extremely strong
indication, is still not a guarantee that
n
is prime. What we would
like to say is that for any number
n
, if we perform the test enough
times and find that
n
always passes the test, then the probability
of error in our primality test can be made as small as we like.

Unfortunately, this assertion is not quite correct. There do exist
numbers that fool the Fermat test: numbers
n
that are not prime and
yet have the property that
a
ⁿ
is congruent to
a
modulo
n
for
all integers
a
<
n
. Such numbers are extremely rare, so the Fermat
test is quite reliable in practice.
⁴⁷
There are variations of the Fermat test that cannot be fooled. In
these tests, as with the Fermat method, one tests the primality of an
integer
n
by choosing a random integer
a
<
n
and checking some
condition that depends upon
n
and
a
. (See
exercise 
1.28
for an example of such a test.) On the
other hand, in contrast to the Fermat test, one can prove that, for
any
n
, the condition does not hold for most of the integers
a
<
n
unless
n
is prime. Thus, if
n
passes the test for some random
choice of 
a
, the chances are better than even that
n
is prime. If
n
passes the test for two random choices of
a
, the chances are better
than 3 out of 4 that
n
is prime. By running the test with more and
more randomly chosen values of
a
we can make the probability of
error as small as we like.

The existence of tests for which one can prove that the chance of
error becomes arbitrarily small has sparked interest in algorithms of
this type, which have come to be known as
probabilistic
algorithms
. There is a great deal of research activity in this area,
and probabilistic algorithms have been fruitfully applied to many
fields.
⁴⁸

Exercise 1.21.
  Use the
smallest-divisor
procedure to find the smallest divisor
of each of the following numbers: 199, 1999, 19999.

Exercise 1.22.
  
Most Lisp implementations include a primitive called
runtime
that returns an integer that specifies the amount of time the system
has been running (measured, for example, in microseconds). The
following
timed-prime-test
procedure, when called with an
integer
n
, prints
n
and checks to see if
n
is prime. If
n
is
prime, the procedure prints three asterisks followed by the amount of time
used in performing the test.

(define (timed-prime-test n)
  (newline)
  (display n)
  (start-prime-test n (runtime)))
(define (start-prime-test n start-time)
  (if (prime? n)
      (report-prime (- (runtime) start-time))))
(define (report-prime elapsed-time)
  (display " *** ")
  (display elapsed-time))