Legendre Polynomials: the Sphere's Natural Basis

Why a sphere gives you Pₙ instead of sin/cos

Separate variables in Laplace's equation ∇²V=0 in spherical coordinates and the polar-angle part comes out as (1−x²)y″ − 2xy′ + n(n+1)y = 0 with x=cos θ — Legendre's equation. Demanding a solution that stays finite at the poles (x=±1) forces n to be a non-negative integer, and the resulting solutions are polynomials: P₀=1, P₁=x, P₂=½(3x²−1), and so on by the recurrence (n+1)P(n+1) = (2n+1)xPn − nP(n−1). They're the angular building blocks behind any gravitational, electric, or magnetic field expanded around a sphere — Earth's gravity field, an atomic orbital, antenna radiation patterns.

The same roots show up in numerical integration

Pull up the roots of P_n for any n and you've just found the sample points for n-point Gaussian quadrature — the most accurate possible way to numerically integrate a polynomial up to degree 2n−1 using only n evaluations. It's not a coincidence: Gauss derived the quadrature rule by requiring it to integrate Legendre polynomials exactly, which forces the nodes to be exactly their roots. So the same special function shows up whether you're expanding a gravitational field or just trying to numerically integrate ∫₋₁¹ f(x) dx as efficiently as possible.

A short history

Adrien-Marie Legendre introduced these polynomials in 1782 while computing the gravitational attraction of spheroids (planets aren't perfect spheres, and Legendre needed a clean way to expand the correction terms). The connection to quadrature came almost a century later, when Carl Friedrich Gauss showed in 1814 that choosing Legendre's roots as integration nodes beats equally-spaced points by a wide margin — turning an obscure planetary-shape tool into one of the most-used algorithms in numerical computing.