viva_math/special
Special mathematical functions.
Implementations focus on functions that complement gleam_community_maths
and are needed by viva_math/distributions for Gamma/Beta/Dirichlet
distributions, by viva_math/free_energy for KL divergences of
exponential-family distributions, and by viva_math/entropy for the
differential entropy of those families.
Algorithms
gamma→ Lanczos approximation (g=7, 9 coefficients). Accurate to ~15 decimal digits across the entire positive reals.lgamma→ Logarithmic form of Lanczos; preferred numerically because the regular gamma overflows past Γ(171.6) ≈ 10³⁰⁸.digamma→ ψ(x) via the asymptotic series for large x with the recurrence ψ(x) = ψ(x+1) - 1/x for small x.beta→ exp(lgamma(x) + lgamma(y) - lgamma(x+y)).factorial→ Stirling-style via gamma for large n, exact for n ≤ 20.
References
- Lanczos (1964) “A precision approximation of the gamma function”
- Press, Teukolsky et al. (1992) “Numerical Recipes” §6.1
- GSL
gsl_sf_lngamma
Values
pub fn beta(x: Float, y: Float) -> FloatBeta function B(x, y) = Γ(x)·Γ(y)/Γ(x+y). Always computed via lgamma
to avoid overflow.
pub fn binomial(n: Int, k: Int) -> Result(Float, Nil)Binomial coefficient C(n, k) = n! / (k!·(n-k)!) via lgamma for stability.
pub fn digamma(x: Float) -> FloatDigamma function ψ(x). Defined for x > 0.
Combines the recurrence ψ(x) = ψ(x+1) - 1/x to push x ≥ 20, then evaluates the asymptotic series ψ(x) ≈ ln x - 1/(2x) - Σ B₂ₖ / (2k·x^(2k)).
pub fn factorial(n: Int) -> Result(Float, Nil)Factorial n! for non-negative integer n.
Exact for n ≤ 20 (fits in i64). For larger n falls back to Γ(n+1) via Lanczos, with attendant ~15-digit relative accuracy.
pub fn gamma(x: Float) -> FloatGamma function Γ(x). Overflows past x ≈ 171.6 — prefer lgamma.