2

This piqued my interest as I was trying to find a larger prime number, but I soon realised that the interpreter of my programming language soon popped an error when trying to insert a number with 24 million digits (2^82,589,933 − 1)to a variable. So I wondered, how do researchers manage to do this? How do they go through the process of finding numbers this large when there are clear computational limitations?

  • 3
    Maybe this question will be more suitable in [Math Stackexchange](https://math.stackexchange.com/). – Yonlif Jun 13 '20 at 20:54
  • 4
    If the question is what algorithms are used to find or test large (possible) primes, that belongs on math. If the question is how do they (or anyone) write programs to implement such algorithms on large numbers, that belongs here and is very basic: some languages have bignums builtin like python, some have them in the standard library like Java, and some need addon libraries like C and C++. See tag [bignum] and https://en.wikipedia.org/wiki/Arbitrary-precision_arithmetic . If you have a question about a specific algorithm and implementation, that would be answerable. – dave_thompson_085 Jun 13 '20 at 21:33
  • 6
    As you can verify from https://en.wikipedia.org/wiki/Largest_known_prime_number the large primes are of special forms that there are special primality tests for. Most have been found by https://www.mersenne.org/. The math for that search is explained at https://www.mersenne.org/various/math.php. – btilly Jun 14 '20 at 00:03
  • https://en.wikipedia.org/wiki/Lucas%E2%80%93Lehmer_primality_test – President James K. Polk Jun 14 '20 at 23:02
  • I converted the LL test to a prime finder. I spend two minutes on this, so thought it was cool to share, so the latest comments on my post to see that it can find primes, even on not so small numbers. Not groundbreaking or earth shattering, just cool since i made a small change just to see what happens and whoah, it finds primes. – oppressionslayer Jun 18 '20 at 05:53

1 Answers1

2

A pow test of 65537 get's you the Mersenne numbers up to it's length at least, i'm not sure what test they use as i'm only able to climb up to 2**4423-1 after a couple of minutes:

for x in range(2,100000): 
   if pow(65537,2**x-2, 2**x-1) == 1: 
       print(f"2**{x}-1 is prime") 

2**2-1 is prime
2**3-1 is prime
2**5-1 is prime
2**7-1 is prime
2**13-1 is prime
2**17-1 is prime
2**19-1 is prime
2**31-1 is prime
2**61-1 is prime
2**89-1 is prime
2**107-1 is prime
2**127-1 is prime
2**521-1 is prime
2**607-1 is prime
2**1279-1 is prime
2**2203-1 is prime
2**2281-1 is prime
2**3217-1 is prime
2**4253-1 is prime
2**4423-1 is prime

Based on President James K. Polk's comment i wrote the Lucas-Lehmer Primality test


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


# 20th Mersenne Prime
In [401]: LucasLehmer(4423)                                                                                                                                                                                 
Out[401]: True

# 21st Mersenne Prime
In [402]: LucasLehmer(9689)                                                                                                                                                                                 
Out[402]: True

# Some non random test:
In [404]: LucasLehmer(2727)                                                                                                                                                                                 
Out[404]: False

# 24th Merseene Prime
In [409]: LucasLehmer(19937)                                                                                                                                                                                
Out[409]: True

I played around with the math in the LL test and came up with a nice prime finder. Slow, not groundbreaking, but not bad:

from sympy import isprime

def PrimeFinderLucasLehmer(N):
   p = 1<<N.bit_length()-1
   if p == 2:
     return True
   s = 4
   M = pow(p, 2) - 1
   for x in range (1, (p-2)+1):
     s = (((s * N) - 2  )) % M
     mgcd = math.gcd(s, N)
     if  mgcd != 1:
        print(s, mgcd)
        if isprime(mgcd) == True:
           break
        else: continue
   return mgcd

From Cracking short RSA keys I was able to get the primes for those numbers rather quickly:

In [149]: PrimeFinderLucasLehmer(8114231289041741)                                                                                                             
15823901318551394674372071332152 1839221
Out[149]: 1839221

In [148]: PrimeFinderLucasLehmer(10142789312725007)                                                                                                            
36241961070961173943814389655700 100711423
Out[148]: 100711423

In [151]: PrimeFinderLucasLehmer(1009732533765203)                                                                                                             
157641275043535187903371634694 1823

Anyways, i thought that was interesting, so wanted to share. That was just with a tweak of the math.

oppressionslayer
  • 6,942
  • 2
  • 7
  • 24
  • 1
    I modified this to be much faster check out: https://github.com/oppressionslayer/primalitytest/blob/master/findprimell.py to see a factorization engine built around Lucas-Lehmer primality test. This isn't earth shattering or groundbreaking, but it's cool, as it can do numbers from the stackoverflow site i posted in seconds. I hope you enjoy, it was a nice find and works great on numbers about that size and under and of course bigger numbers as seen in the comments of the code – oppressionslayer Jun 19 '20 at 05:15