Large Prime Numbers

Date: 12/17/97 at 13:55:22
From: Lynne & Dave Ware
Subject: Prime Numbers

To the Math Swat Team:

I believe that I read somewhere (Ivars Pederson [spelling?], I think)
that it is possible to determine if a very large number is prime using 
a simple procedure. (But not to determine its factors). It didn't say 
how this is done. If this is true, what is the algorithm?

Thanks for your consideration.

Dave Ware

Date: 12/17/97 at 15:50:01
From: Doctor Rob
Subject: Re: Prime Numbers

Good question! The person is probably Ivars Peterson, who is a popular
science writer.

The situation is a bit more complicated than you remember. The fact is
that the simple procedure does not provide a proof of primality, but 
may provide a proof of compositeness. It is based on Fermat's Little 
Theorem, which says:

   Theorem: If p is a prime number, a is an integer, and p is not
            a divisor of a, then p is a divisor of a^(p-1) - 1.

Given a number N we want to test, pick any old a, and find the 
greatest common divisor of a and N, GCD(a,N). If this is not 1, then 
it is a factor of N, and N is composite. If it is 1, then compute the 
remainder r of a^(N-1) when divided by N. If r is not 1, then by 
Fermat's Little Theorem, N cannot be prime, so is composite. This 
gives the proof of compositeness.

By the way, you might as well choose 1 < a < N-1, since a and a-N 
will give the same value of r, and since a = 1 or N-1 will always 
give r = 1.

How to compute the remainder r?  If N is even moderately large, 
computing a^(N-1) and then dividing by N will be a bad idea, since the 
number a^(N-1) will have many, many digits. The trick is to divide by 
N and keep only the remainder at all intermediate steps. It may not be 
obvious that this works, but it does. If N = 67, N-1 = 66, you might 
compute a^66 by doing 65 multiplications. After each one, divide by 67 
and keep only the remainder. 

Better than doing a^66 by 65 multiplications (and 65 divisions by N), 
you can shortcut the computation by the following trick:
a^66 = (a^33)^2, a^33 = a*a^32, a^32 = (a^16)^2, a^16 = (a^8)^2,
a^8 = (a^4)^2, a^4 = (a^2)^2, a^2 = a*a. Working from the last 
equation backwards, you will need only 7 multiplications (and 
7 divisions by N).

If r = 1, then what can we say?  All prime numbers will give r = 1, 
but there are a few composite numbers which will do so, too. For 
example, if a = 2, N = 341 = 11*31 is composite, but r = 1. Such N's 
are called Fermat pseudoprimes with respect to base a. 341 is a Fermat 
pseudoprime with respect to base 2. 

It turns out that Fermat pseudoprimes with respect to any fixed base 
are uncommon.  The chance of picking one at random from some large set 
of N's is very small. Thus, in some sense, which is again a 
complicated matter, if you get r = 1, N is "probably" prime.

A tactic you could use if you get r = 1 is to pick another base a and
repeat the calculation.  If you get r = 1 again, try still another a.
Continue this until you either find a proof that N is composite, or 
you decide that you couldn't possibly be unlucky enough to have an N 
which is a Fermat pseudoprime to all the bases you used.  Then declare 
the number to be prime, and you will be wrong only rarely.

A flaw with the above tactic is that there exist numbers called
Carmichael numbers, for which, for all bases a such that N passes the
GCD test, r will equal 1. The smallest one is 561 = 3*11*17. Every 
base a not divisible by 3, 11, or 17 will give r = 1. These are even 
rarer than Fermat pseudoprimes with respect to a given base a, but 
there are known to be infinitely many. If you happened to choose one 
of these, you might try very many bases, getting r = 1 over and over, 
then declare the number prime, and be wrong.

There are variations of the above method which get rid of the last 
flaw, give you proofs of compositeness, but only probabilistic 
statements about primality.  If you need to know more about them, 
write again.

If you want a true proof of primality, there are more complicated 
methods which will do this, but I am quite sure the above is what you 
remember reading about. The complicated methods can be proven to run 
on a computer in a relatively short time, and produce a certificate of 
primality which can be checked even faster.  They could not, however, 
be termed "simple" in most senses of the word!

-Doctor Rob,  The Math Forum
 Check out our web site!  http://mathforum.org/dr.math/   

Date: 12/18/97 at 10:36:38
From: Lynne & Dave Ware
Subject: Follow up question on testing for primes

Doctor Rob:

Thank you for your quick answer to my question.

I think I was looking for the true test of primality you mentioned in
your last paragraph. I think Ivars Peterson did cover some of the 
theory on Fermat's little theorem. I also remember something about a 
method he only mentioned that could be run on a computer. How can I 
program my computer to do this complicated method? 

I am working on a prime number program (in Pascal) that uses the sieve 
method. I saw an answer on your web site about testing for primes. I 
go further in simplification by only testing with prime numbers less 
than the square of the numbers and only odd numbers. I want to be able 
to cross check on its operation by using a check for primes. The 
"simple method" looks too messy for large numbers as I am trying to 
push my program to work up to 10 billion.

Dave Ware

Date: 12/18/97 at 13:22:02
From: Doctor Rob
Subject: Re: Follow up question on testing for primes

If you think the "simple method" is too messy for large numbers, then 
you will really hate the more sophisticated ones. They are too 
complicated to go into via e-mail, but I will point you to a 
reference. See

Henri Cohen, _A Course in Computational Algebraic Number Theory_,
  Springer-Verlag, Graduate Texts in Mathematics 138, Chapter 9.

You should be able to find a copy in any medium-to-large university
library, and probably other places.  You will need a background in 
either algebraic number theory, or the theory of elliptic curves, to 
understand why the algorithms work, but only a knowledge of 
programming to implement them.

Here is another idea. Since your range is only up to 10^10 (as opposed 
to 10^1000[!!]), here is an alternative. Use the Strong Compositeness 
Test (below) with bases 2, 3, 5, and 7. If it fails all of these, and 
is not equal to 3215031751, then it is prime.

The Strong Compositeness Test with Base a: To test a number N for
compositeness:

1. If GCD(a,N) > 1, stop and declare that N is composite.
2. Write N - 1 = 2^s*d, where d is odd.
3. Compute R, the remainder of a^d when divided by N.  (Do this as
   described in the last message.)
4. If R = 1, stop and declare failure.
5. For i = 0, 1, ..., s-1, do the following:
   a. If R = N-1, stop and declare failure.
   b. Replace R with the remainder of R^2 when divided by N.
   c. If R = 1, stop and declare that N is composite.
6. Stop and declare N composite.

A number which fails the Strong Compositeness Test with Base a, yet is
composite, is called a "strong pseudoprime with respect to base a."  
Every strong pseudoprime is a Fermat pseudoprime to the same base, but 
the reverse is false.  The smallest strong pseudoprime with respect to 
base 2 is 2047.

   Theorem: The only number less than 2.5*10^10 which is a strong 
            pseudoprime with respect to bases 2, 3, 5, and 7, is 
            3215031751.

-Doctor Rob,  The Math Forum
 Check out our web site!  http://mathforum.org/dr.math/