Sturm-Liouville Theory and Special Functions

Confront the special functions of physics all at once — Legendre and Hermite polynomials, Laguerre and Chebyshev, the Bessel functions of the drumhead — and they look like a museum of unrelated curiosities, each with its own weird orthogonality weight, its own three-term recurrence, its own generating function pulled from a hat. Sturm–Liouville theory is the single idea that collapses the whole museum into one exhibit. Every one of those families is nothing more than the set of eigenfunctions of a second-order differential operator that has been made self-adjoint — the continuous analogue of a real symmetric matrix — on a fixed interval with a fixed weight. Once you see the operator, everything else is forced. The central object is $Ly=-\frac{1}{w}[(py')'+qy]$ , engineered so that the eigenvalue equation $Ly=\lambda y$ reproduces the ODE while $L$ is symmetric in the weighted inner product $\langle f,g\rangle_w=\int_a^b\overline{f}\,g\,w\,\mathrm{d}x$ . Symmetry is proved once, by a single integration-by-parts identity — the Lagrange (Green) identity — whose only residue is a boundary term. Make that boundary term vanish and three theorems fall out together: the eigenvalues are real, eigenfunctions for distinct eigenvalues are orthogonal in the weight $w$ , and for a regular problem they form a complete basis of $L^2_w[a,b]$ . The weight is not a decoration; it is the metric that makes the operator symmetric, and it is the same object whether you call it an orthogonality measure or a physical density. This article builds the operator framework, proves reality and orthogonality straight from the boundary identity, shows how any second-order ODE is beaten into self-adjoint form by an integrating factor, and then watches Legendre, Hermite, and Bessel functions drop into the scheme as the eigenfunctions of specific singular problems. Singular endpoints, the limit-point/limit-circle dichotomy, and continuous spectra appear as the natural bridge from this classical machinery to scattering and the full spectral theorem.

Key takeaways

A Sturm–Liouville problem writes a second-order linear ODE in the self-adjoint form $\frac{\mathrm{d}}{\mathrm{d}x}(p\,y')+(q+\lambda w)y=0$ on $[a,b]$ , with boundary conditions chosen so the operator is symmetric.
Self-adjointness forces the eigenvalues $\lambda_n$ to be real, eigenfunctions for distinct $\lambda_n$ to be orthogonal in the weight $w(x)$ , and — for regular problems — the eigenfunctions to be complete in $L^2_w[a,b]$ .
Every classical special function (Legendre, Hermite, Laguerre, Chebyshev, Bessel) is the eigenfunction family of one such operator; their orthogonality, completeness, and recurrences are all consequences, not coincidences.
The Lagrange identity $\int_a^b(uLv-vLu)\,w\,\mathrm{d}x=-[p(uv'-vu')]_a^b$ is the whole engine: making the boundary term vanish is what selects the admissible boundary conditions.
Any $Py''+Qy'+(R+\lambda)y=0$ is cast into SL form by the integrating factor $w=\frac{1}{P}\exp\int\frac{Q}{P}\,\mathrm{d}x$ , which simultaneously reveals the orthogonality weight.
Key idea: Write the ODE as a self-adjoint operator $Ly=-\frac{1}{w}[(py')'+qy]$ ; then reality of eigenvalues, $w$ -weighted orthogonality, and completeness are all forced by one boundary-term identity, and every classical special function is simply the eigenfunction family of one such operator.

The core insight

Write the ODE as a self-adjoint operator $Ly=-\frac{1}{w}[(py')'+qy]$ ; then reality of eigenvalues, $w$ -weighted orthogonality, and completeness are forced by a single boundary-term identity, and every classical special function is simply the eigenfunction family of one such operator.

Where this lives today

The Sturm–Liouville machine is quietly running underneath most of computational physics. Every separation-of-variables solution of Laplace, Helmholtz, heat, wave, or Schrödinger equations decomposes into one SL problem per coordinate, and the completeness theorem is the licence to match arbitrary boundary and initial data by a series — spherical harmonics (Legendre in $\theta$ , periodic in $\phi$ ), cylindrical modes (Bessel in $r$ ), and quantum eigenstates all arrive this way. In numerical practice the same structure powers spectral methods: representing a field in the eigenbasis of an SL operator gives exponential convergence for smooth solutions, which is why global Chebyshev and Legendre expansions dominate high-accuracy fluid and climate codes, and why Gauss–Legendre quadrature — whose nodes are the zeros of $P_\ell$ you counted in the figure — is the default for exact polynomial integration. The limit-point/limit-circle dichotomy is a live tool in quantum mechanics: it decides whether a singular potential (the Coulomb tail, an inverse-square well, the AdS boundary in holography) needs an extra self-adjoint-extension boundary condition to make the Hamiltonian a legitimate observable — an issue that resurfaces whenever someone quantises a field on a space with a boundary or a singularity. And the drumhead of Worked Example 2 is not a toy: the inharmonic Bessel spectrum is exactly what acousticians engineer around when they tune timpani and design loudspeaker diaphragms.

The physical picture

Think of $L$ as a real symmetric matrix with infinitely many rows. A symmetric matrix $A=A^{\mathsf T}$ has real eigenvalues and eigenvectors you can choose orthonormal — that is linear algebra. The Sturm–Liouville operator is the continuum version of exactly that statement: real eigenvalues $\lambda_n$ , eigenfunctions orthogonal in an inner product, a complete basis to expand any state. The one subtlety is that the "dot product" is not the naive $\int \overline f g\,\mathrm{d}x$ but the weighted $\int\overline f g\,w\,\mathrm{d}x$ ; the weight $w(x)$ is the diagonal metric that makes $L$ genuinely symmetric. Pick a weight and an interval and the machine hands you a specific special-function family — the first figure shows how one operator template spawns the entire zoo.

Figure

One self-adjoint operator, many families. Choosing the coefficient p(x), the weight w(x) and the interval selects Legendre p = 1−x² w = 1, Hermite p = w = e^−x², Laguerre p = xe^−x w = e^−x, or Bessel p = w = x as the eigenfunctions of the same template

The second face of the theory is oscillation. On a finite interval the eigenvalues are discrete and ordered, $\lambda_0<\lambda_1<\lambda_2<\cdots$ , and the eigenfunction $y_n$ carries exactly $n$ interior zeros — a wave with more nodes costs more "energy," precisely as for a vibrating string. The classical polynomials wear this theorem on their sleeve: $P_\ell$ has $\ell$ real zeros inside $(-1,1)$ , and those zeros are what make Gauss–Legendre quadrature exact. The second figure plots the first Legendre polynomials so you can literally count nodes and see the family close under a single orthogonality relation.

Figure

The first Legendre polynomials on the interval from −1 to 1. Each P_ℓ has exactly ℓ interior zeros, the oscillation theorem in action, and the family is orthogonal with ∫ P_m P_n dx = 2/(2ℓ+1) δ_mn

The core equations

The starting point is the Sturm–Liouville equation itself: a second-order linear ODE written so that the highest two derivatives collapse into one exact derivative $(py')'$ . The parameter $\lambda$ is the eigenvalue, and $w$ is the weight that will later measure orthogonality. Written this way, the "physics" (the operator) and the "geometry" (the inner product) are already entangled.

\frac{\mathrm{d}}{\mathrm{d}x}\!\left(p(x)\,\frac{\mathrm{d}y}{\mathrm{d}x}\right)+\left(q(x)+\lambda\, w(x)\right)y=0,\qquad a\le x\le b \qquad (1)

where: $p(x)$ — coefficient function, $p>0$ on the open interval (its vanishing at an endpoint makes the problem singular); $q(x)$ — real potential-like term; $w(x)$ — weight (density), $w>0$ , fixes the inner product; $\lambda$ — spectral/eigenvalue parameter; $y(x)$ — eigenfunction, subject to boundary conditions at $a,b$ . Read it as "an operator acting on $y$ equals $\lambda$ times $y$ against the weight $w$ ." Isolating that operator gives the self-adjoint form. The minus sign and the $1/w$ are conventions chosen so that $L$ is a positive, symmetric operator in the weighted inner product — the exact analogue of a real symmetric matrix acting on $\mathbb{R}^n$ .

L y=-\frac{1}{w(x)}\left[\frac{\mathrm{d}}{\mathrm{d}x}\!\left(p\,y'\right)+q\,y\right],\qquad Ly=\lambda y \qquad (2)

where: $L$ — the Sturm–Liouville operator, symmetric in $\langle\cdot,\cdot\rangle_w$ ; $\langle f,g\rangle_w=\int_a^b\overline{f(x)}\,g(x)\,w(x)\,\mathrm{d}x$ — the weighted inner product that turns $L$ into a self-adjoint object. Everything downstream rests on one calculus identity. Integrating $uLv-vLu$ by parts twice, the second-derivative pieces cancel and leave only a boundary term — the Lagrange, or Green, identity. It is the continuum version of $\langle u,Av\rangle-\langle Au,v\rangle=0$ for a symmetric matrix, but with an explicit surface remainder.

\int_a^b\!\big(u\,Lv-v\,Lu\big)\,w\,\mathrm{d}x=-\Big[\,p\,(u v'-v u')\,\Big]_a^b \qquad (3)

where: $u,v$ — any two functions in the domain of $L$ . The bracket $p(uv'-vu')$ is the bilinear concomitant; its value at the endpoints is the sole obstruction to symmetry. Self-adjointness is therefore a single scalar condition: make that boundary bracket vanish. Whether you achieve it by Dirichlet ( $y=0$ ), Neumann ( $y'=0$ ), periodic, or a singular-endpoint regularity condition, the algebraic requirement is identical.

\Big[\,p\,(u v'-v u')\,\Big]_a^b=0\quad\Longrightarrow\quad L=L^\dagger\ \text{on its domain} \qquad (4)

where: $L^\dagger$ — the adjoint of $L$ in $L^2_w$ ; equality of $L$ and $L^\dagger$ with equal domains is self-adjointness, the property that guarantees a real, complete spectrum. The payoff theorem is weighted orthogonality: eigenfunctions belonging to different eigenvalues are perpendicular in the weight $w$ . This is why Hermite polynomials are orthogonal under $e^{-x^2}\,\mathrm{d}x$ and not under $\mathrm{d}x$ — the weight is dictated by the operator, not chosen for convenience.

\langle y_m,y_n\rangle_w=\int_a^b y_m(x)\,y_n(x)\,w(x)\,\mathrm{d}x=0\quad(\lambda_m\neq\lambda_n) \qquad (5)

where: $y_m,y_n$ — eigenfunctions for distinct eigenvalues $\lambda_m\neq\lambda_n$ . Equal eigenvalues (degeneracy) leave the overlap undetermined and must be orthogonalised by hand. Most ODEs of physics do not arrive in self-adjoint form. An integrating factor fixes that: multiply $Py''+Qy'+(R+\lambda)y=0$ by the weight below and the first two terms fuse into $(py')'$ . The very same $w$ is the orthogonality measure — one calculation delivers both.

w(x)=\frac{1}{P(x)}\exp\!\left(\int^x \frac{Q(t)}{P(t)}\,\mathrm{d}t\right),\qquad p(x)=P(x)\,w(x) \qquad (6)

where: $P,Q,R$ — the original (generally non-self-adjoint) coefficients. The construction fails only where $P=0$ , which is exactly where singular endpoints and the classical polynomials live. Completeness is what makes the theory useful: for a regular problem the eigenfunctions span $L^2_w[a,b]$ , so any admissible function is a generalized Fourier series in them, with coefficients read off by weighted projection.

f(x)=\sum_{n=0}^{\infty} c_n\,y_n(x),\qquad c_n=\frac{\langle y_n,f\rangle_w}{\langle y_n,y_n\rangle_w} \qquad (7)

where: $f$ — any $L^2_w[a,b]$ function; the series converges in the $L^2_w$ norm; $c_n$ — generalized Fourier coefficient. Ordinary Fourier series ( $p=w=1$ , periodic) are the simplest special case. Completeness has a distributional shadow: the eigenfunctions resolve the identity, and the resolvent (Green's function) $G(x,\xi;\lambda)=\sum_n y_n(x)y_n(\xi)/[(\lambda_n-\lambda)\langle y_n,y_n\rangle_w]$ has poles exactly at the eigenvalues.

\sum_{n=0}^{\infty}\frac{y_n(x)\,y_n(\xi)}{\langle y_n,y_n\rangle_w}=\frac{\delta(x-\xi)}{w(\xi)} \qquad (8)

where: $\delta$ — the Dirac delta; the factor $1/w$ is what makes the closure relation consistent with the weighted inner product. Test it by hitting both sides with $\int\,\cdot\,w\,\mathrm{d}x$ and it reproduces $f=\sum c_n y_n$ .

Derivation AReality of eigenvalues and orthogonality from one boundary identity

The goal here is to extract the two cornerstone theorems of the whole subject — that Sturm–Liouville eigenvalues are real and that eigenfunctions for distinct eigenvalues are orthogonal in the weight $w$ — from nothing more than the symmetry of $L$ , i.e. the vanishing of the Lagrange boundary term. No properties of any particular special function are used; the argument is pure operator symmetry, and it is the exact continuum echo of the linear-algebra fact that a real symmetric matrix has real eigenvalues and orthogonal eigenvectors. Assume $Ly_m=\lambda_m y_m$ and $Ly_n=\lambda_n y_n$ with boundary conditions making the bilinear concomitant vanish at both ends.

Step 1. Apply Green's identity (3) to $y_m$ and $\overline{y_n}$ . Because $p,q,w$ are all real, $L$ commutes with complex conjugation, so $L\overline{y_n}=\overline{\lambda_n}\,\overline{y_n}$ ; the boundary term vanishes by hypothesis, leaving $\int_a^b(\overline{y_n}\,Ly_m-y_m\,L\overline{y_n})\,w\,\mathrm{d}x=0$ . Why: symmetry of $L$ is the only fact injected — everything else is bookkeeping on this one equation.

Step 2. Insert the eigenvalue equations. The integrand becomes $(\lambda_m-\overline{\lambda_n})\,\overline{y_n}\,y_m\,w$ , so $(\lambda_m-\overline{\lambda_n})\int_a^b\overline{y_n}\,y_m\,w\,\mathrm{d}x=0$ . Why: pulling the eigenvalues out of the operator converts an operator identity into a scalar equation about overlaps.

Step 3. Set $n=m$ . The overlap becomes the squared weighted norm $\|y_m\|_w^2=\int|y_m|^2 w\,\mathrm{d}x>0$ , so $(\lambda_m-\overline{\lambda_m})\|y_m\|_w^2=0$ forces $\lambda_m=\overline{\lambda_m}$ — the eigenvalue is real. Why: a strictly positive norm cannot be killed, so the prefactor must vanish, and that is precisely the reality statement.

Step 4. Take $m\neq n$ with $\lambda_m\neq\lambda_n$ . The eigenvalues are now known to be real, so $(\lambda_m-\lambda_n)\langle y_n,y_m\rangle_w=0$ , and since the prefactor is nonzero the overlap must vanish: $\langle y_n,y_m\rangle_w=0$ . Why: distinct real eigenvalues cannot share a weighted inner product; degenerate eigenspaces are simply Gram–Schmidt-orthogonalised within. *Verified: reproduces the symmetric-matrix theorem in the $w$ -weighted inner product; standard, Arfken–Weber–Harris §8.2, Riley–Hobson–Bence §17.4.*

Go deeper · Derivation B · 1 minLegendre's equation as a singular Sturm–Liouville problem

Derivation BLegendre's equation as a singular Sturm–Liouville problem

The aim is to watch a named special function crystallise out of the operator framework with no extra assumptions: Legendre's equation is already in self-adjoint form, its coefficient $p=1-x^2$ conveniently vanishes at both endpoints of $[-1,1]$ , and that single fact makes the boundary term disappear for free, quantises the spectrum to $\lambda=\ell(\ell+1)$ , and hands back the orthogonal polynomials $P_\ell$ with their textbook normalisation. This is the prototype for every classical-polynomial family born at a singular endpoint. Start from Legendre's equation $(1-x^2)y''-2xy'+\ell(\ell+1)y=0$ .

Step 1. Read off the SL data. Since $\frac{\mathrm{d}}{\mathrm{d}x}[(1-x^2)y']=(1-x^2)y''-2xy'$ matches the first two terms exactly, the equation is $\frac{\mathrm{d}}{\mathrm{d}x}[(1-x^2)y']+\ell(\ell+1)y=0$ , i.e. $p=1-x^2$ , $q=0$ , $w=1$ , $\lambda=\ell(\ell+1)$ . Why: no integrating factor is needed — Legendre's equation is self-adjoint on the nose, which is why its weight is trivially $w=1$ .

Step 2. Kill the boundary term. The concomitant is $[(1-x^2)(uv'-vu')]_{-1}^{1}$ , and because $p(\pm1)=1-x^2=0$ it vanishes automatically for any $u,v$ that stay finite at the ends. Why: at a singular endpoint the vanishing of $p$ does the job a boundary condition would otherwise do — mere boundedness is the whole requirement.

Step 3. Quantise by regularity. The Frobenius solution about $x=\pm1$ is generically singular (the second solution $Q_\ell$ blows up logarithmically); it terminates into a finite polynomial $P_\ell$ only when $\lambda=\ell(\ell+1)$ for integer $\ell=0,1,2,\dots$ . Why: demanding a bounded eigenfunction at the singular ends is exactly what selects the discrete spectrum — regularity, not an imposed value, does the quantising.

Step 4. Collect orthogonality and normalisation. The SL theorem gives $\int_{-1}^{1}P_m P_n\,\mathrm{d}x=0$ for $m\neq n$ , and Rodrigues' formula $P_\ell=\frac{1}{2^\ell \ell!}\frac{\mathrm{d}^\ell}{\mathrm{d}x^\ell}(x^2-1)^\ell$ fixes the constant to $\frac{2}{2\ell+1}$ . Why: orthogonality is free from self-adjointness, but the actual value of $\|P_\ell\|^2$ needs the explicit Rodrigues representation. *Verified: $\int_{-1}^1 P_\ell^2\,\mathrm{d}x=2/(2\ell+1)$ checked numerically for $\ell=0,1,2,3$ ; standard, Arfken §15.1.*

Go deeper · Derivation C · 1 minCasting Hermite and Bessel into Sturm–Liouville form

Derivation CCasting Hermite and Bessel into Sturm–Liouville form

The purpose of this derivation is to run the integrating-factor machine (6) on two equations that do not start self-adjoint, and to show that the weight it produces is never arbitrary: for Hermite it is the Gaussian $e^{-x^2}$ that every physicist already associates with the harmonic oscillator, and for Bessel it is $w=x$ , the radial Jacobian $r\,\mathrm{d}r$ of the plane. In both cases the SL weight turns out to be the physical measure, which is the deep reason these functions are the natural bases for their respective problems. Take the generic form $Py''+Qy'+(R+\lambda)y=0$ .

Step 1. Hermite. From $y''-2xy'+2ny=0$ we have $P=1$ , $Q=-2x$ , so $w=\exp\int(-2x)\,\mathrm{d}x=e^{-x^2}$ and $p=wP=e^{-x^2}$ ; the self-adjoint form is $(e^{-x^2}y')'+2n\,e^{-x^2}y=0$ with eigenvalue $\lambda=2n$ . Why: the integrating factor manufactures precisely the Gaussian weight, so Hermite orthogonality under $e^{-x^2}\,\mathrm{d}x$ is not a definition but a consequence.

Step 2. Hermite orthogonality survives on the whole line. The weight $e^{-x^2}$ decays faster than any polynomial grows, so the concomitant $[e^{-x^2}(uv'-vu')]_{-\infty}^{\infty}$ vanishes and the spectrum stays discrete, giving $\int_{-\infty}^{\infty}H_m H_n\,e^{-x^2}\,\mathrm{d}x=2^n n!\sqrt{\pi}\,\delta_{mn}$ . Why: on an infinite interval it is the decaying weight, not a hard wall, that supplies self-adjointness — the boundary term dies at infinity.

Step 3. Bessel. From $x^2y''+xy'+(k^2x^2-\nu^2)y=0$ , divide by $x$ to get $\frac{\mathrm{d}}{\mathrm{d}x}(xy')+(k^2x-\nu^2/x)y=0$ , so $p=x$ , $q=-\nu^2/x$ , $w=x$ , and the eigenvalue is $\lambda=k^2$ . The point $x=0$ , where $p=0$ , is a singular endpoint. Why: dividing by $x$ is the integrating-factor step in disguise, and it exposes $w=x$ — the same $r$ that appears in the area element of the disk.

Step 4. Bessel orthogonality on a disk of radius $a$ . Imposing $y(a)=0$ (a rim clamp) sets $k_n=\alpha_{\nu,n}/a$ with $\alpha_{\nu,n}$ the $n$ -th positive zero of $J_\nu$ , and regularity at $x=0$ supplies the inner condition, giving $\int_0^a J_\nu(k_m x)J_\nu(k_n x)\,x\,\mathrm{d}x=\frac{a^2}{2}[J_{\nu+1}(k_n a)]^2\,\delta_{mn}$ . Why: the weight $w=x$ is exactly the physical radial measure, so Fourier–Bessel orthogonality is orthogonality in the geometry of the membrane. *Verified: Hermite constant $2^n n!\sqrt{\pi}$ and Bessel normalisation $\tfrac{a^2}{2}J_{\nu+1}(k_na)^2$ are standard (Arfken §14, §18).*

Limiting cases

Whether an SL problem has a clean discrete spectrum or something more exotic is decided by three switches: is the interval finite, does $p$ vanish at an endpoint, and does the domain run to infinity.

Case	Result	Meaning
Regular problem	$p,p',q,w$ continuous, $p>0$ , $w>0$ on finite $[a,b]$	Discrete, simple, real eigenvalues $\lambda_0<\lambda_1<\cdots\to+\infty$ ; eigenfunctions complete in $L^2_w[a,b]$ ; $y_n$ has exactly $n$ interior zeros (oscillation theorem).
Singular: $p$ vanishes at an end	$p(a)=0$ or $p(b)=0$	Boundary condition is replaced by boundedness/regularity; yields the classical polynomials — Legendre ( $p=1-x^2$ at $\pm1$ ), Laguerre ( $p=xe^{-x}$ at $0$ ), Bessel ( $p=x$ at $0$ ). Spectrum stays discrete.
Singular: infinite interval	$[0,\infty)$ or $(-\infty,\infty)$	Hermite ( $p=e^{-x^2}$ on $\mathbb{R}$ ) keeps a discrete spectrum because $w$ decays; a free particle ( $p=w=1$ on $\mathbb{R}$ ) has purely continuous spectrum with $\delta$ -normalised $e^{ikx}$ .
Limit-point vs limit-circle	number of $L^2_w$ solutions near a singular end $=1$ (LP) or $=2$ (LC)	Limit-point: no boundary condition needed there, self-adjointness is automatic. Limit-circle: one extra condition is required — the analogue of choosing a wall at a singular endpoint.

Common trapOrthogonality without the weight

The single most common error is writing $\int y_m y_n\,\mathrm{d}x=0$ and forgetting the weight $w(x)$ . Orthogonality is always in the weighted inner product $\int y_m y_n\,w\,\mathrm{d}x$ , because that is the only product in which $L$ is symmetric — the Lagrange identity carries the factor $w$ explicitly. Drop it and Hermite polynomials are simply not orthogonal (they overlap unless you weight by $e^{-x^2}$ ), Chebyshev polynomials need $(1-x^2)^{-1/2}$ , Laguerre need $e^{-x}$ . A companion trap: eigenfunctions are orthogonal only for distinct eigenvalues. A degenerate eigenspace (two eigenfunctions sharing one $\lambda$ ) is not automatically orthogonal and must be Gram–Schmidt-orthogonalised by hand.

Worked examples

1 · A variational bound on the lowest eigenvalue (fully worked)

Consider the regular SL problem $-y''=\lambda y$ on $[0,1]$ with $y(0)=y(1)=0$ — here $p=w=1$ , $q=0$ . Its exact spectrum is $\lambda_n=(n+1)^2\pi^2$ , so $\lambda_0=\pi^2$ . Use the Rayleigh quotient with the crude trial function $\phi(x)=x(1-x)$ to estimate $\lambda_0$ from above, and check how good the estimate is.

Step 1. Write down the Rayleigh quotient. For this problem it is $R[\phi]=\dfrac{\int_0^1 (\phi')^2\,\mathrm{d}x}{\int_0^1 \phi^2\,\mathrm{d}x}$ , and the min-max principle guarantees $\lambda_0\le R[\phi]$ for any admissible trial function satisfying the boundary conditions. Our $\phi=x(1-x)$ vanishes at both ends, so it qualifies.

Step 2. Do the two integrals. With $\phi'=1-2x$ , the numerator is $\int_0^1(1-2x)^2\,\mathrm{d}x=\int_0^1(1-4x+4x^2)\,\mathrm{d}x=1-2+\tfrac{4}{3}=\tfrac{1}{3}$ . The denominator is $\int_0^1 x^2(1-x)^2\,\mathrm{d}x=\int_0^1(x^2-2x^3+x^4)\,\mathrm{d}x=\tfrac{1}{3}-\tfrac{1}{2}+\tfrac{1}{5}=\tfrac{1}{30}$ .

Step 3. Form the ratio and compare. $R[\phi]=\dfrac{1/3}{1/30}=\mathbf{10}$ , while the exact value is $\lambda_0=\pi^2=9.8696$ . The estimate sits 1.32 % above the truth and, as the theorem promises, never below it — a one-line parabola already pins the ground eigenvalue to better than 2 %. R = 10 vs λ₀ = π² = 9.8696, a 1.3 % overestimate, and rigorously an upper bound. This is the engine behind the Ritz method: enrich the trial function and the bound tightens monotonically.

2 · The circular drumhead and inharmonic overtones

A circular membrane of radius $a$ clamped at its rim vibrates according to the Bessel SL problem of Derivation 3 (with $\nu=0$ for the axially symmetric modes). The mode frequencies are $f_{\nu n}=\dfrac{c\,\alpha_{\nu,n}}{2\pi a}$ , where $c$ is the wave speed and $\alpha_{\nu,n}$ is the $n$ -th zero of $J_\nu$ . Take $a=0.30\ \mathrm{m}$ and $c=100\ \mathrm{m/s}$ . Find the fundamental and the first overtone, and explain why a drum has no clear pitch.

Step 1. Get the Bessel zeros. The eigenvalues are quantised by clamping the rim, $J_\nu(k a)=0$ , so $k_{\nu n}a=\alpha_{\nu,n}$ . The relevant zeros are $\alpha_{0,1}=2.40483$ , $\alpha_{1,1}=3.83171$ , and $\alpha_{0,2}=5.52008$ — these are the eigenvalues $\lambda=k^2$ of the singular SL problem, with regularity at $r=0$ selecting $J_0$ over the divergent $Y_0$ .

Step 2. Convert to frequencies. The fundamental is the $(0,1)$ mode: $f_{01}=\dfrac{c\,\alpha_{0,1}}{2\pi a}=\dfrac{100\times2.40483}{2\pi\times0.30}=\dfrac{240.48}{1.8850}=\mathbf{127.6\ Hz}$ . The next mode is $(1,1)$ : $f_{11}=\dfrac{100\times3.83171}{1.8850}=\mathbf{203.3\ Hz}$ .

Step 3. Take the ratios. $f_{11}/f_{01}=\alpha_{1,1}/\alpha_{0,1}=3.83171/2.40483=\mathbf{1.593}$ , and the first symmetric overtone gives $\alpha_{0,2}/\alpha_{0,1}=5.52008/2.40483=\mathbf{2.295}$ . Neither ratio is a simple fraction like $2{:}1$ or $3{:}2$ , so the overtones are inharmonic: they do not form a harmonic series the way a string's do. The drum's overtones sit at 1.593, 2.295, … times the fundamental — irrational ratios of Bessel zeros — which is exactly why a struck membrane sounds like a "thud" with no definite pitch. Timpani recover a pitch only by loading the membrane and exciting mainly the $(1,1)$ -type modes.

Physics of the portrait

Sturm–Liouville theory is one of the rare corners of mathematical physics named for a genuine collaboration. Two friends in 1830s Paris — a specialist in differential equations and a founder of one of the great mathematical journals — proved together that the eigenfunctions of these operators actually converge to represent arbitrary functions, decades before Hilbert-space language existed to explain why.

Origins · 1 minJoseph Liouville and how this was discovered

Joseph Liouville

1809–1882, French mathematician and founder of the Journal de Mathématiques Pures et Appliquées · the convergence theory that turned Sturm's oscillation results into a working expansion method

Joseph Liouville was one of the most prolific and influential French mathematicians of the nineteenth century, a professor at the Collège de France and the École Polytechnique whose interests ranged from number theory (he was the first to prove that transcendental numbers exist, constructing Liouville's number in 1844) to celestial mechanics, complex analysis (Liouville's theorem — a bounded entire function is constant), and the differential geometry of surfaces. In 1836 he founded the Journal de Mathématiques Pures et Appliquées, still running today as "Liouville's Journal," which became the principal venue for French analysis. His deepest collaboration was with his close friend Charles-François Sturm. Sturm had published, in 1836, a landmark memoir on the qualitative behaviour of solutions of second-order linear ODEs — the comparison, separation, and oscillation theorems that count nodes without solving anything. Over 1836–1837 Sturm and Liouville together attacked the harder analytic question: do the eigenfunctions of such a problem actually form a complete set, so that an arbitrary function can be expanded in them the way Fourier expanded in sines? They proved convergence of these eigenfunction expansions under general conditions, establishing what we now call Sturm–Liouville theory. It was a startlingly modern achievement: a rigorous spectral theorem for a differential operator, complete with orthogonality and completeness, produced two generations before von Neumann and Stone gave self-adjoint operators their definitive Hilbert-space foundation. Every separation-of-variables calculation in physics — every spherical harmonic, every Bessel-mode drumhead, every quantum eigenstate expansion — is a descendant of what Sturm and Liouville proved in those two years.

History

1829

— Sturm begins his investigations into the oscillation of solutions of differential equations

1836

— Sturm publishes his memoir on second-order linear ODEs (comparison, separation, oscillation theorems)

1837

— Sturm and Liouville prove convergence of eigenfunction expansions, founding the theory

1910

— Weyl develops the limit-point/limit-circle classification of singular endpoints

1929

— von Neumann and Stone place self-adjointness on a rigorous Hilbert-space footing

Going further

Texts. Stone & Goldbart, Mathematics for Physics, ch. 8 (Linear Differential Operators) and ch. 9 (Sturm–Liouville) · Arfken, Weber & Harris, Mathematical Methods for Physicists, 7th ed., ch. 8 and ch. 18 · Hassani, Mathematical Physics, ch. 19 (Sturm–Liouville systems) · Byron & Fuller, Mathematics of Classical and Quantum Physics, ch. 5–6 · Riley, Hobson & Bence, Mathematical Methods for Physics and Engineering, 3rd ed., ch. 17 · Courant & Hilbert, Methods of Mathematical Physics, vol. I, ch. V–VI · Coddington & Levinson, Theory of Ordinary Differential Equations, ch. 7–9.

Provenance. The operator framework and the Lagrange-identity proof of reality/orthogonality follow Stone & Goldbart ch. 8–9 and Arfken–Weber–Harris ch. 8. A verifier should confirm: (1) the Lagrange identity (3) by two integrations by parts; (2) the integrating-factor weight $w=P^{-1}\exp\int(Q/P)\,\mathrm{d}x$ reproduces $p=Pw$ and $p'=Qw$ ; (3) the standard normalisations — Legendre $2/(2\ell+1)$ , Hermite $2^n n!\sqrt{\pi}$ , Bessel $\tfrac{a^2}{2}J_{\nu+1}(k_na)^2$ — from Rodrigues formulas and zero conditions; (4) the worked numbers: Rayleigh quotient $R=10$ against $\pi^2=9.8696$ (1.32 % high), and Bessel-zero ratios $\alpha_{1,1}/\alpha_{0,1}=1.593$ , $\alpha_{0,2}/\alpha_{0,1}=2.295$ from $\alpha_{0,1}=2.40483$ , $\alpha_{1,1}=3.83171$ , $\alpha_{0,2}=5.52008$ . No nonstandard constants are introduced.

Prerequisites. Orthogonal Polynomials: General Theory · Hilbert spaces and self-adjoint operators · Second-order linear ODEs. Leads to. Bessel Functions of the First and Second Kind · Legendre Polynomials · Recurrence Relations for Special Functions · Generating Functions for Special Functions.

Ask the tutor: Ask me to take any second-order ODE you name — Chebyshev, Laguerre, associated Legendre, or your own — and reduce it to Sturm–Liouville form, expose its weight, and read off its orthogonality and spectrum step by step.

Why does the weight w(x) appear at all — why not just use ordinary orthogonality?

Because the operator, not the physicist, decides the inner product. The Sturm–Liouville operator $L=-\frac{1}{w}[(py')'+qy]$ is symmetric only in the weighted product $\langle f,g\rangle_w=\int f g\,w\,\mathrm{d}x$ ; that $w$ is baked into $L$ by the $1/w$ out front, and it reappears in the Lagrange identity on the left-hand side. If you insist on ordinary orthogonality $\int f g\,\mathrm{d}x$ you are asking a different, non-symmetric operator to have real eigenvalues and orthogonal eigenfunctions, and it simply will not. The clean example is Hermite: $H_1=2x$ and $H_3=8x^3-12x$ have a nonzero plain overlap, but $\int H_1 H_3\,e^{-x^2}\,\mathrm{d}x=0$ . The weight is the price of self-adjointness.

What is the difference between a regular and a singular Sturm–Liouville problem?

A regular problem lives on a finite interval $[a,b]$ with $p,p',q,w$ continuous and $p,w>0$ everywhere, including the endpoints; you impose honest boundary conditions (Dirichlet, Neumann, or mixed) to kill the boundary term. A singular problem breaks one of those conditions: $p$ vanishes at an endpoint (Legendre at $\pm1$ , Bessel at $0$ ), the interval is infinite (Hermite on $\mathbb{R}$ ), or a coefficient blows up. When $p$ vanishes the boundary term dies on its own, so you replace the boundary condition with a regularity/boundedness requirement — and it is precisely this that selects the physically sensible solution ( $P_\ell$ not $Q_\ell$ , $J_\nu$ not $Y_\nu$ ) and quantises the spectrum. The classical special functions are, almost by definition, the eigenfunctions of singular problems.

How can the same differential equation give different special functions?

By changing the interval and boundary conditions — that is, by changing the domain of the operator rather than the operator itself. The equation $y''+\lambda y=0$ on $[0,\pi]$ with $y=0$ at both ends gives $\sin(nx)$ ; the same equation with periodic conditions on $[0,2\pi]$ gives the full Fourier basis; on $\mathbb{R}$ it gives the continuous spectrum $e^{ikx}$ . The differential expression is identical, but the self-adjoint operator is not defined until you specify its domain, and different admissible domains produce different complete bases. This is why "the SL equation" is a template, not a single problem, and why the same second-order form hosts sines, cosines, and exponentials depending on where you put the walls.

What does the oscillation theorem actually buy you?

It tells you the qualitative shape of every eigenfunction before you compute a single one. For a regular problem the eigenvalues order themselves $\lambda_0<\lambda_1<\cdots$ , and the $n$ -th eigenfunction has exactly $n$ zeros strictly inside the interval — the ground state never crosses zero, the first excited state crosses once, and so on. In quantum mechanics that is the node-counting rule for bound states: the number of nodes labels the energy level. Numerically it is priceless for root-finding, because it guarantees exactly how many eigenvalues lie below any trial value (a fact that underlies shooting and Sturm-sequence methods). And the zeros themselves are useful objects: the zeros of $P_\ell$ are the nodes of Gauss–Legendre quadrature.

What are the limit-point and limit-circle cases, in plain terms?

They classify how a singular endpoint behaves. Near such an endpoint the SL equation has two independent solutions; the question is how many of them are square-integrable (in the weight) close to that point. If only one is, the endpoint is limit-point: the non-integrable solution excludes itself, and you need no boundary condition there — self-adjointness is automatic. If both are square-integrable, the endpoint is limit-circle: you must supply one extra condition to pin down which combination is allowed, exactly as if you were choosing a wall at a place with no physical wall. Weyl proved this dichotomy in 1910. Physically it is the mathematics that decides whether a singular potential (say $1/r^2$ near the origin) needs an extra boundary condition to define a legitimate quantum Hamiltonian.

Is the Schrödinger equation really just a Sturm–Liouville problem?

Yes, and recognising it explains a lot. The time-independent Schrödinger equation $-\psi''+V\psi=E\psi$ (in one dimension, units $\hbar=2m=1$ ) is exactly the SL equation with $p=1$ , $q=-V$ , $w=1$ , and eigenvalue $\lambda=E$ . Self-adjointness of the SL operator is the statement that the Hamiltonian is an observable with a real spectrum; bound states are the discrete SL eigenvalues (a confining potential or a decaying weight keeps them in $L^2$ ), and scattering states are the continuous spectrum you get on an infinite interval. Even the subtle business of whether a singular potential needs an extra boundary condition at the origin is just the limit-point/limit-circle question. So the whole apparatus of orthogonality, completeness, and expansion that physicists use for wavefunctions is Sturm–Liouville theory wearing a quantum hat.