2

I am trying to write a function that should take in a Mersenne number and return whether it is a prime number or not using the Lucas-Lehmer primality test. I am trying to return the last number generated in the Lucas-Lehmer sequence which should be 0 if it is a prime number

Lucas Lehmer sequence

I have written the following function to do the above

def lucas_lehmer_modified(p):
   M=2**p-1
   for i in range (0,p-1):
     if i == 0:
       x = 4
     else : 
        x = (prev**2-2)%(M)
     prev=x
   return x

My problem is that this code works for small numbers such as 127 but not for large numbers like 2305843009213693951 and even for 524287. My Python Jupyter notebook hangs up. Any advice on how I can get a function that takes in a Mersenne prime as an input and returns whether it is a prime number or not using the Lucas Lehmer test. I need to make it work for at least 2^65-1

Paul Brit
  • 5,901
  • 4
  • 22
  • 23
  • Where do you run the Jupyter notebook? Locally or on Google / Internet? – Thomas Weller Apr 22 '20 at 07:06
  • 2
    Unrelated, but according to [wikipedia](https://en.wikipedia.org/wiki/Lucas%E2%80%93Lehmer_primality_test) the LLT applies for Mersenne number of the form `2**p - 1` **where p is an odd prime**. 65 is odd but not prime (5*13). – Serge Ballesta Apr 22 '20 at 07:08
  • What do you actually pass to `lucas_lehmer_modified()`? Do you pass `61` or do you pass `2305843009213693951`? – Thomas Weller Apr 22 '20 at 07:14
  • It works well for me, if I pass `19`, which is the p for `524287` and it also works well for `61` which is the p for `2305843009213693951`. My feeling is that you pass M to the function instead of p. – Thomas Weller Apr 22 '20 at 07:22
  • For p=`19937`, the function already needs 18 seconds to run. Any higher numbers than that will likely cause a hang, so that's another indicator for me that you pass M instead of p. – Thomas Weller Apr 22 '20 at 07:26

1 Answers1

0

I wrote a Lucas Lehmer test from the psuedocode from the wikipedia entry, here is what i came up with ( p being the mersenne power number ):

def LucasLehmer(p):
   if p == 2:
     return True
   s = 4
   M = pow(2, p) - 1
   for x in range (1, (p-2)+1):
      s = ((s * s) - 2) % M
   if s == 0: return True
   else: return False

And some answers:

In [388]: LucasLehmer(9689)                                                                                                                                                                     
Out[388]: True

In [389]: LucasLehmer(19937)                                                                                                                                                                    
Out[389]: True

In [390]: LucasLehmer(500)                                                                                                                                                                      
Out[390]: False



oppressionslayer
  • 6,942
  • 2
  • 7
  • 24