The Rayleigh-Ritz Variational Bound

Derivation

The Rayleigh-Ritz Variational Bound

D-355 Home PU-401 Threads energy Depends on hermitian-operator-spectral-theorem

Statement

For a Hamiltonian \(\hat{H}\) bounded below with ground-state energy \(E_0\), any normalisable trial state \(\lvert\psi\rangle\) in the domain of \(\hat{H}\) satisfies \(\displaystyle E[\psi] \equiv \frac{\langle\psi\lvert\hat{H}\rvert\psi\rangle}{\langle\psi\lvert\psi\rangle} \ge E_0\), with equality if and only if \(\lvert\psi\rangle\) is a true ground state. The functional \(E[\psi]\) is therefore a rigorous upper bound on \(E_0\) for every admissible trial state.

Why it matters

Exact solution of the Schrödinger equation is impossible for almost every system of interest — helium, molecules, nuclei, condensed matter. The variational bound converts the intractable eigenvalue problem into a tractable minimisation: pick a family of trial states, compute \(E[\psi]\), and lower it. Because the result is a genuine bound (never an underestimate), a lower number is always a strictly better estimate of the ground state, giving an unambiguous figure of merit.

It is the theoretical foundation of the Hartree–Fock method, the variational Monte Carlo method, and density-functional theory, and it underlies the practical logic of nearly all of quantum chemistry. The bound is what makes "guess and improve" a controlled scientific procedure rather than uncontrolled fitting.

Assumptions

The spectrum of \(\hat{H}\) is bounded below.If no lowest eigenvalue \(E_0\) exists (e.g. an attractive \(-1/r^2\) potential that "falls to the centre"), there is nothing to bound and \(E[\psi]\) can be driven to \(-\infty\); the theorem is vacuous. \(\hat{H}\) is Hermitian (self-adjoint) with a complete orthonormal eigenbasis.Without the spectral theorem the expansion \(\lvert\psi\rangle=\sum_n c_n\lvert n\rangle\) fails, the eigenvalues need not be real, and "upper bound" loses meaning. The trial state lies in the form domain of \(\hat{H}\), so \(\langle\psi\lvert\hat{H}\rvert\psi\rangle\) is finite.A state with the correct normalisation but divergent kinetic energy (e.g. a cusp too sharp, \(\psi\sim r^{-1}\) at the origin in 3D) gives \(E[\psi]=+\infty\); the bound holds trivially but says nothing, and naive numerics may return a spurious finite value below \(E_0\). The trial state is normalisable, \(0<\langle\psi\lvert\psi\rangle<\infty\).A non-normalisable state (plane wave, unbound Gaussian tail) makes the Rayleigh quotient ill-defined; scattering states require a separate continuum treatment.

Derivation

\[ \hat{H}\lvert n\rangle = E_n\lvert n\rangle,\qquad \langle m\lvert n\rangle=\delta_{mn},\qquad E_0\le E_1\le E_2\le\cdots \]

By the spectral theorem for the self-adjoint operator \(\hat{H}\) (assumed prior result), a complete orthonormal eigenbasis \(\{\lvert n\rangle\}\) exists with real eigenvalues, which we order. A

\[ \lvert\psi\rangle = \sum_n c_n \lvert n\rangle,\qquad c_n=\langle n\lvert\psi\rangle \]

Completeness lets any admissible \(\lvert\psi\rangle\) be expanded in the eigenbasis; the coefficients are the projections. A

\[ \langle\psi\lvert\psi\rangle = \sum_{m,n} c_m^{*}c_n\langle m\lvert n\rangle = \sum_n \lvert c_n\rvert^{2} \]

Insert the expansion into the norm and use orthonormality to collapse the double sum. A

\[ \langle\psi\lvert\hat{H}\rvert\psi\rangle = \sum_{m,n} c_m^{*}c_n\,\langle m\lvert\hat{H}\rvert n\rangle = \sum_{m,n} c_m^{*}c_n E_n\delta_{mn} = \sum_n E_n\lvert c_n\rvert^{2} \]

Act \(\hat{H}\) on \(\lvert n\rangle\) to bring down \(E_n\), then use orthonormality. A

\[ \langle\psi\lvert\hat{H}\rvert\psi\rangle - E_0\langle\psi\lvert\psi\rangle = \sum_n E_n\lvert c_n\rvert^{2} - E_0\sum_n \lvert c_n\rvert^{2} = \sum_n (E_n - E_0)\lvert c_n\rvert^{2} \]

Subtract \(E_0\) times the norm and combine the two sums term by term — a purely algebraic rearrangement, no numbers yet. A

\[ E_n - E_0 \ge 0 \quad\text{and}\quad \lvert c_n\rvert^{2}\ge 0 \;\;\Longrightarrow\;\; \sum_n (E_n-E_0)\lvert c_n\rvert^{2} \ge 0 \]

Each summand is a product of two non-negative factors (ordering from Step 1 and modulus-squared positivity), so the sum is non-negative. This is the crux. B

\[ \langle\psi\lvert\hat{H}\rvert\psi\rangle - E_0\langle\psi\lvert\psi\rangle \ge 0 \;\;\Longrightarrow\;\; \frac{\langle\psi\lvert\hat{H}\rvert\psi\rangle}{\langle\psi\lvert\psi\rangle} \ge E_0 \]

Divide by \(\langle\psi\lvert\psi\rangle>0\) (normalisability assumption preserves the inequality direction). B

\[ E[\psi]=E_0 \iff (E_n-E_0)\lvert c_n\rvert^{2}=0\ \ \forall n \iff c_n=0\ \text{whenever}\ E_n>E_0 \]

A sum of non-negative terms vanishes only if every term vanishes; hence equality forces \(\lvert\psi\rangle\) to be a superposition of ground states alone. This establishes the "if and only if". C

Result

\[ \boxed{\,E[\psi]=\dfrac{\langle\psi\lvert\hat{H}\rvert\psi\rangle}{\langle\psi\lvert\psi\rangle}\;\ge\;E_0\,}\qquad \text{equality} \iff \lvert\psi\rangle \text{ is a ground state} \]

Reading. The energy expectation value of any guess you write down can never dip below the true ground-state energy. Mixing in excited states \((E_n>E_0)\) can only raise the average, weighted by how much of each you admit \((\lvert c_n\rvert^2)\). So minimising \(E[\psi]\) over a trial family drives you toward — and never past — \(E_0\) from above.

Units check. \(\langle\psi\lvert\hat{H}\rvert\psi\rangle\) carries units of (state-norm)²×energy and \(\langle\psi\lvert\psi\rangle\) carries units of (state-norm)²; their ratio is an energy (J or eV), matching \(E_0\). The bound compares like with like — energy \(\ge\) energy — and is invariant under rescaling \(\lvert\psi\rangle\to\lambda\lvert\psi\rangle\).

Limiting cases

Exact ground state: \(\lvert\psi\rangle=\lvert 0\rangle\) gives \(c_0=1\), all other \(c_n=0\), and \(E[\psi]=E_0\) exactly — the bound is saturated.
Pure excited state: \(\lvert\psi\rangle=\lvert k\rangle\) gives \(E[\psi]=E_k\ge E_0\); the bound holds but is loose by \(E_k-E_0\).
Small admixture: if \(\lvert\psi\rangle=\lvert 0\rangle+\epsilon\lvert 1\rangle\), then \(E[\psi]-E_0=\epsilon^2(E_1-E_0)/(1+\epsilon^2)=O(\epsilon^2)\) — the energy error is second order in a first-order wavefunction error. This quadratic insensitivity is why crude trial states give surprisingly good energies.
Orthogonal to the ground state: if \(\langle 0\lvert\psi\rangle=0\) is imposed, the same argument bounds the first excited energy: \(E[\psi]\ge E_1\) (the MacDonald / Hylleraas–Undheim extension).

Breaks when

The Hamiltonian is unbounded below. For a potential steeper than \(-1/r^2\) near the origin (or a "fall to the centre" model) there is no \(E_0\); \(E[\psi]\) has no floor and can be sent to \(-\infty\), so the bound is meaningless.
The trial state leaves the form domain. If \(\langle\psi\lvert\hat{H}\rvert\psi\rangle\) diverges (kinetic energy infinite from an over-sharp cusp or a discontinuity), \(E[\psi]=+\infty\); worse, a truncated numerical estimate can report a finite value spuriously below \(E_0\), a false "variational collapse".
The operator is not self-adjoint. For a non-Hermitian effective \(\hat{H}\) (open systems, complex absorbing potentials, PT-symmetric models above the symmetry-breaking threshold) eigenvalues become complex, "\(\ge\)" is undefined, and the derivation fails at Step 1.
A basis that cannot represent the true kinetic energy. Relativistic Dirac Hamiltonians are unbounded below (negative-energy continuum); naive minimisation dives into the negative sea unless kinetic-balance constraints are imposed — the classic "variational collapse" of relativistic quantum chemistry.

Failure modes

Forgetting to normalise: computing \(\langle\psi\lvert\hat{H}\rvert\psi\rangle\) alone and comparing to \(E_0\) without dividing by \(\langle\psi\lvert\psi\rangle\). The bound is on the quotient, not the numerator.
Believing lower is truer for the wavefunction: a lower \(E[\psi]\) means a better energy, not necessarily a better wavefunction — other observables (e.g. \(\langle r\rangle\)) can be worse even as the energy improves, because the error is \(O(\epsilon^2)\) in energy but \(O(\epsilon)\) in the state.
Minimising over a non-linear parameter and forgetting to check it is a minimum: a stationary point of \(E[\psi]\) can be a saddle; only the global minimum over the trial family is the best bound.
Applying it to an excited state without the orthogonality constraint: \(E[\psi]\ge E_1\) requires \(\langle 0\lvert\psi\rangle=0\) enforced exactly; approximate orthogonality can drop the estimate below \(E_1\).
Using a variational parameter inside the potential that also appears in the measure: double-counting a scale parameter (e.g. in the norm and the exponent) and mis-differentiating \(\partial E/\partial\alpha\).
Trusting a Dirac/relativistic energy that sits in the negative continuum as if it were a bound state — variational collapse mistaken for a good result.

Discussion

The theorem's power lies in what Step 5 isolates: the deviation \(E[\psi]-E_0\) is exactly \(\sum_n(E_n-E_0)\lvert c_n\rvert^2/\langle\psi\lvert\psi\rangle\), a positively weighted average of excitation energies. This is not merely a bound but a quantitative error formula. It shows the bound is tightest when the trial state's spectral weight is concentrated near the ground state, and it explains the celebrated second-order behaviour: an error \(\epsilon\) in the wavefunction produces an error \(\sim\epsilon^2\) in the energy, so even mediocre wavefunctions yield respectable energies. This asymmetry is the engine of practical quantum chemistry.

The Rayleigh–Ritz principle recasts the eigenvalue problem as stationarity of the functional \(E[\psi]\): demanding \(\delta E=0\) subject to normalisation, with a Lagrange multiplier \(\lambda\), reproduces \(\hat{H}\lvert\psi\rangle=\lambda\lvert\psi\rangle\) with \(\lambda=E\). The Schrödinger equation is thus the Euler–Lagrange equation of the energy functional. Expanding \(\lvert\psi\rangle\) in a finite basis turns stationarity into a matrix generalised-eigenvalue problem \(\mathbf{H}\mathbf{c}=E\,\mathbf{S}\mathbf{c}\), the linear variational method; the Hylleraas–Undheim–MacDonald theorem then guarantees that the \(k\)-th ordered root bounds \(E_{k-1}\) from above, and that enlarging the basis can only lower each root — the interleaving theorem.

At the rigorous level the equality condition (Step 8) requires care in a degenerate ground manifold: saturation forces \(\lvert\psi\rangle\) into the degenerate eigenspace but does not select a particular member, so the ground state recovered by minimisation is defined only up to the degeneracy. Domain issues are equally sharp: the algebraic manipulation presumes \(\lvert\psi\rangle\) lies in the form domain \(Q(\hat H)\) (finite \(\langle\psi\lvert\hat H\rvert\psi\rangle\)), which is strictly larger than the operator domain. The min–max (Courant–Fischer) theorem, \(E_k=\min_{\dim V=k+1}\max_{0\neq\psi\in V}E[\psi]\), generalises the bound to the whole discrete spectrum and is the tool used to prove that adding basis functions produces monotone convergence of every level.

Common misconceptions. The bound is on the energy, not on any other expectation value, and not on the wavefunction quality directly; a lower energy does not certify a better density or a better dipole. Nor does the variational estimate ever underestimate \(E_0\) when done correctly — an apparent underestimate signals a domain violation (infinite true energy hidden by truncation) or a non-Hermitian/relativistic Hamiltonian, not a triumph.

Worked examples

\[ \hat{H}=-\frac{\hbar^{2}}{2m}\frac{d^{2}}{dx^{2}}+\frac{1}{2}m\omega^{2}x^{2},\qquad \psi_\alpha(x)=e^{-\alpha x^{2}} \]

1-D harmonic oscillator with a Gaussian trial state of variational width \(\alpha>0\). Exact answer \(E_0=\tfrac12\hbar\omega\) known for comparison. A

\[ \langle\psi\lvert\psi\rangle=\int_{-\infty}^{\infty}e^{-2\alpha x^{2}}dx=\sqrt{\frac{\pi}{2\alpha}},\qquad \langle x^{2}\rangle=\frac{1}{4\alpha} \]

Standard Gaussian integrals; \(\langle x^2\rangle=\int x^2 e^{-2\alpha x^2}dx/\int e^{-2\alpha x^2}dx\). A

\[ \langle T\rangle=\frac{\hbar^{2}}{2m}\frac{\int\lvert\psi'\rvert^{2}dx}{\langle\psi\lvert\psi\rangle}=\frac{\hbar^{2}\alpha}{2m},\qquad \langle V\rangle=\frac12 m\omega^{2}\langle x^{2}\rangle=\frac{m\omega^{2}}{8\alpha} \]

Use \(\psi'=-2\alpha x\,e^{-\alpha x^2}\) and integrate by parts for \(\langle T\rangle\); insert \(\langle x^2\rangle\) for \(\langle V\rangle\). B

\[ E(\alpha)=\frac{\hbar^{2}\alpha}{2m}+\frac{m\omega^{2}}{8\alpha},\qquad \frac{dE}{d\alpha}=\frac{\hbar^{2}}{2m}-\frac{m\omega^{2}}{8\alpha^{2}}=0 \;\Rightarrow\; \alpha_{*}=\frac{m\omega}{2\hbar} \]

Minimise the bound over \(\alpha\); the stationary point is a minimum since \(d^2E/d\alpha^2>0\). B

\[ E(\alpha_{*})=\frac{\hbar^{2}}{2m}\cdot\frac{m\omega}{2\hbar}+\frac{m\omega^{2}}{8}\cdot\frac{2\hbar}{m\omega}=\frac{\hbar\omega}{4}+\frac{\hbar\omega}{4}=\frac{1}{2}\hbar\omega \]

Substitute \(\alpha_*\). Numbers: for an electron with \(\omega=1.0\times10^{15}\,\text{s}^{-1}\), \(\tfrac12\hbar\omega=\tfrac12(1.055\times10^{-34})(10^{15})\,\text{J}=5.3\times10^{-20}\,\text{J}=0.33\,\text{eV}\). B

\[ E_{\min}=\tfrac12\hbar\omega = E_0 \]

Reading. Because the exact ground state is a Gaussian, the trial family contains the true answer and the bound is saturated exactly — a rare case where variation is exact.

Units check. \(\hbar\omega\): \([\text{J·s}][\text{s}^{-1}]=\text{J}\). Correct energy.

\[ \hat{H}=-\frac{\hbar^{2}}{2m}\nabla^{2}-\frac{e^{2}}{4\pi\varepsilon_0 r},\qquad \psi_a(r)=e^{-r/a} \]

Hydrogen atom with a hydrogenic trial state of variational Bohr-like radius \(a\); \(E_0=-13.6\,\text{eV}\) known. A

\[ \langle T\rangle=\frac{\hbar^{2}}{2m a^{2}},\qquad \langle V\rangle=-\frac{e^{2}}{4\pi\varepsilon_0}\,\frac{1}{a} \]

Standard results for the \(1s\)-type exponential: kinetic energy from \(\int\lvert\nabla\psi\rvert^2\), potential from \(\langle 1/r\rangle=1/a\). Let \(k\equiv e^2/4\pi\varepsilon_0\). B

\[ E(a)=\frac{\hbar^{2}}{2m a^{2}}-\frac{k}{a},\qquad \frac{dE}{da}=-\frac{\hbar^{2}}{m a^{3}}+\frac{k}{a^{2}}=0 \;\Rightarrow\; a_{*}=\frac{\hbar^{2}}{m k}=a_0 \]

Minimise; the stationary radius is exactly the Bohr radius \(a_0\). B

\[ E(a_{*})=\frac{\hbar^{2}}{2m a_0^{2}}-\frac{k}{a_0}=\frac{k}{2a_0}-\frac{k}{a_0}=-\frac{k}{2a_0}=-\frac{m k^{2}}{2\hbar^{2}} \]

Substitute \(a_0=\hbar^2/mk\) and simplify — the virial relation \(\langle T\rangle=-\tfrac12\langle V\rangle\) is visible. B

\[ E_{\min}=-\frac{m k^{2}}{2\hbar^{2}} = -\frac{(9.109\times10^{-31})(2.307\times10^{-28})^{2}}{2(1.055\times10^{-34})^{2}}\,\text{J} = -2.18\times10^{-18}\,\text{J} \]

With \(k=e^2/4\pi\varepsilon_0=(1.602\times10^{-19})^2/(1.113\times10^{-10})=2.307\times10^{-28}\,\text{J·m}\). Convert: \(-2.18\times10^{-18}/1.602\times10^{-19}=-13.6\,\text{eV}\). C

\[ E_{\min} = -13.6\ \text{eV} = E_0 \]

Reading. Again the trial family contains the exact ground state, so the bound saturates. Had we instead used a Gaussian \(e^{-\beta r^2}\) (wrong cusp), we would get \(E_{\min}\approx-11.5\,\text{eV}>E_0\) — above the true energy, as the theorem guarantees.

Units check. \(mk^2/\hbar^2\): \([\text{kg}][\text{J·m}]^2/[\text{J·s}]^2=[\text{kg}][\text{J}^2\text{m}^2]/[\text{J}^2\text{s}^2]=[\text{kg}\,\text{m}^2\text{s}^{-2}]=\text{J}\). Correct.

Problems

State precisely the equality condition of the variational theorem and prove that equality \(E[\psi]=E_0\) forces \(\lvert\psi\rangle\) to lie entirely in the ground-state eigenspace.
Solution
From Step 5, \(E[\psi]-E_0=\sum_n(E_n-E_0)\lvert c_n\rvert^2/\sum_n\lvert c_n\rvert^2\). Each numerator term is \(\ge 0\). A sum of non-negative terms equals zero iff every term is zero, so \((E_n-E_0)\lvert c_n\rvert^2=0\) for all \(n\). For any \(n\) with \(E_n>E_0\) this forces \(c_n=0\); only coefficients with \(E_n=E_0\) may be nonzero. Hence \(\lvert\psi\rangle=\sum_{n:E_n=E_0}c_n\lvert n\rangle\) lies in the ground eigenspace. Conversely any such state gives \(E[\psi]=E_0\). \(\square\)
Using the trial state \(\psi(x)=(a^2-x^2)\) on \([-a,a]\) (zero outside) for the infinite square well of width \(2a\), estimate the ground-state energy and compare to the exact \(E_0=\pi^2\hbar^2/(2m(2a)^2)\).
Solution
\(\langle\psi\lvert\psi\rangle=\int_{-a}^{a}(a^2-x^2)^2dx=\tfrac{16}{15}a^5\). \(\psi'=-2x\), so \(\langle T\rangle=\frac{\hbar^2}{2m}\int_{-a}^{a}(2x)^2dx=\frac{\hbar^2}{2m}\cdot\frac{8}{3}a^3\). Then \(E[\psi]=\frac{\hbar^2}{2m}\cdot\frac{(8/3)a^3}{(16/15)a^5}=\frac{\hbar^2}{2m}\cdot\frac{5}{2a^2}=\frac{5\hbar^2}{4ma^2}\). Exact: \(E_0=\frac{\pi^2\hbar^2}{8ma^2}=\frac{1.2337\hbar^2}{ma^2}\). Estimate \(=1.25\,\hbar^2/(ma^2)\). The bound exceeds the exact value by \(0.13\%\) — above \(E_0\), as required, and remarkably tight.
Show that if \(\lvert\psi\rangle=\lvert 0\rangle+\epsilon\lvert 1\rangle\) (unnormalised, \(\epsilon\) real, \(\langle 0\lvert 1\rangle=0\)), the energy error is second order in \(\epsilon\).
Solution
\(\langle\psi\lvert\psi\rangle=1+\epsilon^2\). \(\langle\psi\lvert\hat H\rvert\psi\rangle=E_0+\epsilon^2 E_1\). So \(E[\psi]=\frac{E_0+\epsilon^2E_1}{1+\epsilon^2}=E_0+\frac{\epsilon^2(E_1-E_0)}{1+\epsilon^2}=E_0+\epsilon^2(E_1-E_0)+O(\epsilon^4)\). The leading error is \(\propto\epsilon^2\), i.e. quadratic in the (first-order) wavefunction error, and positive since \(E_1>E_0\).
For the hydrogen atom use a Gaussian trial state \(\psi(r)=e^{-\beta r^2}\). Find the optimal \(\beta\) and the resulting energy in eV, and confirm it exceeds \(-13.6\,\text{eV}\).
Solution
Standard results: \(\langle T\rangle=\frac{3\hbar^2\beta}{2m}\), \(\langle V\rangle=-k\langle 1/r\rangle=-k\sqrt{8\beta/\pi}\) with \(k=e^2/4\pi\varepsilon_0\). \(E(\beta)=\frac{3\hbar^2\beta}{2m}-k\sqrt{\frac{8\beta}{\pi}}\). Set \(dE/d\beta=0\): \(\frac{3\hbar^2}{2m}=k\sqrt{\frac{2}{\pi\beta}}\), giving \(\beta_*=\frac{8m^2k^2}{9\pi\hbar^4}\). Back-substituting yields \(E_{\min}=-\frac{8}{3\pi}\cdot\frac{mk^2}{2\hbar^2}=-\frac{8}{3\pi}\,\text{Ry}\), where \(\text{Ry}=mk^2/2\hbar^2=13.6\,\text{eV}\). Numerically \(E_{\min}=-\frac{8}{3\pi}(13.6)\,\text{eV}=-11.5\,\text{eV}\). Since \(-11.5>-13.6\), the Gaussian (wrong cusp behaviour, \(\psi'(0)=0\) instead of a cusp) gives a looser bound, exactly as the theorem demands.
The Hylleraas–Undheim–MacDonald theorem: prove that if \(\lvert\psi\rangle\) is constrained orthogonal to the exact ground state, \(\langle 0\lvert\psi\rangle=0\), then \(E[\psi]\ge E_1\).
Solution
Orthogonality kills the \(n=0\) coefficient: \(c_0=\langle 0\lvert\psi\rangle=0\), so \(\lvert\psi\rangle=\sum_{n\ge1}c_n\lvert n\rangle\). Then \(E[\psi]-E_1=\frac{\sum_{n\ge1}(E_n-E_1)\lvert c_n\rvert^2}{\sum_{n\ge1}\lvert c_n\rvert^2}\). For \(n\ge1\), \(E_n\ge E_1\), so every term is \(\ge0\) and \(E[\psi]\ge E_1\), with equality iff \(\lvert\psi\rangle\) lies in the first-excited eigenspace. In practice one uses the best available approximation to \(\lvert 0\rangle\); errors in it lower the bound below \(E_1\), which is why the constraint must be exact for a rigorous excited-state bound.