Why naive elimination fails at order (v/c)²

Units are SI; where four-vectors appear the metric is (+, −, −, −). Take A = 0 and a static potential energy V(r), so p = −iħ∇. Quoted from topic 6.01 in one line: in a stationary state with E = mc² + E′ the Dirac rows are (E′ − V)φ = cσ⋅pχ and (2mc² + E′ − V)χ = cσ⋅pφ, both exact. Topic 6.01 set the bracket to 2mc², which captures g = 2 and nothing of relative order (v/c)², where hydrogen's measured 10.97 GHz 2P₃/₂ to 2P₁/₂ interval sits. Expand the inverse one step further: χ = (1/2mc)[1 − (E′ − V)/2mc²]σ⋅pφ, the bracket an operator standing between the two σ⋅p because V(r) does not commute with p. Substituting into the upper row gives, derived, (E′ − V)φ = (p²/2m)φ − (1/4m²c²)σ⋅p(E′ − V)σ⋅pφ. The eigenvalue E′ now sits inside the operator as well: an energy-dependent effective Hamiltonian, Hermitian for each real E′, but not yet an ordinary eigenvalue problem.

To reach an ordinary eigenvalue equation, remove E′ from the correction with the lowest-order relation (E′ − V)φ = (p²/2m)φ. First move E′ − V past σ⋅p: (E′ − V)σ⋅p = σ⋅p(E′ − V) + [σ⋅p, V], with [σ⋅p, V] = −iħσ⋅∇V. Then σ⋅p(E′ − V)σ⋅pφ = (p⁴/2m)φ − iħ(σ⋅p)(σ⋅∇V)φ. The identity (σ⋅a)(σ⋅b) = a⋅b + iσ⋅(a × b) survives only with the order kept, because p differentiates ∇V: (σ⋅p)(σ⋅∇V)φ = [−iħ∇²V + ∇V⋅p − iσ⋅(∇V × p)]φ. Collecting, derived: Hₙ = p²/2m + V − p⁴/8m³c² + (ħ/4m²c²)σ⋅(∇V × p) + (ħ²/4m²c²)∇²V + (iħ/4m²c²)∇V⋅p. The first correction is the next term of √(p²c² + m²c⁴); for central V the second is (1/2m²c²)(1/r)(dV/dr)S⋅L, Thomas half included. The last two are the trouble: (∇V⋅p)† = p⋅∇V = ∇V⋅p − iħ∇²V, so Hₙ − Hₙ† = (iħ/4m²c²)(∇V⋅p + p⋅∇V), nonzero wherever a force acts.

A non-Hermitian operator invites the fear of complex energies, so test it. In any bound eigenstate |b⟩ of H₀ = p²/2m + V, ⟨b|[p², H₀]|b⟩ = 0, and [p², H₀] = [p², V] = −ħ²∇²V − 2iħ∇V⋅p, so ⟨b|iħ∇V⋅p|b⟩ = −½ħ²⟨b|∇²V|b⟩: real, derived. First-order perturbation theory, quoted from Quantum Mechanics II, still applies because H₀ is Hermitian, and the two suspect terms give (ħ²/4m²c²)⟨∇²V⟩ − (ħ²/8m²c²)⟨∇²V⟩ = (ħ²/8m²c²)⟨∇²V⟩, one Darwin energy. For hydrogen ∇²V = 4πke²δ³(r), and with Quantum Mechanics II's |ψ₂₀₀(0)|² = 1/8πa₀³ the Laplacian term gives mc²α⁴/8 = 43.8 GHz on 2s and the gradient term −21.9 GHz, a sum of 21.9 GHz; for ℓ ≥ 1 both vanish. So Hₙ reproduces the n- and j-only Eₙⱼ, with 2S₁/₂ level with 2P₁/₂, which experiment separates by only 1057.8 MHz. Diagonalising the one-dimensional Dirac equation in a harmonic well, in a basis that keeps V far below 2mc² (an electrostatic well binds a Dirac electron only quasi-stably), confirms it numerically: agreement to order c⁻⁴, and a miss of one Darwin shift if the gradient term is dropped.

The Dirac equation conserves ∫ψ†ψ d³r = ∫(φ†φ + χ†χ) d³r. With χ ≈ σ⋅pφ/2mc and (σ⋅p)² = p², the second integral is ⟨φ|p²|φ⟩/4m²c², so the conserved norm written with φ is ⟨φ|1 + p²/4m²c²|φ⟩, derived. Let Hₙ drive time evolution, legitimate for static V because every superposition of its stationary states obeys iħ∂φ/∂t = Hₙφ. Then d⟨φ|φ⟩/dt = (i/ħ)⟨φ|Hₙ† − Hₙ|φ⟩ = (1/4m²c²)⟨∇V⋅p + p⋅∇V⟩, while Ehrenfest's theorem gives d⟨p²⟩/dt = −⟨∇V⋅p + p⋅∇V⟩, the quantum form of d(p²)/dt = 2p⋅F. Their combination ⟨φ|φ⟩ + ⟨p²⟩/4m²c² is constant: Hₙ conserves the Dirac norm, not ⟨φ|φ⟩. As the potential speeds the electron up, probability passes into χ, and an operator on φ alone can record that only by being non-Hermitian. The rule is general: eliminating part of a Hilbert space leaves an energy-dependent H_eff(E′), Hermitian at each E′, whose solutions are orthogonal only with the weight 1 − ∂H_eff/∂E′, here 1 + p²/4m²c².

Define Φ = Ωφ with Ω = 1 + p²/8m²c², Hermitian, so that ⟨Φ|Φ⟩ = ⟨φ|Ω²|φ⟩ = ⟨φ|1 + p²/4m²c²|φ⟩ to this order: Φ carries the whole conserved density, so its Hamiltonian must be Hermitian. From Hₙφ = E′φ follows ΩHₙΩ⁻¹Φ = E′Φ, so H = ΩHₙΩ⁻¹ = Hₙ + [p², Hₙ]/8m²c², and at this order only [p², V] = −ħ²∇²V − 2iħ∇V⋅p survives, derived from p = −iħ∇. The commutator cancels (iħ/4m²c²)∇V⋅p exactly and cuts (ħ²/4m²c²)∇²V to (ħ²/8m²c²)∇²V, giving the Hermitian H = p²/2m + V − p⁴/8m³c² + (ħ/4m²c²)σ⋅(∇V × p) + (ħ²/8m²c²)∇²V, derived. Similar operators share eigenvalues, so not one level moves. Checks: ħ²∇²V/m²c² has units (J s)²(J m⁻²)/(kg m s⁻¹)² = J; with V = 0, H = p²/2m − p⁴/8m³c², the expansion of √(p²c² + m²c⁴) − mc²; and for hydrogen the last term is (πħ²ke²/2m²c²)δ³(r), 21.9 GHz on 2s.

The repair worked because the weight was 1 + p²/4m²c², whose square root could be read off. One order higher it is 1 + p²/4m²c² − σ⋅p(E′ − V)σ⋅p/4m³c⁴, containing V and spin: Ω must then be assembled from non-commuting operators, E′ removed again, and Hermiticity checked by hand. With time-dependent fields the procedure fails before it starts, since E′ becomes iħ∂/∂t and there is no number to trade. Even symmetrising, ½(Hₙ + Hₙ†), which happens to return H here because Hₙ − H ≈ −[p²/8m²c², H₀] is purely anti-Hermitian, is an accident of first order. The stakes are lopsided: in hydrogen's 2s state Ω − 1 averages α²/32 = 1.7 × 10⁻⁶, yet it rearranges operator terms at the 22 GHz scale, because a relative change of order (v/c)² in the wavefunction is the same order as fine structure. What is wanted is a change of representation that is unitary by construction, conserving ψ†ψ at every order and valid for time-dependent fields. That is the Foldy-Wouthuysen transformation of the next topic, which returns exactly this H.