The Variational Method

Why it matters

The Schrödinger equation can be solved exactly for only a handful of systems — the particle in a box, the harmonic oscillator, and the hydrogen atom among them. The moment you add a second electron, as in the helium atom, the electron-electron repulsion term makes the equation analytically intractable. The variational method is the escape hatch: instead of solving H ψ = E ψ, you guess a plausible wavefunction and let the mathematics guarantee that your energy estimate is never too low.

That guarantee is what makes the method so powerful. Many approximation schemes give an answer without telling you which side of the truth you land on. The variational principle always errs on the high side: your computed energy is a rigorous upper bound on the true ground-state energy E₀. So if two different trial functions give −77 eV and −78 eV for helium, you know immediately that −78 eV is the better estimate — no external check required. The lower your bound, the closer you are to the truth.

This single fact scales all the way up to industrial computation. Every time a chemist runs a Hartree-Fock or density-functional calculation to predict a bond length, a reaction barrier, or a spectrum, they are — under the hood — minimizing ⟨H⟩ over a parametrized family of wavefunctions. The variational method is not a niche trick; it is the conceptual engine of computational chemistry and much of condensed-matter physics.

The variational theorem

The central statement is compact. For a Hamiltonian H with ground-state energy E₀ and any (not necessarily normalized) trial state ψ:

E[ψ] = ⟨ψ|H|ψ⟩ / ⟨ψ|ψ⟩ ≥ E₀

where:

• E[ψ] — the energy functional (the Rayleigh quotient), in joules or, more commonly in atomic physics, in Hartree (1 Ha = 27.211 eV).
• ⟨ψ|H|ψ⟩ — the expectation value of the Hamiltonian, ∫ ψ* H ψ dτ, integrated over all space.
• ⟨ψ|ψ⟩ — the normalization integral, ∫ ψ* ψ dτ; dividing by it removes any dependence on the overall scale of ψ.
• E₀ — the exact ground-state (lowest) eigenvalue of H.

Equality holds only when ψ is the true ground-state eigenfunction φ₀. Any other guess overestimates the energy.

How it works, step by step

1. Choose a trial wavefunction. Pick a function ψ(x; α, β, …) that satisfies the boundary conditions and captures the physics — it should decay where the potential is large, and have the right symmetry. Introduce one or more variational parameters (a width, an effective charge, mixing coefficients).

2. Build the energy functional. Compute the two integrals and form E(α, β, …) = ⟨ψ|H|ψ⟩ / ⟨ψ|ψ⟩. The parameters now appear as ordinary variables in a function you can differentiate.

3. Minimize. Set the gradient to zero:

∂E/∂α = 0,   ∂E/∂β = 0,   …

and solve for the optimal parameters. Substituting them back gives the lowest — and therefore best — upper bound your trial family can produce.

4. Improve (optional). Add more parameters or more basis functions. Because the family of accessible states grows, the minimum can only drop (or stay put), so the bound is monotonically improving. In the limit of a complete basis it converges to E₀ exactly.

Why the bound holds — the eigenbasis proof

The proof is a two-line argument that never uses the eigenstates explicitly, only their existence. Expand the trial state in the complete (unknown) energy eigenbasis {φₙ}, with H φₙ = Eₙ φₙ and E₀ ≤ E₁ ≤ E₂ ≤ …:

ψ = Σₙ cₙ φₙ     ⟹     ⟨H⟩ = Σₙ |cₙ|² Eₙ / Σₙ |cₙ|²

Now subtract E₀ from both sides. Because every Eₙ − E₀ ≥ 0 and every |cₙ|² ≥ 0:

⟨H⟩ − E₀ = Σₙ |cₙ|² (Eₙ − E₀) / Σₙ |cₙ|² ≥ 0

so ⟨H⟩ ≥ E₀. The excess energy is a weighted average of the positive gaps (Eₙ − E₀); it vanishes only when all the weight sits on n = 0, i.e. when ψ = φ₀. Crucially, the argument needs only that the eigenstates form a complete set — you never have to find them.

Worked example: the helium atom

Helium is the classic showcase. Its Hamiltonian (in atomic units, ℏ = mₑ = e = 4πε₀ = 1) is

H = −½∇₁² − ½∇₂² − 2/r₁ − 2/r₂ + 1/r₁₂

The first two terms are the electron kinetic energies, the next two are attraction to the Z = 2 nucleus, and the last is the electron-electron repulsion — the term that makes the problem unsolvable in closed form. A natural trial function is a product of two hydrogen-like 1s orbitals, but with the nuclear charge replaced by an effective charge Z that we treat as the variational parameter:

ψ(r₁, r₂) = (Z³/π) e^(−Z r₁) e^(−Z r₂)

Carrying out the integrals gives the energy as a function of Z:

E(Z) = Z² − (27/8) Z    (Hartree)

Minimizing, dE/dZ = 2Z − 27/8 = 0, yields Z = 27/16 = 1.6875. The value is less than the true nuclear charge of 2 because each electron partially screens the nucleus from the other — the variational method discovers screening from first principles. Plugging back in:

E_min = −(27/16)² = −729/256 ≈ −2.848 Ha ≈ −77.5 eV

The experimental ground-state energy of helium is −2.9037 Ha ≈ −79.0 eV, so this one-parameter estimate is within about 2 percent — and, as guaranteed, it lies above the true value. Hylleraas (1929) added terms in the inter-electron distance r₁₂ and drove the agreement to many decimal places; his approach still sets the standard for high-precision helium calculations.

Helium trial wavefunction	Parameters	Energy (Ha)	Error vs −2.9037 Ha
Bare Z = 2 hydrogenic product (no minimization)	0	−2.750	+5.3%
Effective charge Z_eff = 27/16	1	−2.848	+1.9%
Hylleraas 3-term (includes r₁₂)	3	−2.9024	+0.04%
Hylleraas / Pekeris many-term	>1000	−2.903724	<10⁻⁶

Rayleigh-Ritz: the linear variational method

When the trial function is a linear combination of fixed basis functions, ψ = Σᵢ cᵢ χᵢ, the minimization becomes an eigenvalue problem. Substituting into E[ψ] and setting ∂E/∂cᵢ* = 0 gives the secular equation:

Σⱼ (Hᵢⱼ − E Sᵢⱼ) cⱼ = 0     ⟺     H c = E S c

where Hᵢⱼ = ⟨χᵢ|H|χⱼ⟩ is the Hamiltonian matrix and Sᵢⱼ = ⟨χᵢ|χⱼ⟩ is the overlap matrix (S = identity for an orthonormal basis). This is a generalized eigenvalue problem: its lowest root is the best variational bound on E₀. Remarkably, the higher roots are upper bounds on the excited-state energies too — the Hylleraas-Undheim-MacDonald theorem guarantees that an N-function basis brackets the lowest N eigenvalues from above. This is the exact structure that appears in Hartree-Fock's Roothaan-Hall equations and in configuration-interaction calculations.

Aspect	Variational method	Perturbation theory
Requires small parameter?	No	Yes (H = H₀ + λV, small λ)
Bound on energy?	Rigorous upper bound on E₀	None guaranteed (may over- or under-shoot)
Convergence guaranteed?	Monotonic as basis grows	Series may diverge (asymptotic)
Best for	Ground state; lowest state of a symmetry	Small corrections to a solvable H₀
Energy error given wavefunction error δψ	Second order, ∝ (δψ)²	Order of the truncated term

Common misconceptions

• "The variational method gives the exact energy." It gives an upper bound. It equals E₀ only if your trial function happens to be the exact ground state — which, for hard problems, you don't know.
• "A lower energy means a better wavefunction for everything." Not necessarily. The energy error is second order in the wavefunction error, so a state can nail the energy while giving a poor density, dipole moment, or ⟨r²⟩. Judge the wavefunction by the observable you care about.
• "It only works for the ground state." It works cleanly for the ground state, but with an orthogonality (or symmetry) constraint it bounds excited states too: a trial function orthogonal to φ₀ satisfies ⟨H⟩ ≥ E₁.
• "You must normalize the trial function first." You don't — dividing by ⟨ψ|ψ⟩ in the Rayleigh quotient handles normalization automatically, and the bound is independent of the overall scale of ψ.
• "More parameters can make the bound worse." Adding parameters (while keeping the old ones reachable) can only lower or keep the minimum, never raise it — the accessible family of states strictly grows.
• "Getting E(Z) = Z² − (27/8)Z wrong means the method failed." A wrong energy from a sensible trial function is still a valid upper bound; the method never gives you a number below E₀. If you ever compute ⟨H⟩ < E₀, you have made an algebra or normalization error.

A little history

The Rayleigh quotient dates to Lord Rayleigh's Theory of Sound (1877), where he bounded the frequencies of vibrating systems. Walther Ritz turned it into a systematic computational recipe in 1908, expanding trial functions in a basis and minimizing — the linear variational method now universally called Rayleigh-Ritz. When quantum mechanics arrived, the principle transferred immediately: the eigenvalue problem for a vibrating membrane and for the Schrödinger equation are mathematically identical. Egil Hylleraas's 1929 helium calculation was an early triumph, and by the 1950s the Roothaan-Hall formulation had made variational optimization the backbone of quantum chemistry. Today the largest supercomputers on Earth still spend much of their time doing exactly what Ritz proposed: minimizing ⟨H⟩.