mathhook_core/functions/special/
beta.rs

1use super::{gamma, lanczos_gamma};
2use crate::{Expression, Number};
3/// Numerically evaluates the beta function B(a, b) = Γ(a)·Γ(b)/Γ(a+b)
4///
5/// Uses Lanczos gamma approximation for high accuracy.
6///
7/// # Arguments
8///
9/// * `a` - First parameter
10/// * `b` - Second parameter
11///
12/// # Examples
13///
14/// ```rust
15/// use mathhook_core::functions::special::beta::beta_numerical;
16///
17/// let result = beta_numerical(2.0, 3.0);
18/// assert!((result - 1.0/12.0).abs() < 1e-14);
19/// ```
20pub fn beta_numerical(a: f64, b: f64) -> f64 {
21    if a.is_nan() || b.is_nan() || a.is_infinite() || b.is_infinite() {
22        return f64::NAN;
23    }
24    let gamma_a = lanczos_gamma(a);
25    let gamma_b = lanczos_gamma(b);
26    let gamma_ab = lanczos_gamma(a + b);
27    (gamma_a * gamma_b) / gamma_ab
28}
29
30/// Beta function B(a, b)
31///
32/// The Beta function is defined as:
33/// B(a, b) = Γ(a)·Γ(b) / Γ(a+b)
34///
35/// # Mathematical Properties
36///
37/// - B(a, b) = B(b, a) (symmetric)
38/// - B(a, b) = ∫₀¹ t^(a-1)·(1-t)^(b-1) dt
39///
40/// # Numerical Evaluation
41///
42/// Float inputs are evaluated numerically using Lanczos gamma approximation.
43/// Mixed Float/Integer inputs are converted to numerical evaluation.
44///
45/// # Arguments
46///
47/// * `a` - First parameter
48/// * `b` - Second parameter
49///
50/// # Examples
51///
52/// ```rust
53/// use mathhook_core::{Expression, Number};
54/// use mathhook_core::functions::special::beta;
55///
56/// let a = Expression::Number(Number::Integer(2));
57/// let b = Expression::Number(Number::Integer(3));
58/// let result = beta(&a, &b);
59///
60/// let a_float = Expression::Number(Number::Float(2.5));
61/// let b_float = Expression::Number(Number::Float(3.7));
62/// let result_num = beta(&a_float, &b_float);
63/// ```
64pub fn beta(a: &Expression, b: &Expression) -> Expression {
65    match (a, b) {
66        (Expression::Number(Number::Float(x)), Expression::Number(Number::Float(y))) => {
67            let result = beta_numerical(*x, *y);
68            Expression::Number(Number::Float(result))
69        }
70        (Expression::Number(Number::Float(x)), Expression::Number(Number::Integer(n))) => {
71            let result = beta_numerical(*x, *n as f64);
72            Expression::Number(Number::Float(result))
73        }
74        (Expression::Number(Number::Integer(n)), Expression::Number(Number::Float(y))) => {
75            let result = beta_numerical(*n as f64, *y);
76            Expression::Number(Number::Float(result))
77        }
78        _ => {
79            let gamma_a = gamma(a);
80            let gamma_b = gamma(b);
81            let sum = Expression::add(vec![a.clone(), b.clone()]);
82            let gamma_sum = gamma(&sum);
83            Expression::mul(vec![
84                gamma_a,
85                gamma_b,
86                Expression::pow(gamma_sum, Expression::Number(Number::Integer(-1))),
87            ])
88        }
89    }
90}
91
92#[cfg(test)]
93mod tests {
94    use super::*;
95
96    #[test]
97    fn test_beta_mixed_evaluation() {
98        let a = Expression::Number(Number::Float(2.5));
99        let b = Expression::Number(Number::Integer(3));
100        let result = beta(&a, &b);
101        match result {
102            Expression::Number(Number::Float(_)) => {}
103            _ => panic!("Beta with mixed inputs should return numerical result"),
104        }
105    }
106
107    #[test]
108    fn test_beta_symmetry() {
109        let result_ab = beta_numerical(2.5, 3.7);
110        let result_ba = beta_numerical(3.7, 2.5);
111        assert!(
112            (result_ab - result_ba).abs() < 1e-14,
113            "Beta symmetry: B(a,b) = B(b,a)"
114        );
115    }
116
117    #[test]
118    fn test_beta_numerical_evaluation() {
119        let result = beta_numerical(2.0, 3.0);
120        assert!((result - 1.0 / 12.0).abs() < 1e-14, "B(2,3) = 1/12");
121        let result_2_5 = beta_numerical(2.0, 5.0);
122        assert!((result_2_5 - 1.0 / 30.0).abs() < 1e-14, "B(2,5) = 1/30");
123    }
124
125    #[test]
126    fn test_beta_float_evaluation() {
127        let a = Expression::Number(Number::Float(2.5));
128        let b = Expression::Number(Number::Float(3.7));
129        let result = beta(&a, &b);
130        match result {
131            Expression::Number(Number::Float(_)) => {}
132            _ => panic!("Beta with float inputs should return numerical result"),
133        }
134    }
135}