The standard shortcut overstates photon emission by a factor of 4.16. Six validation gates, two analytic limits the code was never tuned to, and one number I had to retract.
My last few articles were the same problem in different clothes: a metric that looks healthy while the thing underneath is broken. Word Error Rate passing mangled speech. A loss curve falling 62% on pure noise. An objective with a hole in it, and two training curves you cannot tell apart.
This one is that problem in its hardest form — when there is no experiment to check against at all.
If you write physics code, nobody hands you a validation set. The black hole is not going to email you. So I picked a calculation with a known answer and tried to reach it without looking: the photon emission spectrum of a Schwarzschild black hole, first computed by Don Page in 1976. The question was not really whether I could get his number. It was what it actually takes to believe your own output.
The shortcut, and where it breaks
Hawking radiation is not a blackbody spectrum. The curvature around the hole acts as a potential barrier that reflects most of what is emitted; the fraction escaping is the greybody factor, and getting it right means solving a wave equation mode by mode.
Almost everyone starts with the shortcut — the geometric-optics limit, which treats the hole as swallowing anything aimed inside the photon capture cross-section. It is clean, analytic, and valid only when the emitted wavelength is small compared with the hole.
That condition has an exact form. With u = E/(k_B·T_H):
E · r_s / (hbar·c) = 2GME/(hbar·c^3) = u / 4pi
The thermal spectrum peaks at u = 1.594. Substitute, and you get 0.127.
The typical emitted photon has a wavelength about eight times the Schwarzschild radius. The approximation is not marginally strained here — it fails across the entire range where the hole actually radiates, and the exponential cutoff kills the spectrum long before you reach the regime where it would become valid. Every textbook says the limit requires short wavelengths. Nobody makes you check.
Two details that decided the answer
A massless spin-1 field on a Schwarzschild background separates into a one-dimensional scattering problem against a barrier peaking at the photon sphere:
(d^2/dr*^2 + w^2 - V_l) Psi = 0, r* = r + 2 ln(r/2 - 1)
V_(S,l) = f(r) · [ l(l+1)/r^2 + (1 - S^2)·2M/r^3 ], f = 1 - 2M/r
For S = 1 the spin term vanishes identically. The greybody factor is the transmission probability through that barrier, summed over modes and over the photon's two polarization states — and a detail worth flagging: for Maxwell on Schwarzschild both field parities obey this same radial equation, and those parities are the two polarizations. That is what pins the degeneracy factor at exactly 2, with no double counting.
Two implementation choices determined whether the answer came out right at all.
Don't shoot from the horizon — expand around it. The horizon is a singular point, so integrating outward from beside it is numerically hostile. Substituting Psi = exp(-i·w·r*)·u(r) makes the equation merely regular-singular, and in x = r - 2 it becomes polynomial with a clean recursion. The nearest other singularity sits at x = -2, so the series converges comfortably and you step over the stiff region instead of fighting it.
The long-range tail, which nearly ruined it. The potential falls off as l(l+1)/r^2 — long-range. Matching the numerical solution to plain plane waves at large radius leaves an error of order L/(w·r), which at low frequency is ruinous.
This was the dominant error in my first working version, and how I found it is the point: not by looking at the answer. The spectrum looked entirely reasonable. It was caught by a flux-conservation monitor that had no business agreeing unless everything else was right. Expanding the asymptotic solution to second order fixed it and dropped the truncation to O((L/2wr)^3), which lets you match at modest radius instead of integrating to absurd distances.
The correction, across the spectrum
Geometric optics (dashed) overstates low-energy photon emission by up to six orders of magnitude. The curves agree only above ~7 TeV, where the wavelength finally becomes small compared with the hole.
The ripple in the lower panel, as the ratio approaches one, is not numerical noise — it is the absorption cross-section ringing around its geometric value, caused by surface waves trapped near the photon sphere. I did not put that in. It arrived on its own, which is the kind of thing that makes you trust a solver.
Six gates, two of which the code could not have anticipated
Agreement with one published number proves very little. Anyone can tune until a single value matches. The checks that mean something are the ones with an independently known answer the implementation had no way to see coming.
| Gate | What it tests | Residual |
|---|---|---|
| A | Convergence under start radius, step size, match radius, series length | 3.1e-08 |
| B | Low-frequency scaling law; measured slopes 4.005, 6.005, 8.005 | 5.5e-05 |
| C | High-frequency limit returns to the analytic geometric-optics form | 2.2% |
| D | Flux conservation, never imposed anywhere in the solver | 2.5e-09 |
| E | Feeding the old approximation through the same pipeline reproduces the earlier frozen result | 4.1e-07 |
| F | Hawking temperature matches the frozen value | 1.9e-07 |
Gate B is the one I would show a skeptic. At low frequency the greybody factor scales as a specific power law whose leading coefficient has a known closed form. My measurement extrapolates to 0.44441988 against the analytic 4/9 = 0.44444444 — agreement to 5.5e-05. Nothing in the code knows about 4/9. It cannot be fitted, because it is not an input.
The fit also has to be done at genuinely small frequency to see it. At w ~ 0.005 the measured slope still reads 4.04, and only settles to 4.005 an order of magnitude further down. Stopping early would have produced a number that looked wrong and wasn't.
Gate D nearly failed, and I didn't loosen it
Flux conservation — what goes in comes out — is never imposed anywhere in the solver, so it is a free error monitor. It also failed my first threshold, and the instinct is to widen the threshold until it passes.
Instead I measured why. The residual falls as h^5 with step size, while the quantity I actually care about moves only in its ninth significant figure. And when the transmission probability is ~1e-16, the identity requires sixteen digits of cancellation — more than double precision has.
So the test is not failing. It is inapplicable in that regime, for a reason you can write down. The honest fix was to apply it only where it is meaningful and document exactly why — not to tune the number until the dashboard went green. That distinction is the whole job.
The result, and the retraction
Photon luminosity written as P = alpha · hbar·c^6/(G^2·M^2):
| Source | alpha |
|---|---|
| Geometric optics | 1.399e-04 |
| This work (Regge–Wheeler, s = 1) | 3.364e-05 |
| Page 1976 (from his 17% photon share) | 3.400e-05 |
Independent implementation, 1.06% from the published value. The shortcut overstates the luminosity by 4.16×.
Now the part I would rather leave out.
An earlier draft compared against a more precise-looking figure for Page's total power — four significant figures instead of one. It came from a secondary summary, not from Page. When I finally pulled the published abstract, that number was not in it. I had propagated a value I never verified, and the tighter agreement it implied was an illusion.
It is retracted, and the comparison now uses Page's own published figures. Which forces a second correction: 1.06% is not a 1% validation. Page's abstract quotes the total to one significant figure. The honest claim is consistency at the precision the source supports — no tighter. Inverted, this calculation implies a photon share of 16.82% against his stated 17%.
What it still rests on
- Page's abstract is verified; his paper is not. The comparison values were read from the published abstract. The body is paywalled and I have not read it.
- Photons only, Schwarzschild only, primary spectrum only, fixed mass. No neutrinos, no gravitons, no spin, no charge, no secondary emission, no evaporation.
- This is reproduction, not discovery. Page did it fifty years ago and the physics has been standard ever since. The point was to own the number instead of citing it.
The transferable part
The physics here is fifty years old. The method is not, and it is the same one behind the linters:
- Test against limits you cannot fit. An analytic value the code has never seen is worth more than a hundred plausible-looking outputs.
- Keep a monitor you never optimise. Flux conservation found a real bug that the answer itself concealed.
- When a check fails, diagnose before you retune. Measure the scaling. Find out whether it is failing or inapplicable. Those are different, and only one is fixed by changing a threshold.
- Read the primary source. A number you have not checked is not evidence, no matter how many digits it carries.
One file, one command, no external physics code and no scipy — quadrature, root-finding and the series recursion are all hand-rolled. It runs in about three minutes and prints every number above.
Related: loss curves lie and one integer that deleted my stop token — the same lesson where a validation set does exist. Primary reference: D. N. Page, Phys. Rev. D 13, 198 (1976).