University Physics V · The Hydrogen Atom · 11.3
Radial Hamiltonian & Effective Potential
Every central potential ends in the same place: one operator on a half-line, one unknown u = rR, and a boundary condition at the origin you have to argue for rather than assume. This is where you learn to name the domain before you solve, because the differential equation on its own will not stop you inventing states.
Build the model
Connect the measurement to the mechanism.
The model is one operator: hₗ = −(ℏ²/2μ)d²/dr² + V(r) + ℏ²l(l+1)/(2μr²), acting on u(r) = rR(r) in L²(0,∞). It comes from separation of variables, which is really the statement that (H, L², Lz) commute, so a shared eigenbasis exists and the angular factor Yₗₘ hands the radial problem one number, ℏ²l(l+1). The substitution u = rR then absorbs the r² of the volume element and cancels the first-derivative term exactly, leaving something that reads like a one-dimensional Schrödinger equation with the centrifugal term stacked on top of V.
What it costs is the endpoint. The half-line is not the line: −d²/dr² on (0,∞) is symmetric on many domains and self-adjoint on almost none of them, and the eigenvalue problem is not defined until one is chosen. Frobenius analysis at r = 0 returns two exponents, l+1 and −l, and only for l ≥ 1 does normalisation kill the second; for l = 0 the singular solution R ~ 1/r is square-integrable and has to be excluded by a domain condition, u(0) = 0, whose warrant is that ∇²(1/r) = −4πδ³(r) would otherwise plant a point source at the nucleus.
Impose it and the Coulomb spectrum follows; leave it out and the same differential equation cheerfully supplies states hydrogen does not have.
- Simple definition
- The radial Hamiltonian is the operator −(ℏ²/2μ)d²/dr² + V(r) + ℏ²l(l+1)/(2μr²) acting on u = rR on the half-line r > 0, made self-adjoint by the domain condition u(0) = 0.
- Example
- For hydrogen with l = 1 it adds a +27.21 eV centrifugal term to the −27.21 eV Coulomb term at r = a₀, so Veff(a₀) = 0 exactly; the well bottoms out at −6.803 eV at 2a₀, and the 2p level at −3.401 eV floats halfway up it.
Cancels the first-derivative term exactly and flattens the measure, so the problem reads like a 1-D one on r > 0.
μ is the reduced mass in kg, r in m, u in m⁻¹⁄²; the r² of the volume element is now inside u
Every l is a different one-dimensional problem, and the same centrifugal term appears for any central V(r).
l is the L² eigenvalue label, a number here and not an operator; Veff in J with r in m
For l = 1 the floor is −6.803 eV at 2a₀, so the 2p level at −3.401 eV sits halfway between floor and threshold.
a₀ = 0.052918 nm; minimum at r = l(l+1)a₀/Z, depth −13.6057 Z²/[l(l+1)] eV
The next order fixes u = r(l+1)[1 − Zr/((l+1)a₀) + …], the same coefficient for every state of that l.
from the leading balance u″ = l(l+1)u/r²; the Coulomb term is subleading at r = 0
Dirichlet u(0) = 0 is one self-adjoint extension; a point Coulomb nucleus is what picks it out of the Robin family.
must vanish for every u, v in the domain; normalisability handles the r → ∞ end
Tridiagonal and O(h²): at h = 0.05a₀ the 1s comes out −13.5972 eV, and Richardson lifts it to −13.6057 eV.
rᵢ = ih for i = 1…N in a box R = (N+1)h; the row never built at i = 0 is u(0) = 0
u = rR flattens the measure
Write ψ = R(r)Yₗₘ(θ, φ). The radial part of the Laplacian is R″ + (2/r)R′, and the norm drags the Jacobian along: ∫|R|²r²dr. One substitution disposes of both. Put R = u/r, so R′ = u′/r − u/r² and R″ = u″/r − 2u′/r² + 2u/r³; adding (2/r)R′ = 2u′/r² − 2u/r³ cancels those two correction terms in pairs and leaves R″ + (2/r)R′ = u″/r exactly. The norm follows suit, ∫|R|²r²dr = ∫|u|²dr, so u lives in L²(0,∞) with the flat measure and obeys −(ℏ²/2μ)u″ + Veff u = Eu. The word to keep is lives. The line has been cut in half: r ∈ (0,∞) has an endpoint, the substitution has made that endpoint singular, and everything difficult about this topic is concentrated there.
The centrifugal term is an eigenvalue, not a force
L² never acts as a differential operator on the radial factor: on Yₗₘ it returns the number ℏ²l(l+1), and that number is what enters the radial equation as ℏ²l(l+1)/(2μr²). Nothing is rotating; the barrier is bookkeeping for angular momentum the state already carries. It is also universal — the same term appears for a square well, a harmonic trap or a Yukawa potential, so it is not a Coulomb effect. Written in eV with r in Bohr radii, hydrogen's radial problem reads Veff = (−27.2114/r + 13.6057 l(l+1)/r²) eV. Put l = 1 and r = a₀: the centrifugal coefficient 2 × 13.6057 = 27.2114 exactly matches the Coulomb one, so Veff(a₀) = 0 identically, not approximately. The floor sits at r = l(l+1)a₀ = 2a₀ with depth −13.6057/[l(l+1)] = −6.803 eV, and the 2p level at −3.401 eV floats halfway between that floor and threshold.
Turning points, and where l ≤ n − 1 comes from
Set Eₙ = Veff(r) and solve. In Bohr units the quadratic Eₙ r² + 27.2114r − 13.6057 l(l+1) = 0 has roots r± = n²[1 ± √(1 − l(l+1)/n²)] a₀. For 2p that is 1.17 a₀ and 6.83 a₀, a band 5.66 a₀ wide holding ⟨r⟩ = 5a₀; for 3d it runs from 3.80 a₀ to 14.20 a₀. The square root is real only if l(l+1) ≤ n², and for integers that is exactly l ≤ n − 1. Read it off the well: raise l at fixed n and the floor −13.6057/[l(l+1)] eV rises until the level has nowhere to sit. A 2d state would have to live at −3.401 eV in a well whose floor is −2.268 eV. Truncation of the Laguerre series proves the rule; the effective potential is what makes it obvious.
Frobenius at the origin: two exponents
Near r = 0 the centrifugal term beats everything else in the equation: it is O(r⁻²) while the Coulomb term is O(r⁻¹) and E is O(1). The leading balance is therefore u″ = l(l+1)u/r². Substituting u = rs gives s(s−1)r(s−2) = l(l+1)r(s−2), so the indicial equation is s(s−1) = l(l+1), with roots s = l+1 and s = −l. The regular root gives R = u/r ~ rl, the familiar behaviour of a p or d orbital near the nucleus; the singular root gives R ~ r(−l−1). Carrying the series one term further fixes the first correction, u = r(l+1)[1 − Zr/((l+1)a₀) + …], and because it comes from the recursion — and because E first enters two orders later — it is the same for every state of that l: for l = 2 it is 1 − r/3a₀, whether the state is 3d, 4d or 5d.
Two different reasons to drop the singular root
For l ≥ 1 the argument is normalisation: u ~ r(−l) makes ∫|u|²dr ~ ∫r(−2l)dr diverge at the origin, and nothing more need be said — in Weyl's language the operator is limit point at r = 0 and needs no boundary condition at all. The dividing line is sharp: the limit-point test on the half-line asks whether l(l+1) ≥ 3/4, which l ≥ 1 clears and l = 0 fails. For l = 0 the singular root is u → constant, R ~ 1/r, and ∫|R|²r²dr = ∫dr converges. Normalisation cannot reject it. What rejects it is ∇²(1/r) = −4πδ³(r), which means R ~ 1/r satisfies the eigenvalue equation only with an unadvertised point source at the nucleus. Formally the operator is limit circle at r = 0 and admits a one-parameter family of self-adjoint extensions, the Robin conditions cos α u(0) + sin α u′(0) = 0, each of which kills the boundary term of formula five. Dirichlet is the member a point Coulomb nucleus selects: give the nucleus a radius RN, solve inside and out, then let RN → 0, and u(0) = 0 is what survives.
The domain, on a grid
Numerically the domain condition is free if the grid is built correctly. Take rᵢ = ih for i = 1…N in a box of radius R = (N+1)h and never create a node at i = 0: the missing row is Dirichlet. In atomic units the matrix is tridiagonal — diagonal 1/h² − Z/rᵢ + l(l+1)/(2rᵢ²), off-diagonal −1/(2h²) — and one symmetric tridiagonal eigensolver call returns the spectrum. Two errors then compete. The three-point second difference is O(h²): for l = 0 in a 40a₀ box the ground state comes out −0.498756 Ha at h = 0.1a₀ and −0.499688 Ha at h = 0.05a₀, and Richardson extrapolation [4E(h/2) − E(h)]/3 = −0.499998 Ha lands on −13.6057 eV. The box is the other error and it bites the top of the spectrum: shrink the box to 10a₀ and 1s is untouched at −0.499687 Ha while 2s reads −0.1128 Ha instead of −0.125 Ha and 3s is pushed above threshold as a pseudo-state.
Change one variable at a time
Make the relationship visible.
With n = 2, drag l from 0 to 1: the well floor lifts to −6.80 eV at 2 a₀ and the inner turning point leaves the origin for 1.17 a₀. Push l to 2 and the dots vanish, the level having dropped below the floor. Then set s = l + 1 to drive the indicial residual to zero.
BARRIER AT r = a₀27.21 eV
LEVEL Eₙ-3.401 eV
ALLOWED BAND r₊ − r₋5.66 a₀
INDICIAL RESIDUAL0.00
Live interpretationBARRIER AT r = a₀: 27.21 eV. LEVEL Eₙ: −3.401 eV. ALLOWED BAND r₊ − r₋: 5.66 a₀. INDICIAL RESIDUAL: 0.00
Catch the common trap
Explain before calculating.
For l = 0 the two Frobenius roots of the radial equation at r = 0 are u ~ r and u ~ constant. The second gives R = u/r ~ 1/r. On what grounds is it discarded?
Choose an answer to test the model.
Practice & worked examples
Reason from the model, then test the result.
EasyHydrogen with l = 2. Write the effective potential in eV with r in Bohr radii, find the radius and the depth of its minimum, and use them to decide which of a 2d and a 3d level can exist.
- Veff(r) = −27.2114/r + 13.6057 l(l+1)/r² eV with r in a₀. For l = 2, l(l+1) = 6, and 13.6057 × 6 = 81.6342, so Veff = −27.2114/r + 81.6342/r².
- dVeff/dr = 27.2114/r² − 163.2684/r³ = 0 gives rₘᵢₙ = 163.2684/27.2114 = 6.000 a₀ = 0.3175 nm, the general result rₘᵢₙ = l(l+1)a₀/Z.
- Depth: Veff(6a₀) = −27.2114/6 + 81.6342/36 = −4.5352 + 2.2676 = −2.2676 eV, the general result −13.6057/[l(l+1)] eV.
- A level can only sit above the floor of its own well. E₃ = −13.6057/9 = −1.5117 eV, which clears −2.2676 eV by 0.7559 eV, so 3d is a legitimate state.
- E₂ = −13.6057/4 = −3.4014 eV lies 1.1338 eV below that same floor, so no 2d state exists — this is l ≤ n − 1 seen as geometry rather than as a series truncation.
AnswerThe minimum is at 6a₀ = 0.3175 nm with a floor of −2.2676 eV. The 3d level at −1.5117 eV is bound in that well; a 2d level would have to sit 1.1338 eV below the floor, so it does not exist.
MediumTake l = 2 in hydrogen and work at the origin. Find the two Frobenius exponents of the radial equation, say which survives and on what grounds, then carry the series one term further and check it against the exact 3d radial function.
- In atomic units the radial equation is u″ = [l(l+1)/r² − 2/r − 2E]u. As r → 0 the first term is O(r⁻²) and dominates the O(r⁻¹) and O(1) terms, so the leading balance for l = 2 is u″ = 6u/r².
- Substituting u = rs gives s(s−1)r(s−2) = 6r(s−2), so s(s−1) = 6 and s = 3 or s = −2 — the general roots l+1 and −l.
- The singular root fails on the half-line: u ~ r⁻² makes ∫|u|²dr ~ ∫r⁻⁴dr diverge at the origin, so for l ≥ 1 normalisation alone kills it. Keep s = 3.
- Write u = r³(1 + a₁r + …) and match the r² coefficients. The −2Eu term acts on r³ and so first appears at order r³; to order r² the equation is u″ = (6/r² − 2/r)u. The left side gives 4⋅3⋅a₁ = 12a₁, the right side 6a₁ − 2, so 6a₁ = −2 and a₁ = −1/3 per a₀ — the general result −Z/[(l+1)a₀].
- Check against the nodeless l = 2 state, 3d: u ∝ r³e(−r/3a₀) = r³(1 − r/3a₀ + …). The coefficients agree, and since a₁ comes from the recursion and never saw E, 4d and 5d share it and differ only at the next order.
Answers = 3 or s = −2; the singular root is not square-integrable, so u = r³(1 − r/3a₀ + …), matching u ∝ r³e(−r/3a₀) for 3d.
HardThe l = 0 radial Hamiltonian is discretised in atomic units on a uniform grid rᵢ = ih, i = 1…N, inside a box of radius 40a₀, with u(0) = 0 imposed by omitting the i = 0 node. The lowest eigenvalue comes out −0.498756 Ha at h = 0.1a₀ and −0.499688 Ha at h = 0.05a₀. Confirm the order of the method, extrapolate, and explain why the same grid does far better for l = 1.
- In atomic units the matrix is tridiagonal: diagonal 1/h² − 1/rᵢ, off-diagonal −1/(2h²). The boundary condition is not a penalty term or an extra equation — it is the row that was never built, at i = 0.
- Errors against the exact −0.5 Ha are +1.244 × 10⁻³ Ha and +3.121 × 10⁻⁴ Ha. Their ratio is 3.99, so the three-point second difference is second order in h, and the convergence is from above.
- Richardson: [4E(h/2) − E(h)]/3 = [4(−0.499688) − (−0.498756)]/3 = −0.499998 Ha = −13.6057 eV, within 4 × 10⁻⁵ eV of the exact value.
- Unextrapolated, the h = 0.05a₀ grid gives −13.5972 eV, which is 0.0085 eV high — nearly 200 times the 4.5 × 10⁻⁵ eV fine-structure splitting of 2p, so an unconverged radial grid would bury any relativistic correction added on top of it.
- At l = 1 the same box and h = 0.05a₀ give −0.125007 Ha, an error of −6.5 × 10⁻⁶ Ha: 48 times smaller and of the opposite sign. The l = 0 solution meets the wall with u′(0) ≠ 0 and the 1/r singularity sitting on the first grid point, while an l = 1 solution leaves the origin as r² and gives the difference formula almost nothing to spoil.
AnswerSecond order confirmed by the 3.99 error ratio; Richardson gives −0.499998 Ha = −13.6057 eV. The l = 0 error is the origin's doing, and the same grid errs 48 times less at l = 1.