Primes and Polynomials March 20, 2021

To put it mildly, the prime numbers are the subject of some mathematical interest, and of great mathematical mystery. Intriguingly, it’s possible to come up with a polynomial $P : \mathbb{Z}^n \to \mathbb{Z}$ in $n$ variables such that the range of positive values taken on by $P$ (on nonnegative inputs) enumerates the primes. To put it more precisely, $P(\mathbb{N}^n) \cap \mathbb{N}$ is the set of prime numbers. This result is due to Yuri Matiyasevich.

The cleanest such construction is probably the Jones-Sato-Wada-Wiens polynomial, which has $n=26$ variables and is of degree $25$. (This is convenient, because the English alphabet happens to have $26$ letters.) I won’t try to strain my poor website’s mathematical typesetting system by writing out the polynomial here, but you can see it in their original paper¹. The main idea is to encode the relationship $x = y!$ in a system of Diophantine equations, and then characterize the prime numbers using Wilson’s theorem. You can then combine the system of equations into one polynomial.

One day, this polynomial came up during lecture, and the professor jokingly said that we could come up with an assignment to the $26$ variables witnessing a small prime like $5$ for “extra credit.” I thought it’d be a fun challenge, so I decided to try my hand at it. I ended up working out the numbers (or at least a closed-form expression for the numbers) for the prime $2$, but anything higher would be completely computationally intractable, at least the way I did it. A more mathematically detailed description of the solution than what I could fit in a blog post is available as a PDF—this post is more of a quick play-by-play of the process that I went through.

Wisdom of the Ancients

(Relevant XKCD: #979)

After some online searching, the only thing that I could find was a Google Groups discussion from nearly twenty years ago, with one poster giving a purported assignment for the prime $2$:

a = 20
b = 0
c = 0
d = 0
e = 3
f = 1
g = 0
h = 1
i = 0
j = 0
k = 0
l = 0
m = 1
n = 0
o = 244
p = 1
q = 1
r = 0
s = 0
t = 0
u = 1
v = 0
w = 0
x = 1
y = 0
z = 1

This very nearly works², except that $$16(k+1)^3(k+2)(n+1)^2 + 1 - f^2 = 32 \neq 0.$$ I figured that it couldn’t be too hard to fix, so I messed around with some numbers for a bit in Julia. But I found that the behavior of the polynomial was rather unpredictable, and with a search space consisting of $26$ parameters that I had no bound for, it was clear that I’d need a better approach.

Computers Are Our Friends

After skimming the original paper, it became apparent that there was a lot of structure to the way the authors had designed the polynomial. In particular, it all seemed to revolve around the Pell equation $$x^2 = (a^2-1)y^2 + 1.$$ By reading the authors’ original proof, I was able to reverse engineer the values of a few of the smaller parameters, like $k$ (obviously) and $n$. Given this foothold, I then needed to solve a chain of these Pell equations; I hashed out a quick Python program to search for a solution, confident that my work would be over soon.

Unfortunately, more than an hour later, my program still had not found a solution. Perhaps this would require some more refined intellectual firepower; I had ruled out pretty much any values small enough to find with brute force. I tried reading a textbook on solving Diophantine equations, but it just confused me, so I did the next best thing: I plugged the system into Mathematica and hoped it would work. To my surprise, Mathematica immediately identified a viable solution. To my dismay, it was far too large; while I could solve for some of the parameters, the remaining ones would be computationally out of reach.

So I did the first thing that all good mathematicians do: if you can’t solve a problem, try to solve an easier version. Rather than try the ambitiously large prime number $5$, I decided to try to find a solution for the much more modest prime $2$. This initially worked well; most of the numbers were pretty small, like $f = 17$ and $e = 32$. Unfortunately, I then ran into a bit of a problem:

FindInstance[e^3 (e + 2) (a + 1)^2 + 1 == o^2, {a, o}, Integers]

This returned:

a = 7901690358098896161685556879749949186326380713409290912
o = 8340353015645794683299462704812268882126086134656108363777

Clearly, things were getting a little out of hand. Perhaps Mathematica wasn’t returning the optimal solution? Unfortunately, the authors proved that this condition would force $$a \ge e - 1 + e^{e-2} = 1427247692705959881058285969449495136382746655,$$ so it seemed unlikely that a significantly smaller solution existed. I had gotten pretty far, and with some more careful reading of this paper and a related one³, I was able to get all of the parameters besides $c$, $d$, $r$, and $u$.

Numbers Are Hard

At this point, I was pretty sure that I was on the right track. By reverse engineering the original proofs and applying some computational power, I had arrived at a partial assignment satisfying almost every constraint I needed—except for two. Unfortunately, these remaining two constraints required a little bit of (fairly elementary) number theory to resolve, and I am pretty bad at that. So I did the second thing that all good mathematicians do: if you’re stumped, put the problem down for a month and do something else.

So I went and did some other stuff for a month. Later, after my first round of midterms, I decided to come back to this problem for a bit. I wasn’t expecting much, but to my surprise, everything seemed to make much more sense a month later. I still wasn’t able to figure out numeric values for the remaining parameters $c$, $d$, $r$, and $u$, but I did get closed form expressions that should (in theory) yield the correct answer. So yay, I guess?

If you’re interested, I’ll again plug the more mathematically detailed writeup, and I’ll mention a short program that verifies most of the constraints. There isn’t really a moral to this story; it was just kind of fun.

As a side note, I had originally titled this post “Prime Polynomials,” but I realized that this choice of title was unfortunate, because the term “prime polynomial” is also used to refer to irreducible polynomials over a unique factorization domain. This has little to do with the sense in which I was thinking of the term when I originally titled this post, and so I’ve changed the title to be less ambiguous.