Orthogonal Polynomials Example

Demonstrates three families of classical orthogonal polynomials using specialfunctions.legendre_poly(), specialfunctions.chebyshev_t(), and specialfunctions.hermite_poly(). Degrees 0 through 4 are shown for each family on their natural domains.

Example Code

"""Example: IMSL Special Functions — Orthogonal Polynomials."""
from __future__ import annotations

import numpy as np
import matplotlib
matplotlib.use("Agg")
import matplotlib.pyplot as plt
from specialfunctions import legendre_poly, chebyshev_t, hermite_poly

# ---------------------------------------------------------------------------
# Known values
# ---------------------------------------------------------------------------
print(f"P3(0.5)  = {legendre_poly(3, 0.5):.10f}  "
      f"(expected: {(5*0.5**3 - 3*0.5)/2:.10f})")
print(f"T3(0.5)  = {chebyshev_t(3, 0.5):.10f}  "
      f"(expected: {4*0.5**3 - 3*0.5:.10f}  = 4·(0.5)³ − 3·(0.5))")
print(f"He3(1.0) = {hermite_poly(3, 1.0):.10f}  "
      f"(expected: {1.0**3 - 3*1.0:.10f}  = x³ − 3x at x=1)")

# ---------------------------------------------------------------------------
# Plots — 3 panels, one per polynomial family
# ---------------------------------------------------------------------------
fig, axes = plt.subplots(1, 3, figsize=(16, 5))
fig.suptitle("Orthogonal Polynomials", fontsize=14)

colors = ["tab:blue", "tab:orange", "tab:green", "tab:red", "tab:purple"]

# --- Legendre P0 .. P4 on [-1, 1] ---
x_leg = np.linspace(-1.0, 1.0, 400)
for n in range(5):
    y = legendre_poly(n, x_leg)
    axes[0].plot(x_leg, y, color=colors[n], lw=1.5, label=f"$P_{n}$")
axes[0].axhline(0, color="k", lw=0.5, ls="--")
axes[0].set_title("Legendre Polynomials $P_n(x)$")
axes[0].set_xlabel("x")
axes[0].set_xlim(-1, 1)
axes[0].set_ylim(-1.1, 1.1)
axes[0].legend(fontsize=8)
axes[0].grid(True, alpha=0.3)

# --- Chebyshev T0 .. T4 on [-1, 1] ---
x_cheb = np.linspace(-1.0, 1.0, 400)
for n in range(5):
    y = chebyshev_t(n, x_cheb)
    axes[1].plot(x_cheb, y, color=colors[n], lw=1.5, label=f"$T_{n}$")
axes[1].axhline(0, color="k", lw=0.5, ls="--")
axes[1].set_title("Chebyshev Polynomials $T_n(x)$")
axes[1].set_xlabel("x")
axes[1].set_xlim(-1, 1)
axes[1].set_ylim(-1.1, 1.1)
axes[1].legend(fontsize=8)
axes[1].grid(True, alpha=0.3)

# --- Hermite He0 .. He4 on [-4, 4] ---
x_herm = np.linspace(-4.0, 4.0, 400)
for n in range(5):
    y = hermite_poly(n, x_herm)
    y_clipped = np.where(np.abs(y) > 20, np.nan, y)
    axes[2].plot(x_herm, y_clipped, color=colors[n], lw=1.5, label=f"$He_{n}$")
axes[2].axhline(0, color="k", lw=0.5, ls="--")
axes[2].set_title(r"Hermite Polynomials $He_n(x)$")
axes[2].set_xlabel("x")
axes[2].set_xlim(-4, 4)
axes[2].set_ylim(-20, 20)
axes[2].legend(fontsize=8)
axes[2].grid(True, alpha=0.3)

plt.tight_layout()
plt.savefig("test_output/example_imsl_orthogonal_polynomials.svg", bbox_inches="tight")
print("\nSaved: test_output/example_imsl_orthogonal_polynomials.svg")

Plot Output

Left: Legendre polynomials \(P_n(x)\) on \([-1,1]\). Centre: Chebyshev polynomials of the first kind \(T_n(x)\) on \([-1,1]\), bounded in \([-1,1]\) for all \(n\). Right: probabilist’s Hermite polynomials \(He_n(x)\) on \([-4,4]\) (clipped at \(\pm 20\) for visibility).

Console Output

P3(0.5)  = -0.4375000000  (expected: -0.4375000000)
T3(0.5)  = -1.0000000000  (expected: -1.0000000000  = 4·(0.5)³ − 3·(0.5))
He3(1.0) = -2.0000000000  (expected: -2.0000000000  = x³ − 3x at x=1)

Saved: test_output/example_imsl_orthogonal_polynomials.svg