{"id":246235,"date":"2025-05-28T04:37:33","date_gmt":"2025-05-28T09:37:33","guid":{"rendered":"https:\/\/www.johndcook.com\/blog\/?p=246235"},"modified":"2025-05-28T04:38:13","modified_gmt":"2025-05-28T09:38:13","slug":"lucas-lehmer-test","status":"publish","type":"post","link":"https:\/\/www.johndcook.com\/blog\/2025\/05\/28\/lucas-lehmer-test\/","title":{"rendered":"The bad version of a good test"},"content":{"rendered":"

Ever since 1952, the largest known primes have all had the form 2n<\/em><\/sup> \u2212 1, with one exception in 1989. The reason the largest known primes have this form is that it is easier to test whether these numbers are prime than other numbers.<\/p>\n

A number of the form 2n<\/em><\/sup> \u2212 1 is called a Mersenne number, denoted M<\/em>n<\/em><\/sub>. A Mersenne prime is a Mersenne number that is also prime. There is a theorem that says that if M<\/em>n<\/em><\/sub>\u00a0is prime then n<\/em> must be prime.<\/p>\n

Lehmer’s theorem<\/h2>\n

Derrick Henry Lehmer proved in 1930 that M<\/em>p<\/em> <\/sub>is prime if\u00a0p<\/em> is an odd prime and M<\/em>p<\/em><\/sub> divides\u00a0s<\/em>p<\/em>\u22122<\/sub> where the\u00a0s<\/em> terms are defined recursively as follows. Define s<\/em>0<\/sub> = 4 and s<\/em>n<\/em><\/sub> = s<\/em>n<\/em>\u22121<\/sub>\u00b2 \u2212 2 for n<\/em> > 0. This is the Lucas-Lehmer primality test for Mersenne numbers, the test alluded to above.<\/p>\n

Let’s try this out on M<\/em>7<\/sub> = 27<\/sup> \u2212 1 = 127. We need to test whether 127 divides s<\/em>5<\/sub>. The first six values of s<\/em>n<\/em><\/sub> are<\/p>\n

4, 14, 194, 37634, 1416317954, 2005956546822746114<\/p>\n

and indeed 127 does divide 2005956546822746114. As you may have noticed, s<\/em>5<\/sub> is a big number, and s<\/em>6<\/sub> will be a lot bigger. It doesn’t look like Lehmer’s theorem is going to be very useful.<\/p>\n

The missing piece<\/h2>\n

Here’s the idea that makes the Lucas-Lehmer [1] test practical. We don’t need to know \u00a0s<\/em>p<\/em>\u22122<\/sub> per se; we only need to know whether it is divisible by M<\/em>p<\/em><\/sub>. So we can carry out all our calculations mod M<\/em>p<\/em><\/sub>. In our example, we only need to calculate the values of s<\/em>n<\/em><\/sub> mod 127:<\/p>\n

4, 14, 67, 42, 111, 0<\/p>\n

Much better.<\/p>\n

Lucas-Lehmer test takes the remainder by M<\/em>p<\/em><\/sub> at every step along the way when computing the recursion for the\u00a0s<\/em> terms. The naive version does not, but computes the s<\/em> terms directly.<\/p>\n

To make the two versions of the test explicit, here are Python implementations.<\/p>\n

Naive code<\/h3>\n
def lucas_lehmer_naive(p):\r\n    M = 2**p - 1\r\n    s = 4\r\n    for _ in range(p-2):\r\n        s = (s*s - 2) \r\n    return s % M  == 0\r\n<\/pre>\n
Good code<\/h3>\n
def lucas_lehmer(p):\r\n    M = 2**p - 1\r\n    s = 4\r\n    for _ in range(p-2):\r\n        s = (s*s - 2) % M\r\n    return s == 0\r\n<\/pre>\n
The bad version<\/h2>\n
How bad is the naive version of the Lucas-Lehmer test?<\/p>\n
The OEIS sequence A003010<\/a> consists of the s<\/em>n<\/em><\/sub>. OEIS gives the following formula:<\/p>\n
\"s_n<\/p>\n
This shows that the s<\/em>n<\/em><\/sub> grow extremely rapidly, which suggests the naive algorithm will hit a wall fairly quickly<\/p>\n
When I used the two versions of the test above to verify that the first few Mersenne primes M<\/em>p<\/em><\/sub> really are prime, the naive version gets stuck after p<\/em> = 19. The better version immediately works through\u00a0p<\/em> = 9689 before slowing down noticeably.<\/p>\n
Historical example<\/h3>\n
Let’s look at M<\/em>521<\/sub>. This was the first Mersenne number verified by a computer to be prime. This was done on the vacuum tube-based SWAC computer in 1952 using the smart version of the Lucas-Lehmer test.<\/p>\n
Carrying out the Lucas-Lehmer test to show that this number is prime requires doing modular arithmetic with 521-bit numbers, which was well within the 9472 bits of memory in the vacuum tube-based SWAC computer that proved M<\/em>521<\/sub> was prime.<\/p>\n
But it would not be possible now, or ever, to verify that M<\/em>521<\/sub> is prime using the naive version of the Lucas-Lehmer test because we could not store, much less compute, s<\/em>519<\/sub>.<\/p>\n
Since<\/p>\n
\"\\log_{2}<\/p>\n
we’d need around 2520<\/sup> bits to store s<\/em>519<\/sub>. The number of bits we’d need is itself a 157-digit number. We’d need more bits than there are particles in the universe, which has been estimated to be on the order of 1080<\/sup>.<\/p>\n
The largest Mersenne prime that the SWAC could verify using the naive Lucas-Lehmer test would have been M<\/em>13 <\/sub>= 8191, which was found to be prime in the Middle Ages.<\/p>\n
Related posts<\/h2>\n