University Quantum Mechanics II · Deriving the Spherical Harmonics · 1.4
Generating function and multipole expansion of 1/|r − r′|
Helium's repulsion, the screening of an outer electron, the quadrupole moment of a molecule: each needs 1/|r − r′| integrated over the positions of the charges, and the square root will not factor. This lesson rewrites it as a series whose terms separate the two radii from the angle between them, and reads the Legendre polynomials off the same root.
Build the model
Connect the measurement to the mechanism.
Take the larger radius out of the distance between two points: |r − r′| = r_>(1 − 2t cos γ + t²)¹⁄², with t = r_</r_> ≤ 1 and γ the angle between r and r′. So G(x, t) = (1 − 2xt + t²)(−1/2) holds the whole geometry of the Coulomb interaction, and it has a second life. Expanded in powers of t, its coefficients are polynomials in x that satisfy Legendre's equation with λ = l(l + 1) and equal 1 at x = 1, so they are the Pₗ(x) of the previous topic, derived here without Rodrigues' formula. One t-derivative of G gives Bonnet's recursion, the working tool for computing Pₗ and the reason cos θ links l only to l ± 1.
Together they give 1/|r − r′| = Σ (r_<ˡ/r_>(l+1)) Pₗ(cos γ). This is an exact identity for every r ≠ r′, not a far-field approximation. The radii and the angle sit in separate factors, so an angular integral can remove whole terms: against a spherical charge every l ≥ 1 term vanishes. The terms fall as tˡ, so convergence is geometric.
The model has three edges. As r′ → r convergence slows, and at r = r′ with γ = 0 the series diverges because the points coincide. When neither vector lies on the polar axis, Pₗ(cos γ) must still be split into spherical harmonics by the addition theorem of UQMII-01.09. And 1/|r − r′| is instantaneous electrostatics: its magnetic and retardation corrections enter at relative order α², a result quoted from relativistic quantum mechanics.
- Simple definition
- The Legendre generating function G(x, t) = (1 − 2xt + t²)(−1/2) is the function whose coefficient of tˡ is Pₗ(x); with x = cos γ and t = r_</r_> it equals r_> times 1/|r − r′|, which makes its series the multipole expansion of the Coulomb interaction.
- Example
- Charges at r = 0.50 nm and r′ = 0.10 nm on the same ray have t = 0.20 and γ = 0. The first four terms give (1/0.50 nm)(1 + 0.2 + 0.04 + 0.008) = 2.496 nm⁻¹ against the exact 1/0.40 nm = 2.500 nm⁻¹, short by t⁴ = 0.16%.
The binomial series of the root produces P₀ = 1, P₁ = x and P₂ = (3x² − 1)/2 with no Rodrigues' formula and no integral, and x = 1 gives Pₗ(1) = 1 for every l.
x = cos γ and t are dimensionless; converges for |t| < 1 when −1 ≤ x ≤ 1
One t-derivative of the generating function gives it. It computes Pₗ(x) in l steps, and it says x Pₗ contains only Pₗ₊₁ and Pₗ₋₁.
l = 1, 2, 3, …; x and Pₗ dimensionless; seeds P₀ = 1 and P₁ = x
Radii and angle separate term by term, so an integral over two positions becomes a sum of radial integrals times angular integrals.
r_< and r_> are the smaller and larger of r and r′, in m; γ is the angle between r and r′; every term is in m⁻¹
Relative to the exact value at γ = 0 the error is exactly t(L+1): 0.1% needs L = 2 at t = 0.1 but L = 9 at t = 0.5.
RL = exact 1/|r − r′| minus the terms l = 0 to L, in m⁻¹; t = r_</r_>; equality at γ = 0
∫ Pₗ(cos γ) dΩ′ = 4π for l = 0 and zero otherwise, so a spherical charge acts only through its l = 0 term: the shell theorem.
dΩ′ = sin θ′ dθ′ dφ′ over all directions of r′ at fixed r and r′; the result is in m⁻¹
Q₀ is the net charge and Q₁ the dipole moment. The first nonzero Qₗ sets the far-field fall-off, 1/r(l+1), and the angular pattern, Pₗ(cos θ).
charges qᵢ in C at heights zᵢ in m; r beyond every |zᵢ|; Qₗ in C mˡ; φ in V
Read the square root as a power series in t
Expand G(x, t) = (1 − 2xt + t²)(−1/2) with the binomial series (1 − u)(−1/2) = 1 + u/2 + 3u²/8 + 5u³/16 + …, where u = 2xt − t². Collecting powers of t, the t¹ coefficient is x, the t² coefficient is 3x²/2 − 1/2 = (3x² − 1)/2, and the t³ coefficient is (5x³ − 3x)/2. These are P₁, P₂ and P₃, and nothing was differentiated or integrated to get them. Each power of u costs at least one power of t and brings at most one power of x, so the tˡ coefficient is a polynomial of degree l. Two facts then hold for every l at once. At x = 1, G = 1/(1 − t) = 1 + t + t² + …, so Pₗ(1) = 1. G is unchanged by x → −x together with t → −t, so Pₗ(−x) = (−1)ˡ Pₗ(x), the parity the Yₗᵐ will inherit. The series converges for |t| < 1 because, with x = cos γ, the roots of 1 − 2xt + t² sit at t = cos γ ± i sin γ, on the unit circle.
Prove that the coefficients are the Legendre polynomials
Matching three coefficients is evidence, not proof. Differentiate G directly and a few lines give (1 − x²) ∂²G/∂x² − 2x ∂G/∂x + t ∂²(tG)/∂t² = 0 for all x and t. Put G = Σ Qₗ(x) tˡ. The last term becomes Σ l(l + 1) Qₗ tˡ, and since different powers of t are independent, each coefficient obeys (1 − x²)Qₗ″ − 2xQₗ′ + l(l + 1)Qₗ = 0: Legendre's equation with λ = l(l + 1). Qₗ is a polynomial, so it is the terminating solution of UQMII-01.02, fixed up to a constant, and Qₗ(1) = 1 fixes the constant. So Qₗ = Pₗ, and this lesson derives it rather than quoting it. An independent check reuses UQMII-01.03. An elementary integral of G(x, t) G(x, s) over −1 ≤ x ≤ 1 gives (ts)(−1/2) ln[(1 + √(ts))/(1 − √(ts))] = Σ 2(ts)ˡ/(2l + 1). It depends on t and s only through their product, so no tˡsᵏ with l ≠ k appears, which is orthogonality, and the tˡsˡ coefficient is the norm 2/(2l + 1).
Differentiate once in t to get the recursion
Since ∂G/∂t = (x − t)G/(1 − 2xt + t²), the generating function obeys (1 − 2xt + t²) ∂G/∂t = (x − t)G. Substitute the series and collect the coefficient of tˡ: (l + 1)Pₗ₊₁ − 2lxPₗ + (l − 1)Pₗ₋₁ = xPₗ − Pₗ₋₁, which rearranges to Bonnet's recursion, (2l + 1)xPₗ = (l + 1)Pₗ₊₁ + lPₗ₋₁. Seeded with P₀ = 1 and P₁ = x, it gives any Pₗ(x) in l steps without expanding a polynomial: at x = 0.5 it returns P₂ = −0.125, P₃ = −0.4375 and P₄ = −0.2891. Read the other way, it says that multiplying Pₗ by x produces only Pₗ₊₁ and Pₗ₋₁. With x = cos θ, a perturbation proportional to z = r cos θ can therefore connect l only to l ± 1, the seed of the dipole selection rule that UQMII-01.08 extends to the Yₗᵐ. The x-derivative of G gives a companion, P′ₗ₊₁ − P′ₗ₋₁ = (2l + 1)Pₗ, after one elimination.
Factor out the larger radius, never the smaller
The law of cosines gives |r − r′|² = r² + r′² − 2rr′ cos γ. Take out the larger radius: |r − r′|² = r_>²(1 − 2t cos γ + t²) with t = r_</r_>, so 1/|r − r′| = G(cos γ, t)/r_> = Σ (r_<ˡ/r_>(l+1)) Pₗ(cos γ). Taking out the smaller radius gives the same root with t > 1, outside the region of convergence, and a series that diverges. Check the result three ways. Dimensions: r_<ˡ/r_>(l+1) is a length⁻¹ for every l. Limit: as r′ → 0 only 1/r survives, and the l = 1 term r′ cos γ/r² is the dipole correction for a charge moved off the origin by r′. Symmetry: exchanging r and r′ leaves every term unchanged, as a distance must. The l = 0 term 1/r_> has a kink at r = r′, where its slope in r jumps from 0 to −1/r′². That jump is Gauss's law for a unit charge spread evenly over the sphere of radius r′.
Count the terms, and find where the series gives out
Since |Pₗ(x)| ≤ 1 on −1 ≤ x ≤ 1, a bound quoted here rather than derived, the terms after l = L add up to at most t(L+1)/[(1 − t)r_>]. At γ = 0 every Pₗ equals 1 and the bound is reached, so aligned points converge most slowly, and there the error relative to the exact 1/(r_> − r_<) is exactly t(L+1). Reaching 0.1% therefore takes L = 2 at t = 0.1, L = 9 at t = 0.5 and L = 65 at t = 0.9: geometric convergence, at a rate set by how nearly equal the radii are. At t = 1 it still converges for 0 < γ < π, slowly and only conditionally, to 1/(2 sin(γ/2)); it fails at γ = π, where the terms are ±1, and diverges at γ = 0, where the two points coincide. For two electrons in one atom r′ runs through r, so a repulsion integral splits its radial range at r′ = r and uses r_</r_> on each side. The identity is exact on both sides; what makes it practical is that the angular integrals usually stop the sum after a few terms.
Average over a sphere and only l = 0 survives
Average 1/|r − r′| over every direction of r′ at fixed r and r′. Put the polar axis along r, so γ = θ′ and ∫ Pₗ(cos θ′) dΩ′ = 2π ∫ Pₗ(x) dx from −1 to 1. That last integral is orthogonality against P₀ = 1: it is 2 for l = 0 and zero for every l ≥ 1. The directional average of 1/|r − r′| is therefore exactly 1/r_>, a derived result. For a spherically symmetric charge it is the shell theorem: charge inside radius r acts from the centre, and each outer shell adds a constant potential. It predicts a measured number. In helium's 1s3d state the outer electron lies almost entirely outside the 1s cloud, so it sees charge 2 − 1 = 1 and should be bound like hydrogen's n = 3, by 13.606 eV/9 = 1.512 eV. The tabulated binding of the 1s3d ¹D₂ level is 1.513 eV. Both are of order α²mc², and the same reduction to l = 0 is what evaluates the 1s² repulsion integral of helium's ground state in Unit 10.
Change one variable at a time
Make the relationship visible.
Hold L = 2 and slide t from 0.1 to 0.8: the solid curve peels away at γ = 0 first, where the error is t³ of the peak. At t = 0.8 raise L one step at a time and watch the bar shrink slowly. Then set γ = 90° and step L from 2 to 3: the open circle stays put, because every odd Pₗ vanishes there.
EXACT AT γ1.1547 × 1/r_>
WORST ERROR, AT γ = 00.2500 × 1/r_>
WORST ERROR / EXACT PEAK12.50 %
L NEEDED FOR 0.1%9
Live interpretationEXACT AT γ: 1.1547 × 1/r_>. WORST ERROR, AT γ = 0: 0.2500 × 1/r_>. WORST ERROR / EXACT PEAK: 12.50 %. L NEEDED FOR 0.1%: 9
Catch the common trap
Explain before calculating.
Two electrons sit at r₁ = 2.0 a₀ and r₂ = 5.0 a₀ from a nucleus, with 60° between their position vectors. Keeping the l = 0 and l = 1 terms of the multipole expansion, what is 1/|r₁ − r₂|?
Choose an answer to test the model.
Practice & worked examples
Reason from the model, then test the result.
EasyUsing only Bonnet's recursion seeded with P₀ = 1 and P₁ = x, find P₂, P₃ and P₄ at x = 0.5. Then sum Pₗ(0.5)tˡ from l = 0 to 4 at t = 0.3, and compare the result with the closed-form generating function and with the truncation bound.
- l = 1: 2P₂ = 3xP₁ − P₀ = 3(0.5)(0.5) − 1 = −0.25, so P₂ = −0.125, which the closed form (3x² − 1)/2 = −0.125 confirms.
- l = 2: 3P₃ = 5xP₂ − 2P₁ = 5(0.5)(−0.125) − 2(0.5) = −1.3125, so P₃ = −0.4375. l = 3: 4P₄ = 7xP₃ − 3P₂ = −1.53125 + 0.375 = −1.15625, so P₄ = −0.2890625, which is −37/128.
- Series at t = 0.3: 1 + 0.5(0.3) − 0.125(0.09) − 0.4375(0.027) − 0.2891(0.0081) = 1 + 0.15 − 0.01125 − 0.01181 − 0.00234 = 1.12460.
- Closed form: (1 − 2(0.5)(0.3) + 0.09)(−1/2) = 0.79(−1/2) = 1.12509. The series falls short by 4.9 × 10⁻⁴.
- Bound on everything omitted: t⁵/(1 − t) = 0.00243/0.7 = 3.5 × 10⁻³. The actual shortfall is seven times smaller because |Pₗ(0.5)| stays well below 1; the bound is reached only at x = 1, where every Pₗ equals 1.
AnswerP₂ = −0.125, P₃ = −0.4375, P₄ = −0.2891. The series to t⁴ gives 1.1246 against the exact 1.1251, short by 4.9 × 10⁻⁴ and inside the bound 3.5 × 10⁻³.
MediumA point-charge model of CO₂ puts +0.80e on the carbon at the origin and −0.40e on each oxygen, at z = ±0.116 nm. Find its first nonzero axial moment Qₗ, then the potential at r = 0.348 nm both on the axis and broadside (θ = 90°), from that term alone and exactly. Take e/4πε₀ = 1.440 V nm.
- Moments Qₗ = Σ qᵢ zᵢˡ: Q₀ = 0.80e − 2(0.40e) = 0. Q₁ = −0.40e(0.116 nm) − 0.40e(−0.116 nm) = 0, and every odd Qₗ vanishes because the oxygens form a ±z pair. Q₂ = −2(0.40e)(0.116 nm)² = −0.01076e nm² = −1.72 × 10⁻³⁹ C m².
- Leading term φ₂ = (1/4πε₀) Q₂ P₂(cos θ)/r³, with r³ = 0.04214 nm³. On the axis P₂(1) = 1: φ₂ = 1.440 V nm × (−0.01076 nm²)/0.04214 nm³ = −0.368 V.
- Exact on the axis: φ = 1.440 V nm × (0.80/0.348 − 0.40/0.232 − 0.40/0.464) nm⁻¹ = −0.414 V. The quadrupole term carries 8/9 of it: at r = 3d each even term is (d/r)² = 1/9 of the one before.
- Broadside, P₂(0) = −1/2 flips the sign: φ₂ = +0.184 V. Exactly, φ = 1.440 V nm × (0.80/0.348 − 0.80/0.3668) nm⁻¹ = +0.170 V, so the leading term is 8.3% high. Adding l = 4, with P₄(0) = 3/8, gives +0.169 V, within 0.8%.
- Checks: Q₂/(4πε₀r³) is C m² × V m C⁻¹ × m⁻³ = V. As the bond length d → 0 the potential vanishes, as three charges summing to zero at one point must. The measured CO₂ quadrupole moment, about −1.4 × 10⁻³⁹ C m² from tables, has the model's sign and is about a fifth smaller.
AnswerQ₀ = Q₁ = 0 and Q₂ = −1.72 × 10⁻³⁹ C m². At 0.348 nm: on the axis −0.368 V from l = 2 against −0.414 V exact; broadside +0.184 V against +0.170 V exact.
HardIn helium's 1s3d state the inner electron has the He⁺ ground-state density |ψ₁ₛ|² = (8/πa₀³)e(−4r/a₀). Treat that cloud as a fixed charge, use the multipole expansion to find the potential energy U(r) of the outer electron, and predict the 3d binding energy below He⁺(1s), including the first-order penetration correction. Compare with the tabulated 1.5133 eV of the 1s3d ¹D₂ level and say what the residual is.
- The cloud's density depends on r′ only, so the angular integral of every Pₗ(cos γ) with l ≥ 1 is zero and 1/|r − r′| averages to 1/r_>. With k = e²/4πε₀, U(r) = −2k/r + k Nᵢₙ(r)/r + k∫ |ψ₁ₛ|² 4πr′ dr′ over r′ > r, where Nᵢₙ(r) is the probability that the 1s electron lies inside radius r.
- The two radial integrals give U(r) = −k/r − k e(−4r/a₀)(1/r + 2/a₀). Limits check it: as r → ∞ the outer electron sees net charge 2 − 1 = 1; as r → 0, U → −2k/r + 2k/a₀, the full nucleus plus the 54.4 eV that the cloud contributes at its centre, which is k⟨1/r′⟩ for Z = 2.
- Relative to −k/r the correction is e(−4r/a₀)(1 + 2r/a₀): 5.5% at r = a₀, 0.17% at 2a₀, and about 10⁻¹⁷ at 10.5a₀, the mean radius of a hydrogenic 3d orbit. The 3d electron has only 5 × 10⁻⁴ of its probability inside 2a₀.
- Zeroth order, U = −k/r, so the level is hydrogen's n = 3 (quoted from QM I): 13.6057 eV/9 = 1.5117 eV, or 1.5115 eV with the ⁴He reduced mass, μ/m = 1 − 1.37 × 10⁻⁴. First-order penetration: ⟨ΔU⟩ = −k⟨e(−4r/a₀)(1/r + 2/a₀)⟩ over the 3d density ∝ r⁶e(−2r/3a₀); the two gamma integrals give −3.37 × 10⁻⁶ k/a₀ = −9.2 × 10⁻⁵ eV.
- Prediction: 1.5115 eV + 0.0001 eV = 1.5116 eV of binding. The table gives 1.5133 eV for 1s3d ¹D₂, 0.0017 eV or 0.11% more; the ³D level lies another 0.4 meV deeper, the singlet-triplet exchange splitting, which a frozen Hartree cloud also leaves out. Fine structure cannot close the gap: it is at most of relative order α² ≈ 5 × 10⁻⁵, under 0.1 meV here, twenty times too small.
- The residual is mainly core polarisation: the outer electron's field distorts the 1s cloud, a second-order effect of the l = 1 term that a frozen spherical cloud leaves out. A standard dipole-polarisation estimate, quoted here, gives about 2 meV, the right size.
AnswerU(r) = −k/r − k e(−4r/a₀)(1/r + 2/a₀). Predicted binding 1.5116 eV (1.5115 eV hydrogenic plus 0.09 meV of penetration) against 1.5133 eV measured for 1s3d ¹D₂; the 0.11% residual is mainly core polarisation.