The Gamma-Integral (Mellin Transform) Representation of ζ(s) — A Deep Dive

$$\zeta(s)\,\Gamma(s) = \int_0^\infty \frac{x^{s-1}}{e^x - 1}\,dx, \qquad \operatorname{Re}(s) = \sigma > 1$$

1. What This Equation Is, Structurally

This identity is a Mellin transform. The Mellin transform of a function f(x) is defined as

$$\mathcal{M}[f](s) = \int_0^\infty x^{s-1} f(x)\,dx,$$

so the equation says precisely:

ζ(s)Γ(s) is the Mellin transform of the Bose–Einstein kernel 1/(eˣ − 1).

This viewpoint is not cosmetic — it is the engine of modern analytic number theory. Mellin transforms convert Dirichlet series (sums over n) into integrals, and integrals can be deformed, shifted, and analyzed by contour methods in ways that raw infinite sums cannot. Nearly every analytic continuation and functional equation in the field (zeta, L-functions, modular forms via theta series) is proven through some variant of this maneuver.

2. History and Attribution

Who came up with it? The story has three layers:

• Leonhard Euler (1729–1730) introduced the Gamma function itself as the integral $\Gamma(s)=\int_0^\infty t^{s-1}e^{-t}dt$ (in a letter to Goldbach, and in his paper De progressionibus transcendentibus). Euler also studied ζ(s) for real s and proved the Euler product (1737), but he did not write this particular integral identity for ζ.

• Niels Henrik Abel (1823) derived closely related integral formulas connecting sums like Σn⁻ˢ to integrals of x^{s−1}/(eˣ−1) type in his early work on functional equations and series (Solution de quelques problèmes à l'aide d'intégrales définies). The identity was essentially "in the air" among analysts of the 1820s–1840s (Abel, Plana, Cauchy, Malmsten).

• Bernhard Riemann (1859) is the one who made it matter. In his eight-page memoir — the only paper he ever wrote on number theory —

"Ueber die Anzahl der Primzahlen unter einer gegebenen Grösse" ("On the Number of Primes Less Than a Given Magnitude"), Monatsberichte der Berliner Akademie, November 1859,

Riemann opens with exactly this equation (his equation appears in the second displayed formula of the paper), then immediately deforms the integration path into what we now call a Hankel contour (named after Hermann Hankel, who used the same contour for Γ in 1864 — Riemann actually used it first here). This yields

$$\zeta(s) = \frac{\Gamma(1-s)}{2\pi i}\oint_{\mathcal{H}} \frac{(-z)^{s-1}}{e^z - 1}\,dz,$$

valid for all s ∈ ℂ except s = 1, giving the first proof that ζ extends meromorphically to the whole plane. Summing the residues at z = 2πin then yields the functional equation

$$\zeta(s) = 2^s \pi^{s-1} \sin\!\left(\frac{\pi s}{2}\right)\Gamma(1-s)\,\zeta(1-s).$$

It is on this foundation that Riemann states, almost in passing, the conjecture we now call the Riemann Hypothesis: all non-trivial zeros have real part ½.

• Hjalmar Mellin (1890s–1900s) later systematized the transform pair that bears his name, making Riemann's trick into a general-purpose tool.

3. Full Derivation (Expanded, With Rigor Notes)

Step 1 — Scale the Gamma integral. For σ > 0 and integer n ≥ 1, substitute t = nx in Γ(s):

$$\Gamma(s) = \int_0^\infty t^{s-1}e^{-t}\,dt = n^s\int_0^\infty x^{s-1}e^{-nx}\,dx \;\Longrightarrow\; \frac{1}{n^s} = \frac{1}{\Gamma(s)}\int_0^\infty x^{s-1}e^{-nx}\,dx.$$

This "unfolds" each term of the zeta series into an integral.

Step 2 — Sum over n and interchange Σ and ∫. For real s > 1 every term x^{s−1}e^{−nx} is non-negative, so Tonelli's theorem permits the swap unconditionally:

$$\sum_{n=1}^{\infty}\frac{1}{n^s} = \frac{1}{\Gamma(s)}\int_0^\infty x^{s-1}\sum_{n=1}^\infty e^{-nx}\,dx.$$

For complex s with σ > 1, one checks absolute convergence (|x^{s−1}| = x^{σ−1}) and invokes Fubini; alternatively, both sides are analytic in σ > 1 and agree on (1, ∞), so they agree everywhere there by the identity theorem.

Step 3 — Geometric series. For x > 0, e^{−x} < 1, so

$$\sum_{n=1}^\infty e^{-nx} = \frac{e^{-x}}{1-e^{-x}} = \frac{1}{e^x - 1}. \qquad\blacksquare$$

Convergence check at the endpoints (why σ > 1 is needed):

• Near x = 0: 1/(eˣ − 1) ~ 1/x, so the integrand ~ x^{σ−2}, integrable iff σ > 1. This is the sharp constraint.
• Near x = ∞: the integrand decays like x^{σ−1}e^{−x} — always fine.

The pole of ζ at s = 1 is thus visible in the integral: as σ → 1⁺, the divergence at x = 0 produces exactly the simple pole with residue 1.

4. Why It Matters for RH Research

1. Gateway to analytic continuation. RH is a statement about zeros in 0 < σ < 1, where the Dirichlet series diverges. Without continuation there is no hypothesis to state. This integral is the historical (and still the cleanest) route.

2. Template for the "modern" continuation. Replacing 1/(eˣ−1) with the theta function θ(x) = Σe^{−πn²x} in the same Mellin machinery gives Riemann's second proof and the completed zeta function ξ(s) — the symmetric object ξ(s) = ξ(1−s) at the heart of RH.

3. Physics connection. The kernel 1/(eˣ−1) is the Bose–Einstein distribution. Setting s = 4 gives ∫₀^∞ x³/(eˣ−1)dx = Γ(4)ζ(4) = π⁴/15, which is exactly the integral in the Stefan–Boltzmann law for blackbody radiation. This zeta–quantum-statistics link feeds into the Hilbert–Pólya program (seeking a self-adjoint operator whose spectrum gives the zeta zeros) — see Berry & Keating's work.

5. Key References (with links)

• Riemann's original 1859 paper — English translation by D. R. Wilkins, hosted at Trinity College Dublin: https://www.maths.tcd.ie/pub/HistMath/People/Riemann/Zeta/EZeta.pdf (German original also there: https://www.maths.tcd.ie/pub/HistMath/People/Riemann/Zeta/)
• Clay Mathematics Institute — RH official problem description by E. Bombieri: https://www.claymath.org/millennium/riemann-hypothesis/
• H. M. Edwards, Riemann's Zeta Function (Academic Press 1974; Dover reprint) — Chapter 1 walks through the 1859 paper line by line, including this integral and the Hankel contour. The single best companion for your paper.
• E. C. Titchmarsh (rev. D. R. Heath-Brown), The Theory of the Riemann Zeta-Function, 2nd ed., Oxford, 1986 — §2.4 gives this representation and its contour version rigorously.
• M. V. Berry & J. P. Keating, "The Riemann zeros and eigenvalue asymptotics," SIAM Review 41 (1999) 236–266: https://epubs.siam.org/doi/10.1137/S0036144598347497
• Abel's related integral work: Œuvres complètes de Niels Henrik Abel, available at https://archive.org/details/oeuvrescompltesd01abel

(Note: URLs are correct as of my training data; if a link has moved, searching the title will locate it.)

6. Computational Examples in Python

Example 1 — Verify the identity numerically for real and complex s

import numpy as np
from scipy.integrate import quad
from mpmath import mp, mpf, mpc, zeta, gamma, quad as mpquad, exp, pi

mp.dps = 30  # 30 digits of precision

def zeta_via_integral(s):
    """Compute zeta(s) via the Gamma-integral representation, Re(s) > 1."""
    integrand = lambda x: x**(s - 1) / (exp(x) - 1)
    integral = mpquad(integrand, [0, mp.inf])
    return integral / gamma(s)

# --- Real s: check against known closed forms ---
tests = {
    2: pi**2 / 6,        # Basel problem (Euler, 1734)
    4: pi**4 / 90,
    6: pi**6 / 945,
    3: zeta(3),          # Apery's constant, no closed form
}
print("Real s:")
for s, exact in tests.items():
    approx = zeta_via_integral(mpf(s))
    print(f"  zeta({s}): integral = {approx},  exact = {exact},  "
          f"err = {abs(approx - exact):.2e}")

# --- Complex s (still needs sigma > 1) ---
print("\nComplex s:")
for s in [mpc(2, 5), mpc(1.5, 14.134725), mpc(3, -7)]:
    approx = zeta_via_integral(s)
    exact  = zeta(s)
    print(f"  zeta({s}): err = {abs(approx - exact):.2e}")

Typical output shows agreement to ~28–30 digits — the identity holds on the nose.

Example 2 — Watch the pole at s = 1 emerge from the x → 0 divergence

from mpmath import mp, mpf, zeta, gamma, quad as mpquad, exp
mp.dps = 25

def zeta_via_integral(s):
    return mpquad(lambda x: x**(s-1)/(exp(x)-1), [0, mp.inf]) / gamma(s)

print("s -> 1+ : (s-1)*zeta(s) should approach 1 (residue of the pole)")
for eps in [mpf('0.1'), mpf('0.01'), mpf('0.001'), mpf('0.0001')]:
    s = 1 + eps
    z = zeta_via_integral(s)
    print(f"  s = 1 + {eps}:  zeta(s) = {float(z):14.4f},  (s-1)*zeta(s) = {float(eps*z):.6f}")

You will see (s−1)ζ(s) → 1, confirming the simple pole with residue 1 directly from the integral.

Example 3 — The physics payoff: Stefan–Boltzmann integral

import numpy as np
from scipy.integrate import quad
from math import pi

# Planck/Bose-Einstein integral: ∫ x^3/(e^x - 1) dx = Γ(4)ζ(4) = 6 * π^4/90 = π^4/15
val, err = quad(lambda x: x**3 / np.expm1(x), 0, np.inf)
print(f"Integral  = {val:.12f}")
print(f"pi^4 / 15 = {pi**4/15:.12f}")   # matches — this constant sits inside σ_SB

Example 4 — Riemann's Hankel-contour continuation, numerically

The contour version $\zeta(s) = \frac{\Gamma(1-s)}{2\pi i}\oint_{\mathcal H}\frac{(-z)^{s-1}}{e^z-1}dz$ works for all s ≠ 1. We approximate the Hankel contour (in from +∞ below the axis, circle the origin, back out above):

from mpmath import mp, mpc, mpf, exp, gamma, zeta, pi, mpmathify, quad as mpquad
mp.dps = 25

def zeta_hankel(s, R=1.0, X=60):
    """Riemann's contour formula. (-z)^{s-1} uses the branch cut along positive real axis."""
    j = mpc(0, 1)
    f = lambda z: (-z)**(s - 1) / (exp(z) - 1)

    # Leg 1: from X - i*eps down... parametrize just below axis: z = t*exp(-i*pi_almost)
    # Cleaner: use log form. Below axis: (-z)^{s-1} = exp((s-1)*(log t - i*pi)) for z = t (t>0)
    below = mpquad(lambda t: exp((s-1)*(mp.log(t) - j*pi)) / (exp(t) - 1), [X, R])
    above = mpquad(lambda t: exp((s-1)*(mp.log(t) + j*pi)) / (exp(t) - 1), [R, X])
    # Circle of radius R around origin, theta from -pi to pi
    circ  = mpquad(lambda th: f(R*exp(j*th)) * j*R*exp(j*th), [-pi, pi])

    contour = below + circ + above
    return gamma(1 - s) * contour / (2*pi*j)

# Test at points where the ORIGINAL integral diverges (sigma <= 1):
for s in [mpc('0.5','14.134725141734693'),  # near the first nontrivial zero!
          mpf('-1'),                        # zeta(-1) = -1/12
          mpf('0'),                         # zeta(0)  = -1/2
          mpc('0.5','0')]:
    print(f"s = {s}:  Hankel = {zeta_hankel(s)}")
    print(f"           mpmath = {zeta(s)}\n")

Note the highlights: it recovers ζ(−1) = −1/12, ζ(0) = −1/2, and evaluates ζ on the critical line σ = ½ — the home of the Riemann Hypothesis — where it should return ≈ 0 at t ≈ 14.1347251417…, the first nontrivial zero (Gram, 1903).

7. Summary Table for Your Paper