Mersenne primes are prime numbers of the form 2^n – 1. They were named after Marin Mersenne, a 17th century French monk who compiled a table of the first dozen (although it turns out he got the last few wrong).

The first few are:

3 = 2^2 – 1

7 = 2^3 – 1

31 = 2^5 – 1

127 = 2^7 – 1

8191 = 2^13 – 1

131071 = 2^17 – 1

…

**Lucas-Lehmer**

Because the numbers grow large rapidly, simple primality tests such as trial division are not feasible for large values of n. The Lucas-Lehmer primality test, developed by mathematicians Edouard Lucas and Derrick Lehmer, provides a straight-forward means for quickly testing large primes.

For testing the primarity of p = 2^n – 1, they define the sequence S(i) as follows:

S(0) = 4

S(i) = ((S(i-1))^2 – 2) mod p

Then S(n-2) = 0 if and only if p = 2^n – 1 is an odd prime.

Because arithmetic is conducted modulo p, the size of terms in S(i) remains relatively small, i.e. less than p.

For more details, Terry Tao has posted an excellent description of the Lucas-Lehmer test and the proof for its correctness.

**Fast modulus**

Traditionally, modulus is one of the slower arithmetic operations, especially on numbers with potentially thousands of digits. However, it turns out that in this case, because of the relationship between Mersenne primes and binary representation of numbers in the computer, we can implement the modulus part of the equation wholly in terms of shifts and adds.

Here’s a quick base 10 example. Suppose we want to compute 1234567 mod 10^2 – 1. Now we can break this number into a sum of products of 100.

1234567 mod 99 = (1*100*100*100 + 23*100*100 + 45*100 + 67) mod 99

A basic property of modulus is that (a*b) mod n = a mod n * b mod n. For this example, this means we can mod all of the 100 terms by 99 and get 1. For example, (1*100*100*100) mod 99 = 1 mod 99 * 100 mod 99 * 100 mod 99 * 100 mod 99 = 1 * 1 * 1 * 1 = 1.

Applying this to the example above:

(1*100*100*100 + 23*100*100 + 45*100 + 67) mod 99 = (1*1*1*1 + 23*1*1 + 45*1 + 67) mod 99

= (1 + 23 + 45 + 67) mod 99

Effectively, dropping of 100 terms equates to shifting the digits to the right.

= 136 mod 99.

We can then repeat the process:

136 mod 99 = (1*100 + 36) mod 99 = 1*100 mod 99 + 36 mod 99 = 37 mod 99 = 37.

In binary, the process is the same, except that instead of base 10, we’re working in base 2. Consider 1234 mod (2^3 – 1). 1234 = 10011010010 and 7 = 111 in binary.

10011010010 mod 111 = (1*1000*1000*1000 + 001*1000*1000 + 010*1000 + 010) mod 111

1000 = (111 + 1), so 1000 mod 111 is 1, thus like in the decimal example, the intermediate product terms drop out.

10011010010 mod 11 = (10*1000*1000*1000 + 011*1000*1000 + 010*1000 + 010) mod 111

= (10*1*1*1 + 11*1*1 + 10*1 + 10) mod 111

= (10 + 11 + 10 + 10) mod 111

= (1001) mod 111

= (1*1000 + 1) mod 111

= (1*1 + 1) mod 111

= (1 + 1) mod 111

= 10 mod 111

= 10

10 is 2 in decimal, so we find that 1234 mod 7 = 5

In practice, it may be easier to deal with one term at a time, starting at the left. The process converges reasonably quickly, and you don’t have as many temporary results to keep track of. Using the above example, we would do:

10011010010 mod 111 = ((10*1000*1000*1000) + 011010010) mod 111

= (10 + 011010010) mod 111

= (011010100) mod 111

= (011*1000*1000 + 010100) mod 111

= (11*1*1 + 010100) mod 111

= (11 + 010100) mod 111

= (010111) mod 111

= (010 * 1000 + 111) mod 111

= (10 + 111) mod 111

= (1001) mod 111

= (1*1000 + 001) mod 111

= (1 + 1) mod 111

= 10 mod 111

Again, this amounts to shifting and adding. The only place where this method doesn’t work is where the you have n mod n, but this case can easily be checked.

**Multiplication**

Multiplication ends up being the main bottleneck in the Lucas-Lehmer test. The standard algorithm for multiplication require multiplying each digit of one operand by each digit of the other, and adding up intermediate products. If you’re multiplying two n-digits numbers, you have to perform n^2 digit multiplications. More complicated methods such as Karatsuba multiplication and Fast Fourier Transform multiplication reduce the number of multiplications at the expense of increasing other operations.

More on that in future postings.

**Example**

gmp, the GNU Multiple Precision arithmetic package, takes care of arithmetic operations on arbitrarily large numbers.