With all the recent math AI breakthroughs in my feed, I got reintroduced to some problems i remember reading about at university. One of which is the famous Riemann Hypothesis. I never studied it deeply, but remember hearing about the famous $\zeta(-1) = 1 + 2 + 3 + 4 … = -\frac{1}{12}$ and that the Riemann Zeta functions are somehow related to the primes.
Fast forward to today and I saw a simulation of the Riemann zeta zeros adding up to the prime counting function in the video What Does the Riemann Hypothesis ACTUALLY Mean? by Physics Explained, and got very excited. Like any physicist schooled in signal processing would, when they see a sum of waves approximating a stepped function.
The step function in question is the prime counting function $\pi(x)$, the number of primes up to $x$. It jumps by one at every prime. Here it is (orange), next to what you get from the formula below when you add the Riemann zeta zeros one at a time (blue).
1. The Riemann Formula
If you look up the standard Riemann Formula you find the below equation. Where $\pi_0(s)$ is the estimated prime counting function and $\mu(m)$ is the Möbius function. The final atypical notation is $\operatorname{Li}(x)$ which corresponds to the integral of the natural logarithm, that Gauss showed has related to the prime counting function. For a full explanation of how Riemann arrived at this result I highly recommend the video I linked above.
$$\pi_0(x)=\sum_{m\ge1}\frac{\mu(m)}{m}J\left(x^{1/m}\right),\qquad J(x)=\operatorname{Li}(x)-\sum_{\rho}\operatorname{Li}(x^{\rho})-\ln 2+\int_x^\infty\frac{dt}{t(t^2-1)\ln t}$$
The $\rho$ is relating to a point in the complex plane where $\zeta(\rho) = 0$, and the Riemann hypthosesis states that $\rho = \frac{1}{2} + i\gamma$, so all such solution are on one line. The millenium problem is proving that there are no other no-trivial zeros for $\zeta(s)$. What made we curious is that there’s no obvious place where the wave pattern is coming from, so I set out to find where it appears. But seeing something raised to a complex number may give us a clue.
2. Deriving the waves
We first note that $\operatorname{Li}(u) = \operatorname{Ei}(\ln {u})$, and swapping the variable $u = x^\rho$, we can tell that
$$\operatorname{Li}(x^\rho) = \operatorname{Ei}(\rho \ln x)$$
As seen in DLMF §6.12, you can expand $\operatorname{Ei}(x)$ when $x \rightarrow \infty$ as
$$\operatorname{Ei}(y)\approx\frac{e^{y}}{y}\left(1+\frac{1}{y}+\frac{2!}{y^{2}}+\frac{3!}{y^{3}}+\cdots\right)$$
We note the size of $y$ is $|y| = |\rho| \ln x$, and the table shows it for a few values of $x$
| Zero | $\lvert\rho\rvert$ | $\lvert y\rvert$ at $x=2$ | $\lvert y\rvert$ at $x=10$ | $\lvert y\rvert$ at $x=100$ |
|---|---|---|---|---|
| $\rho_1 \approx 1/2 + 14.134725i$ | 14.144 | 9.8 | 32.6 | 65.1 |
| $\rho_2 \approx 1/2 + 21.022039i$ | 21.028 | 14.6 | 48.4 | 96.8 |
| $\rho_3 \approx 1/2 + 25.010857i$ | 25.016 | 17.3 | 57.6 | 115.2 |
The size of $y$ grows quite fast, so the first term of the expansion is a good estimate, meaning
$$\operatorname{Ei}(y) = \operatorname{Ei}(\rho \ln x) = \operatorname{Li}(x^{\rho}) \approx\dfrac{e^{y}}{y}=\dfrac{e^{\rho\ln x}}{\rho\ln x} = \frac{x^{\rho}}{\rho\ln x}.$$
Now $\rho=\tfrac12+i\gamma$ gives $x^{\rho}=x^{1/2}x^{i\gamma}=\sqrt{x} e^{i\gamma\ln x}$, so
$$\operatorname{Li}(x^{\rho})\approx\frac{\sqrt{x} e^{i\gamma\ln x}}{\rho\ln x}.$$
And we have found our wave. $e^{i\gamma\ln x}$ oscillates with frequency $\gamma$ in $\ln x$, and the factor $\sqrt{x}/\ln x$ sets its size. To further drive this home, we note that each $\rho$ comes with a complex conjugate $\bar{\rho}$. The two terms are then complex conjugates of each other, so their sum is twice the real part:
$$\operatorname{Li}(x^{\rho})+\operatorname{Li}(x^{\bar\rho})\approx2\operatorname{Re}\left[\frac{\sqrt{x}e^{i\gamma\ln x}}{\rho\ln x}\right]$$
Writing $\rho=|\rho|e^{i\theta}$ and using $\operatorname{Re}e^{i\phi}=\cos\phi$, this is a real wave:
$$\operatorname{Li}(x^{\rho})+\operatorname{Li}(x^{\bar\rho})\approx\frac{2\sqrt{x}}{|\rho|\ln x}\cos\big(\gamma\ln x-\theta\big).$$
3. Back to the formula
Expanding the original series we get
$$\sum_{\rho} \operatorname{Li}(x^{\rho}) = (\operatorname{Li}(x^{\rho_1}) + \operatorname{Li}(x^{\bar{\rho_1}})) + (\operatorname{Li}(x^{\rho_2}) + \operatorname{Li}(x^{\bar{\rho_2}})) + … $$
and by substitution we see that
$$\sum_{\rho} \operatorname{Li}(x^{\rho}) \approx \sum_{k\ge1}\frac{2\sqrt{x}}{|\rho_k|\ln x}\cos\big(\gamma_k\ln x-\theta_k\big),$$
where $\rho_k=\tfrac12+i\gamma_k$ and $\theta_k=\arg\rho_k$.
Putting this back into $J(x)$, where the sum over zeros carries a minus sign, we get
$$J(x)\approx\operatorname{Li}(x)-\sum_{k\ge1}\frac{2\sqrt{x}}{|\rho_k|\ln x}\cos\big(\gamma_k\ln x-\theta_k\big)-\ln 2+\int_x^\infty\frac{dt}{t(t^2-1)\ln t}.$$
So, the wave pattern in the approximation is coming from a sum of waves with an increasing frequency and a changing amplitude.
Plotting the first few
Playing the wave with $x$ as time, both the pitch and the amplitude change, and we can hear the primes! For this one I used the 100th zero. It starts at 4000 Hz and falls to 20 Hz:
I hope this post offers some intuition into this quite exotic branch of mathematics. Thanks for reading!