## Mersenne Primes

### June 3, 2011

Here is our version of the Lucas-Lehmer test for Mersenne primes:

`(define (mersenne-prime? p)`

(let ((m (- (expt 2 p) 1)))

(if (not (prime? m)) #f

(let loop ((i 3) (s 4))

(if (< p i) (zero? s)

(loop (+ i 1)

(modulo (- (* s s) 2) m)))))))

We made two changes to the algorithm stated on the previous page: First, the S sequence is accumulated modulo 2^{p}−1, which keeps the numbers smaller. Second, the divisibility test has become a test for zero-congruence modulo 2^{p}−1, which is necessary because we no longer have the true value of S. The indices of the Mersenne primes less than 256 are shown below; we had to handle M_{2} specially, since the Lucas-Lehmer test specifies only odd primes:

`> (cons 2 (filter mersenne-prime? (range 256)))`

(2 3 5 7 13 17 19 31 61 89 107 127)

We used range and filter from the Standard Prelude, and `prime?`

from a previous exercise. You can run the program at http://programmingpraxis.codepad.org/VEUVtatQ.

Pages: 1 2

My Haskell solution (see http://bonsaicode.wordpress.com/2011/06/03/programming-praxis-mersenne-primes/ for a version with comments):

Perhaps it’s too early for me to try these without coffee, but I think there may be some mistakes. M_256 = 2^256 – 1, while not prime, is a 78 digit number. Moreover, 2 is not one less than a power of 2, so it’s not a Mersenne Prime.

@Graham: For your first point, both version give the exponents. The actual primes are trivially obtained by calculating 2^p – 1. For your second point, since we’re calculating the exponents, the actual value is 2^2 – 1 = 3, which is prime.

@Remco: Thanks, but it seems I’m still not getting it. For the second point, yes, 2^2-1 = 3 is prime, but that’s M_2. However, there is no n such that M_n = 2, since 2^n – 1 is an odd number. Hence, 2 is not a Mersenne Number, and therefore not a Mersenne Prime. I’m not quite sure what you mean for the first point. Basically, it seems to me that the question as written is not “Find all primes p less than 256 such that 2^p – 1 is prime,” but “Find all numbers of the form M_n = 2^n – 1 where M_n is prime for 2 < n < 256."

@Graham: From the wikipedia article on Mersenne primes: “Some definitions of Mersenne numbers require that the exponent p be prime, since the associated number must be composite if p is.” Testing only primes is an optimization since if p is not prime, 2^p – 1 won’t be. The definitions “Find all primes p less than 256 such that 2^p – 1 is prime” and “Find all numbers of the form M_n = 2^n – 1 where M_n is prime for 2 < n < 256." are therefore equivalent, save for the fact that the former returns p and the second 2^p – 1, which is the trivial transformation I mentioned.

@Remco: Thanks for the point about only testing primes from 2 to 256 instead of all integers. I hadn’t seen that. Appreciate you setting me straight! (and sorry for being so obtuse)

Graham: If you were confused about definitions, that’s my fault, since I set the terms of the exercise.

Mersenne numbers are numbers of the form 2^n-1, with n>0: 1, 3, 7, 15, 31, …. (Sometimes the sequence is extended to n=0, so the first term is 0, but we’ll stick with the standard definition.) A Mersenne number is either prime or composite; thus, the Mersenne number 7 is prime, and the Mersenne number 15=3*5 is composite.

Since Mersenne numbers grow so quickly, they are generally referred to by their index. For instance, 162259276829213363391578010288127 is a Mersenne number, but it’s easier to say M107. And it’s also prime, so it’s a Mersenne prime. And, as a special added bonus, we also know that since M107 is prime, so is 107, since it is a necessary (but not sufficient) condition that the index of a Mersenne prime is also a prime.

Sometimes when mathematicians talk about Mersenne numbers they really mean Mersenne numbers such that the index is prime, whether or not the number itself is prime, which is a ready source of confusion; that sequence is 3, 7, 31, 127, 2047, 8191, …, and it’s usually written 2^p-1 instead of 2^n-1. You may notice that Wikipedia uses the first definition, while Sloane uses the second. There is no “official” definition, you just have to sort out for yourself which meaning is intended.

The indices of the Mersenne primes are 2, 3, 5, 7, 13, 17, 19, 31, 61, 89, 107, 127, … and the Mersenne primes themselves are 3, 7, 31, 127, 8191, 131071, 524287, 2147483647, 2305843009213693951, 618970019642690137449562111, 162259276829213363391578010288127, 170141183460469231731687303715884105727, …; notice that M11=2047 is not a Mersenne prime, even though it is a Mersenne number under either definition.

If you want to know more about Mersenne numbers, of all kinds, a good place to start is OEIS (which used to be called “Sloane’s” after its founder N. J. A. Sloane, but he is now retired); just type “mersenne” in the search bar and click “Search.” After you’re finished being amazed, you can read some of the entries there, and follow the references given in each entry, and soon enough your entire evening will be gone. Sloane’s is one of my favorite web sites to sit and browse.

The program you were asked to write was to find the indices of the Mersenne primes to M256. The solution first tests the index, which must be prime, and then performs a Lucas-Lehmer test to determine if the Mersenne number is also prime. The list of indices and the corresponding numbers are shown in the prior paragraph.

Hope this helps,

Phil

Thanks Phil, I must have just been particularly slow yesterday. I appreciate the background, though, and look forward to learning more. Interestingly, the math software I use (Sage) has an interface to Sloane’s built in.

Thanks again to you and Remco for putting up with me and setting me straight, and thank you for this great blog!

Graham

I reused a revised version of the “Python Cookbook’s” Sieve of Eratosthenes recipe for

`primes()`

. Also, my original code computed the function`s`

on its own, but the code was incredibly slow (lots of recursive calls with very large numbers?). I opted for a translation of the Scheme solution.Factor code. the sequence S1, S2, etc. is a lazily computed list. primes word is sieve of eratosthenes from previous exercise.

( scratchpad ) 256 mersennes .

V{ 2 3 5 7 13 17 19 31 61 89 107 127 }

I managed to find 2^44497 – 1 was prime before giving up :)

Revisit in Clojure…