A few weeks ago I heard someone casually refer to ‘that formula of Euler’s that generates primes’. I hadn’t heard of this, but it turns out that in 1772 Euler produced this formula:
\[ f(x) = x^2 + x + 41\text{.} \]
Using this, \(f(0)=41\), which is prime. \(f(1)=43\), which is also prime. \(f(2)=47\) is another prime. In fact this sequence of primes continues for an incredible forty integer inputs until \(f(40)=41^2\). It might generate more primes for higher inputs, but what’s interesting here is the uninterrupted sequence of forty primes.
This got me wondering. Clearly \(f(0)\) is prime because 41 is prime, so that much will work for any function
\[ f(x) = x^2 + x + p \]
for prime \(p\), since \(f(0)=0^2+0+p=p\). Are there other values of \(p\) that generate a sequence of primes? Are there any values of \(p\) that generate longer sequences of primes?
I wrote some code to investigate this. Lately, I’ve taken to writing C++ when I need a bit of code, for practice, so I wrote this in C++.
I figured the cases where \(f(0)\) is prime but \(f(1)\) isn’t weren’t that interesting, since \(f(0)\) is trivially prime. In fact, \(f(x)=x g(x)+p=p\) when \(x=0\) for any prime \(p\), but saying so doesn’t seem worth the effort.
So I kept track of the primes \(p\) whose functions \(f(x)=x^2+x+p\) generate more than one prime, and the lengths of the sequences of primes generated by each of these. This produced a pair of integer sequences.
I put the primes that work into the OEIS and saw that I had generated a list of the smaller twin in each pair of twin primes. I was momentarily spooked by this, until I realised it was obvious. Since \(f(0)=p\) and \(f(1)=1^2+1+p=p+2\), any prime this works for will generate at least a twin prime pair \(p,p+2\).
What about the lengths of the sequences of consecutive primes generated? The table below shows the sequences of consecutive primes generated for small values of \(p\). Most primes that generate a sequence produce just two, and \(p=41\) definitely stands out by generating forty.
I was pleased to see this sequence of lengths of primes generated was not in the OEIS. So I submitted it, and it is now, along with the code I wrote. (I discovered along the way that the version where sequences of length one are included was already in the database.)
Anyway, I amused myself by having some C++ code published, and by citing Euler in a mathematical work. Enjoy: A371896.
In the last Finite Group livestream, Katie told us about emirps. If a number p is prime, and reversing its digits is also prime, the reversal is an emirp (‘prime’ backwards, geddit?).
For example, 13, 3541 and 9999713 are prime. Reversing their digits we get the primes 31, 1453 and 3179999, so these are all emirps. It doesn’t work for all primes – for example, 19 is prime, but 91 is \(7 \times 13 \).
In the livestream chat the concept of primemirp emerged. This would be a concatenation of a prime with its emirp. There’s a niggle here: just like in the word ‘primemirp’ the ‘e’ is both the end of ‘prime’ and the start of ’emirp’, so too in the number the middle digit is end of the prime and the start of its emirp.
Why? Say the digits of a prime number are \( a_1 a_2 \dots a_n \), and its reversal \( a_n \dots a_2 a_1 \) is also a prime. Then the straight concatenation would be \( a_1 a_2 \dots a_n a_n \dots a_2 a_1 \). Each number \(a_i\) is in an even numbered place and an odd numbered place. Now, since
\[ 10^k \pmod{11} = \begin{cases} 10, & \text{if } k \text{ is even;}\\ 1, & \text{otherwise,} \end{cases} \]
it follows that each \(a_i \) contributes a multiple of eleven to the concatenation. A mismatched central digit breaks this pattern, allowing for the possibility of a prime.
I wrote some code to search for primemirps by finding primes, reversing them and checking whether they were emirps, then concatenating them and checking the concatenation. I found a few! Then I did what is perfectly natural to do when a sequence of integers appears in front of you – I put it into the OEIS search box.
Imagine my surprise to learn that the concept exists and is already included in the OEIS! It was added by Patrick De Geest in February 2000, based on an idea from G. L. Honaker, Jr. But there was no program code to find these primes and only the first 32 examples were given. I edited the entry to include a Python program to search for primemirps and added entries up to the 8,668th, which I believe is all primemirps where the underlying prime is less than ten million. My edits to the entry just went live at A054218: Palindromic primes of the form ‘primemirp’.
Here’s a game I’ve been trying to make for a while.
For a while I’ve had a hunch that there’s fun to be had in moving between numbers by using something related to the prime numbers.
Over the years I’ve tried out a few different ideas, but none of them ever worked out – they were either too easy, too hard, or just not interesting. This time, I think I’ve found something close enough to the sweet spot that I’m happy to publish it.
Prime Run is a game about adding and subtracting prime numbers. You start at a random number, with a random target. Your goal is to reach the target, by adding or removing any prime factor of your current number.
Following on from the series of ‘Pascal’s Triangle and its Secrets‘ posts, guest author David Benjamin shares another delightful piece of mathematics – this time relating to prime numbers.
At the time of writing the largest known prime number has $24862048$ digits. The number of digits does not reflect the true size of this prime but if we were to type it out at Times New Roman font size 12, it would reach approximately $51.5$ km, or about $32$ miles. Astonishing!
Patrick Laroche from Ocala, Florida discovered this Mersenne prime on December 7, 2018. I was surprised to discover that it’s exponent $82589933$ is the length of the hypotenuse of a primitive Pythagorean triple where $82589933^{2} = 30120165^{2} + 76901708^{2}$ as indeed are 8 of the exponents of those currently ranked from 1 to 10.
The Greek mathematician Euclid of Alexandria ($\sim$325 BC-265 BC) was arguably the first to prove that there are an infinite number of primes – and since then, people have been searching for new ones. Some do it for kudos, for the prize money, to test the power of computers and the need to find more of the large primes used to help protect the massive amount of data which is being moved around the internet.
Mersenne primes, named after the French monk Marin Mersenne, are of the form $2^{p} -1$, where the exponent $p$ is also prime. Mersenne primes are easier to test for primality, which is one reason we find so many large ones (all but one of the top ten known primes are Mersenne). When Mersenne primes are converted to binary they become a string of $1$s, which makes them suitable for computer algorithms and an excellent starting point for any search.
Marin Mersenne
Since generally testing numbers for primality is slow, some have tried to find methods to produce primes using a formula. Euler’s quadratic polynomial $n^2+n+41$ produces this set of $40$ primes for $n = 0$ to $39$. When $n=40$, the polynomial produces the square number $1681$. Other prime-generating polynomials are listed in this Wolfram Mathworld entry.
The French mathematician LejeuneDirichlet proved that the linear polynomial $a+nb$ will produce an infinite set of primes if $a$ and $b$ are coprime for $n=0,1,2,3,4,…$. Then again, it also produces an infinite number of composite numbers! However, this gem: $224584605939537911 + 1813569659748930n$ produces 27 consecutive primes for $n=0$ to $n=26$ – and of course, all the primes are in arithmetic progression.
14 fruitful fractions
The primes are unpredictable, and become less common as they get larger. Consequently there is no formula that will generate all the prime numbers. However, there is a finite sequence of fractions, that – given an infinite amount of time – would generate all the primes, and in sequential order.
They are the fruitful fractions, created by the brilliant Liverpool-born mathematician, John Horton Conway (1937–2020) who, until his untimely death from complications related to COVID-19, was the John von Neumann Emeritus Professor in Applied and Computational Mathematics at Princeton University, New Jersey, USA.
John Horton Conway (Photo: Denise Applewhite, Office of Communications)
The fruitful fractions are
$\frac{17}{91}$
$\frac{78}{85}$
$\frac{19}{51}$
$\frac{23}{38}$
$\frac{29}{33}$
$\frac{77}{29}$
$\frac{95}{23}$
$\frac{77}{19}$
$\frac{1}{17}$
$\frac{11}{13}$
$\frac{13}{11}$
$\frac{15}{44}$
$\frac{15}{2}$
$\frac{55}{1}$
A
B
C
D
E
F
G
H
I
J
K
L
M
N
The first time I encountered this set of fractions was in the wonderful book, The Book of Numbers, by Conway and Guy. I was so intrigued as to how Conway came up with his idea, I emailed him to ask. I was delighted to receive an outline of an explanation and even a second set of fractions, neither of which I can now find – it was 1996 and pre-cloud storage! But no worries… Conway explains everything in this lecture, which also demonstrates his passion for mathematics and his ability to express his ideas in a relaxed and humorous way, even when he searches for an error in his proof on 26 minutes. The lecture also includes an introduction to Conway’s computer language, FRACTRAN, which includes the statement:
‘It should now be obvious to you that you can write a one line fraction program that does almost anything, or one and a half lines if you want to be precise‘.
Using the fractions to find prime numbers
Here’s how the fractions are used to generate primes.
Start with the number $2$
Multiply by each of the fourteen fractions until you find one for which the product is an integer
Starting with this new integer, continue multiplying through the fractions until another integer is produced. (If this process reaches fraction $N=\frac{55}{1}$, the integer’s product with N is guaranteed to be another integer as N has a denominator of $1$; the process continues with this new integer being multiplied by fraction A)
Continue multiplying through the set to create more integers
When the integer is a power of $2$, its exponent will be a prime number.
The 19 steps needed to produce the first prime number are:
The successive primes are produced almost like magic – but the number of multiplications needed to produce each new prime becomes larger and larger, and so the method, though wonderfully inventive, is not at all efficient.
Edit: Since this article was first published, the exponent $82589933$ of the Laroche prime has been accepted as the next term in the sequence http://oeis.org/A112634
There are many sequences of numbers to be found in Pascal’s triangle. The Natural numbers occur in the second diagonal, running in either direction, and the next two diagonals after that contain other important sequences:
Sequences in the diagonals
There are many sequences of numbers to be found in Pascal’s triangle. The Natural numbers occur in the second diagonal, running in either direction, and the next two diagonals after that contain other important sequences:
This article originally attributed the last Mersenne prime found by hand to Leonhard Euler. It was Édouard Lucas. So, that’s what the title is about. We also only realised after publication that there’s a wordplay in the title that only works if you say Euler wrong. The proofreaders apologise.
There is an ongoing quest in mathematics to find larger and larger prime numbers. Since the primes go on forever, there is no single largest one, but every day around the world computers are churning to find a new record holder.
This project was inspired by a Numberphile video about how mathematicians find enormous prime numbers. The goal is to find a large prime number ourselves, one that took mathematicians over 2,000 years to find. The record for the largest known prime has since been broken since that video was released; the largest known prime as of May 2021 is
\[ M_{82,589,933} = 2^{82,589,933} – 1 \]
If you were to write the number down it would be longer than the entire text of the Harry Potter series. We won’t find any prime numbers nearly that big, but we can try to find one too big to fit in a single tweet.
I don’t believe anyone has ever read the largest prime number from beginning to end.
2,000 years of math in 10 lines of Python
Mathematicians have been trying to find progressively larger prime numbers for several thousand years – the University of Tennessee at Martin has a great introduction to the topic. While the current record holders are found with computers, these techniques were once done by hand. The last time a human set the record for calculating the largest known prime was done by Édouard Luca in 1876 when he proved that $M_{127}$ is prime. Lucas’ record stood for 76 years until 1952 when a computer finally found a larger prime number. Since then all records have been set with the help of computers. In 10 lines of Python code, we will prove Lucas was right and find a bigger prime, surpassing more than 2,000 years of mathematical history.
With pen and paper, Lucas showed that $M_{127}$ is prime.
$M_{82,589,933}$ was found by a worldwide computer search called the Greater Internet Mersenne Prime Search, known as GIMPS. GIMPS works by picking a number that might be prime and testing it using a massive network of volunteered computers. By breaking the rather difficult question “Is this number prime?” into smaller pieces each machine in the network can work on a chunk of the problem they can handle, then report its results. This massively parallelized structure allows their network of more than two million computers to tackle prime numbers far larger than any single machine could handle.
GIMPS runs several tests on each candidate number before running the most conclusive, and slowest, test to see if a candidate number is indeed prime. These initial screens, which include trial division and Pollard’s P-1 method, can quickly show if a candidate is composite, so they can eliminate candidates, but cannot prove a number is prime. Once a candidate passes all the early rounds it’s time for the conclusive test. Until 2018, their conclusive test was the Lucas-Lehmer Primality test. The Lucas-Lehmer test is a definitive method to prove whether or not a number is prime and the one we’ll be focusing on. In 2018, GIMPs did change their final test to the Fermat Probable Prime test since it could be implemented with a lower chance of hardware error. They improved the method again in 2020 when new work by Krzysztof Pietrzak eliminated the need to double-check results. It’s exciting to see new techniques still being created to solve such an old problem, but we’ll stick with the Lucas-Lehmer test for now.
The Lucas-Lehmer Primality Test
Computers are helpful but the Lucas-Lehmer test was originally developed and done by hand. The method is so easy to write down it fits in 10 lines of Python code. Proving it works is not as straightforward, so we’ll understand a related piece of mathematics to get a feel for how the Lucas-Lehmer test works. Let’s focus on Fermat’s Little Theorem as a way to understand how the behavior of prime numbers can help us build a test to find them. We won’t prove anything here but I hope this provides a sense of how mathematicians can turn a property of a number into a way to detect it.
Fermat’s Little Theorem is a statement that shows one way in which prime numbers tie together exponents and modular arithmetic. It gives us an exponent that can be used to make any base equal to $1 \mod p$.
Fermat’s Little Theorem:
\[ a^{p-1} \equiv 1 \mod p \]
where $a$ is an integer and $p$ is a prime.
This feels very abstract, so let’s try the example from Wikipedia where $a=2$ and $p=7$.
\[ a^{p-1} \equiv 1 \mod p \implies a^{p-1} – 1 \equiv 0 \mod p \]
Keep in mind that $\equiv 0 \mod p$ is the same thing as saying a number is a multiple of $p$.
\[ 2^{7-1} – 1= 2^6 -1 = 64 -1= 63 \text{ is a multiple of 7, so } 63 \equiv 0 \mod 7 \]
The takeaway from this is to understand that raising any base to the power $p-1$ and having the result be equivalent to $1 \mod p$ is a property of every prime number, no matter how big. This means we know something about the prime numbers that have not yet been discovered. We use this knowledge to our advantage by changing how we use the equation. If we have a number that we hope is prime but aren’t sure, we can use Fermat’s Little Theorem to see if the candidate behaveslike a prime.
Let’s consider a candidate, $n$, and see how we can use Fermat’s Little Theorem to test if the candidate behaves like a prime. If we have a candidate, $n$, and check that $a^{n-1} \equiv 1 \mod n$ for many values of $a$, then $n$ is behaving like a prime number. We have to be careful here because the converse of Fermat’s Little Theorem here is not true. In our example, $n$ is behaving like a prime number but we don’t know for sure if $n$ is a prime number. The Carmichael Numbers live in the murky gap of behaving like something without being that thing, but we won’t venture there today.
We introduced Fermat’s Little Theorem to show how true statements about prime numbers can be used to detect them. We took a statement that is true for all primes and then used that relationship to test if an integer is behaving the same way prime numbers do. Unlike our approach with Fermat’s Little Theorem, the Lucas-Lehmer test is definitive and does not have the same “walks like a duck, quacks like a duck, but isn’t a duck” problem that the Carmichael Numbers create. Lucas-Lehmer’s proof closes that gap in two steps. First, Fermat’s Little Theorem is true for all primes, while the Lucas-Lehmer test only applies to primes of a very specific form. Second, the Lucas-Lehmer test uses iteration alongside exponentiation and modular arithmetic to create a sequence of numbers that end in zero if and only if the candidate is prime. This iteration is also why the test can be slow and why GIMPS runs it last.
With the background that the Lucas-Lehmer Test uses iteration to definitively prove a number prime, we can finally write it down. The entire test fits in 10 lines of Python. This code isn’t good enough to set new records, but it’ll prove, and beat, Lucas!
def LucasLehmer(p):
'''Returns True if the number 2^(p - 1) is prime, otherwise returns False'''
s = 4
M = (2**p) - 1
for i in range(p - 2):
s = (s**2 - 2) % M
if s == 0:
return True
else:
return False
p = 127
lucas_prime = (2 ** p) - 1
if LucasLehmer(p):
print(lucas_prime, 'is prime!')
else:
print(lucas_prime, 'is not prime.')
Running the code above gives us the output below.
170141183460469231731687303715884105727 is prime!
The number that took the record from Lucas, and from humans once and for all, was $M_{521} = 2^{521} – 1$. We can test that one too.
p = 521
M_p = (2 ** p) - 1
if LucasLehmer(p):
print(M_p, 'is prime!')
else:
print(M_p, 'is not prime.')
6864797660130609714981900799081393217269435300143305409394463459185543183397656052122559640661454554977296311391480858037121987999716643812574028291115057151 is prime!
We’ve beaten Lucas! To be fair, we did not really find a number bigger than Lucas’ – we looked one up and checked it. Next, we show how to discover $M_{521}$ for ourselves and then find even bigger primes. We use the same ideas GIMPS does when deciding which numbers to test next.
How do we know which numbers to test?
The Lucas-Lehmer test we’re using is powerful and easy to write down, but like many good things, these advantages come with a trade-off. The trade-off is our test only works in very specific settings. This test can only determine if the number $M_p = 2^p – 1$ is prime. It is not guaranteed that $M_p$ is prime, but $p$ must be prime for $M_p$ to even stand a chance. In short, we need a small prime number, like 2, 3, 127, or 521, to try to make a bigger one.
The proof that $p$ must be prime is challenging but we can try a few values to see how this idea could work.
$p=2$ is prime and $2^2 – 1 = 4 – 1 = 3$ is also prime.
$p=3$ is prime and $2^3 – 1 = 8 – 1 = 7$ is also prime.
$p=4$ is not prime and $2^4 – 1 = 16 – 1 = 15 = 3 \times 5$ is also not prime.
Every time I focus on abstractions for too long I’m always a bit shocked when examples like this work out. For completeness, we’ll see that if $p$ is prime $M_p$ doesn’t have to be prime, so one last example.
$p=11$ is prime and $2^{11} – 1 = 2048 – 1 = 2047 = 23 \times 89$ is not prime.
The point is that before we find a big prime we need a lot of small ones. There are several ways to find small primes, but an easy method that is fast enough for us is the Sieve of Eratosthenes. This wonderful animation from Wikipedia shows how it works, but you don’t need to understand the details to run the code.
By SKopp at German Wikipedia, used under the CC-BY-SA 3.0 licence.
def sieveOfEratosthenes(n):
'''Returns a list of all prime numbers up to n. Fairly quick for n up to 1,000,000.'''
candidates = {i:True for i in range(2, n)}
for i in candidates:
if candidates[i]:
for j in range(2*i, n, i):
candidates[j] = False
primes = [i for i in candidates if candidates[i]]
return primes
With this function we can generate lots of small primes. Let’s get all the primes up to 1,000 and beat Lucas without cheating.
N = 1000
primes = sieveOfEratosthenes(N)
print(f'Testing the {len(primes)} primes less than {N}. The largest prime less than {N} is {max(primes)}.')
mersenne_primes = []
for p in primes:
if LucasLehmer(p):
print(f'2^{p} - 1 = {2**p -1} is prime')
mersenne_primes.append(p)
print(f'\nFound {len(mersenne_primes)} Mersenne Primes - {mersenne_primes}')
print(f'The largest is 2^{max(mersenne_primes)} - 1 and it has {len(str(2**max(mersenne_primes)-1))} digits.')
The block above will generate a list of small prime numbers, each less than 1,000, and test each Mersenne Prime $M_p$ using the Lucas-Lehmer Primality Test. The code prints out each Mersenne Prime as it finds it and then tells us how big the longest one is. Its output is below.
Testing the 168 primes less than 1000. The largest prime less than 1000 is 997.
2^3 - 1 = 7 is prime
2^5 - 1 = 31 is prime
2^7 - 1 = 127 is prime
2^13 - 1 = 8191 is prime
2^17 - 1 = 131071 is prime
2^19 - 1 = 524287 is prime
2^31 - 1 = 2147483647 is prime
2^61 - 1 = 2305843009213693951 is prime
2^89 - 1 = 618970019642690137449562111 is prime
2^107 - 1 = 162259276829213363391578010288127 is prime
2^127 - 1 = 170141183460469231731687303715884105727 is prime
2^521 - 1 = 6864797660130609714981900799081393217269435300143305409394463459185543183397656052122559640661454554977296311391480858037121987999716643812574028291115057151 is prime
2^607 - 1 = 531137992816767098689588206552468627329593117727031923199444138200403559860852242739162502265229285668889329486246501015346579337652707239409519978766587351943831270835393219031728127 is prime
Found 13 Mersenne Primes - [3, 5, 7, 13, 17, 19, 31, 61, 89, 107, 127, 521, 607]
The largest is 2^607 - 1 and it has 183 digits.
The largest prime we found is $M_{607}$, which is even bigger than the 1952 recorder holder $M_{527}$. That means your computer has caught up with more than 2,000 years of mathematicians. Your computer just beat every human who has ever tried to find a large prime number by hand and you know something the famous Lucas did not. Try turning up the $N$ from $1,000$ to $2,000$ and see if your computer can find a prime so big it can’t fit in a tweet – that is, find a prime with more than 280 digits.
Going further
My laptop was able to check all numbers $2^p -1 $ for $p < 10,000$ and found the following $p$ generate Mersenne Primes.
That matches the OEIS list and Wikipedia page for Mersenne primes with $p < 10,000$. The biggest of these, $M_{9941} = 2^{9941} – 1$, held the record in 1963 and is mind-bogglingly big.
M_9941 has 2993 digits!
34608828249085121524296039576741331672262866890023854779048928344500622080983411446436437554415370753366448674763505018641470709332373970608376690404229265789647993709760358469552319045484910050304149809818540283507159683562232941968059762281334544739720849260904855192770626054911793590389060795981163838721432994278763633095377438194844866471124967685798888172212033000821469684464956146997194126921284336206463313859537577200462442029064681326087558257488470489384243989270236884978643063093004422939603370010546595386302009073043944482202559097406700597330570799507832963130938739885080198416258635194522913042562936679859587495721031173747796418895060701941717506001937152430032363631934265798516236047451209089864707430780362298307038193445486493756647991804258775574973833903315735082891029392359352758617185019942554834671861074548772439880729606244911940066680112823824095816458261761861746604034802056466823143718255492784779380991749580255263323326536457743894150848953969902818530057870876229329803338285735419228259022169602665532210834789602051686546011466737981306056247480055071718250333737502267307344178512950738594330684340802698228963986562732597175372087295649072830289749771358330867951508710859216743218522918811670637448496498549094430541277444079407989539857469452772132166580885754360477408842913327292948696897496141614919739845432835894324473601387609643750514699215032683744527071718684091832170948369396280061184593746143589068811190253101873595319156107319196071150598488070027088705842749605203063194191166922106176157609367241948160625989032127984748081075324382632093913796444665700601391278360323002267434295194325607280661260119378719405151497555187549252134264394645963853964913309697776533329401822158003182889278072368602128982710306618115118964131893657845400296860012420391376964670183983594954112484565597312460737798777092071706710824503707457220155015899591766244957768006802482976673920392995410164224776445671222149803657927708412925555542817045572430846389988129960519227313987291200902060882060733762075892299473666405897427035811786879875694315078654420055603469625309399653955932310466430039146465805452965014040019423897552675534768248624631951431493188170905972588780111850281190559073677771187432814088678674286302108275149258477101296451833651979717375170900505673645964696355331369819296000267389583289299126738345726980325998955997501176664201042888546085699446442834195232948787488410595750197438786353119204210855804692460582533832967771946911459901921324984968810021189968284941331573164056304725480868921823442538199590383852412786840833479611419970101792978355653650755329138298654246225346827207503606740745956958127383748717825918527473164970582095181312905519242710280573023145554793628499010509296055849712377978984921839997037415897674154830708629145484724536724572622450131479992681684310464449439022250504859250834761894788889552527898400988196200014868575640233136509145628127191354858275083907891469979019426224883789463551 is prime.
To catch up with the current record holder you would need a much more powerful computer and an improved version of our test. Currently, GIMPS takes about a month to test a single number and you can contribute directly to that project by running their code on your computer. You can learn more about the Greater Internet Mersenne Prime Search on their “How GIMPS Works” page.
A conversation about mathematics and education inspired by a hundred square. Presented by Katie Steckles and Peter Rowlett, with special guest Susan Okereke.
