Orthogonal polynomials and angular momentum

Legendre polynomial P_n(x) (DLMF 18.3; P_n(1) = 1).

With an integer degree n (a Python int or an integer Expr) and expand=True (the default) the explicit polynomial is returned, with exact rational coefficients; P_{-n-1} = P_n. A symbolic or non-integer degree, or expand=False, gives the unevaluated node legendre_p(n, x), which differentiates in x (d/dx P_n = C_{n-1}^{(3/2)}), evaluates numerically (f64 and rigorous balls, any real degree) and simplifies at special values (P_n(1) = 1, parity for integer n).

Raises SpecialFunctionError (E-SPECFN-004) for an explicit expansion past degree 1000 — keep it symbolic with expand=False.

>>> import alkahest as ak
>>> p = ak.ExprPool(); x = p.symbol("x")
>>> ak.legendre_p(2, x)            # 3/2*x^2 - 1/2

Associated Legendre (Ferrers) function P_n^m(x), Condon–Shortley phase included: P_n^m(x) = (-1)^m (1-x^2)^{m/2} d^m P_n/dx^m for m >= 0 and P_n^{-m} = (-1)^m (n-m)!/(n+m)! P_n^m (DLMF 14.6.1, 14.9.3; SymPy's assoc_legendre). |m| > n gives 0.

Integer n and m expand ((1-x^2)^{m/2} times a polynomial); otherwise the node assoc_legendre_p(n, m, x). For odd m the function is real only on [-1, 1]; numeric evaluation outside it declines rather than return a real part.

Physicists' Hermite polynomial H_n(x) (weight exp(-x^2), H_n = 2^n x^n + ...; SymPy's hermite). Integer n >= 0 expands; d/dx H_n = 2n H_{n-1}.

Probabilists' Hermite polynomial He_n(x) (monic, weight exp(-x^2/2); SymPy's hermite_prob). He_n(x) = 2^{-n/2} H_n(x/sqrt(2))`; `d/dx He_n = n He_{n-1}.

Laguerre polynomial L_n(x) = L_n^{(0)}(x) (weight exp(-x) on [0, oo); SymPy's laguerre). d/dx L_n = -L_{n-1}^{(1)}.

Generalised (associated) Laguerre polynomial L_n^{(alpha)}(x) (SymPy's assoc_laguerre). alpha may be symbolic: the expansion is then polynomial in it. d/dx L_n^{(a)} = -L_{n-1}^{(a+1)}. Hydrogen radial functions use L_{n-l-1}^{(2l+1)}.

Chebyshev polynomial of the first kind T_n(x) (T_n(cos t) = cos(n t)`; `T_{-n} = T_n`). `d/dx T_n = n U_{n-1}.

Chebyshev polynomial of the second kind U_n(x) (U_n(cos t) = sin((n+1)t)/sin t`; `U_{-1} = 0`, `U_{-n} = -U_{n-2}).

Gegenbauer (ultraspherical) polynomial C_n^{(a)}(x) (SymPy's gegenbauer; C_n^{(0)} = 0 for n >= 1). d/dx C_n^{(a)} = 2a C_{n-1}^{(a+1)}.

Jacobi polynomial P_n^{(a,b)}(x) (SymPy's jacobi). d/dx P_n^{(a,b)} = (n+a+b+1)/2 P_{n-1}^{(a+1,b+1)}; parity swaps a and b.

Rewrite every orthogonal-polynomial node in expr whose degree (and Legendre order) is now an integer literal into its explicit polynomial — e.g. after subs(legendre_p(n, x), {n: 5}). Nodes with a symbolic degree are left alone.

Complex spherical harmonic Y_l^m(theta, phi) as an explicit expression.

Y_l^m = sqrt((2l+1)/(4 pi) (l-m)!/(l+m)!) P_l^m(cos theta) exp(i m phi) with the Condon–Shortley phase inside P_l^m — SymPy's Ynm and Mathematica's SphericalHarmonicY; orthonormal on the unit sphere, Y_l^{-m} = (-1)^m conj(Y_l^m). (1 - cos^2 theta)^{|m|/2} is written sin(theta)^|m|, exact for the polar range 0 <= theta <= pi. |m| > l gives 0. The imaginary unit is the pool's I.

Raises SpecialFunctionError: E-SPECFN-003 for non-integer l or m, E-SPECFN-002 for l < 0, E-SPECFN-004 for l > 1000.

Real spherical harmonic Y_{l,m}(theta, phi) as an explicit expression.

The tesseral (real) basis, Wikipedia's convention: for m > 0 sqrt(2)*N*P_l^m(cos theta)*cos(m phi), for m < 0 sqrt(2)*N*P_l^{|m|}(cos theta)*sin(|m| phi), and Y_l^0 for m = 0, where N = sqrt((2l+1)/(4 pi) (l-|m|)!/(l+|m|)!) and P_l^m is the Legendre function without the Condon-Shortley phase (the phase of the complex Y_l^m cancels the (-1)^m in the definition sqrt(2) (-1)^m Re/Im Y_l^{|m|}). So Y_{1,1} = sqrt(3/(4 pi)) x/r, Y_{1,-1} = sqrt(3/(4 pi)) y/r, Y_{1,0} = sqrt(3/(4 pi)) z/r. The functions are real and orthonormal on the unit sphere. sin(theta)**|m| stands for (1 - cos(theta)**2)**(|m|/2) (exact for 0 <= theta <= pi). |m| > l gives 0.

Raises SpecialFunctionError: E-SPECFN-003 for non-integer l or m, E-SPECFN-002 for l < 0, E-SPECFN-004 for l > 1000.

>>> import alkahest as ak
>>> p = ak.ExprPool(); t = p.symbol("t"); f = p.symbol("f")
>>> ak.real_spherical_harmonic_y(1, 1, t, f)   # sqrt(3/(4 pi)) sin(t) cos(f)

Wigner 3j symbol (j1 j2 j3; m1 m2 m3), exactly.

Arguments are integers or half-integers: int, Fraction, a float that is exactly k/2 (0.5), or a numeric Expr. The result is an Expr of the form q or q*r**(1/2) (q rational, r square-free), computed by Racah's formula over big integers — the same normalisation and phase as SymPy's wigner_3j and Mathematica's ThreeJSymbol. Selection-rule violations (triangle, |m| > j, m1+m2+m3 != 0, j-m not an integer) give 0.

The result lives in the pool of any Expr argument, else pool, else a module-wide default pool.

Raises SpecialFunctionError: E-SPECFN-001 for an argument that is not an integer or half-integer, E-SPECFN-002 for a negative j, E-SPECFN-004 past 2j = 20000.

>>> import alkahest as ak
>>> ak.wigner_3j(1, 1, 0, 0, 0, 0)      # -1/3*3^(1/2)

Wigner 6j symbol {j1 j2 j3; j4 j5 j6}, exactly (Racah's formula; SymPy's wigner_6j). A triad (j1 j2 j3), (j1 j5 j6), (j4 j2 j6), (j4 j5 j3) that fails the triangle condition gives 0. Inputs, pool and refusals as wigner_3j().

Clebsch–Gordan coefficient <j1 m1; j2 m2 | j3 m3> (Condon–Shortley phase), exactly. Argument order is SymPy's: clebsch_gordan(j1, j2, j3, m1, m2, m3).

<j1 m1 j2 m2|J M> = (-1)^(j1-j2+M) sqrt(2J+1) (j1 j2 J; m1 m2 -M). Inputs, pool and refusals as wigner_3j().

>>> import alkahest as ak
>>> from fractions import Fraction as F
>>> ak.clebsch_gordan(F(1, 2), F(1, 2), 1, F(1, 2), F(-1, 2), 0)   # 1/2*2^(1/2)

:no-index:

E-SPECFN-001 … E-SPECFN-005; see errors.