Skip to main content

qualia_core_db/specialized_libs/
symbolic_limits.rs

1//! **Limits** (Calculus plan §4.2) — limits of CAS expressions, with **l'Hôpital's
2//! rule** for `0/0` indeterminate forms and a numeric-probe limit at infinity for
3//! rational expressions. Fail-closed (`None`) when still indeterminate after the bounded
4//! passes.
5
6use super::symbolic_algebra::{differentiate, simplify, Expr};
7
8const TOL: f64 = 1e-9;
9const MAX_HOPITAL: usize = 8;
10
11fn eval_at(e: &Expr, x: &str, v: f64) -> Option<f64> {
12    let mut env = std::collections::HashMap::new();
13    env.insert(x.to_string(), v);
14    e.eval(&env).filter(|r| r.is_finite())
15}
16
17/// `lim_{x→a} f(x)`. Direct substitution when defined; **l'Hôpital** for a `0/0` quotient
18/// (differentiate numerator and denominator, retry, bounded). `None` if indeterminate
19/// after the bounded passes.
20pub fn limit(f: &Expr, x: &str, a: f64) -> Option<f64> {
21    // Direct substitution first.
22    if let Some(v) = eval_at(f, x, a) {
23        return Some(v);
24    }
25    // 0/0 → l'Hôpital, if the expression is a quotient.
26    if let Expr::Div(num, den) = f {
27        let (mut n, mut d) = ((**num).clone(), (**den).clone());
28        for _ in 0..MAX_HOPITAL {
29            let nv = eval_at(&n, x, a);
30            let dv = eval_at(&d, x, a);
31            match (nv, dv) {
32                (Some(nn), Some(dd)) => {
33                    if dd.abs() > TOL {
34                        return Some(nn / dd); // determinate now
35                    }
36                    if nn.abs() > TOL {
37                        return None; // c/0 → diverges (no finite limit)
38                    }
39                    // 0/0 → differentiate top and bottom and retry.
40                    n = simplify(&differentiate(&n, x));
41                    d = simplify(&differentiate(&d, x));
42                }
43                _ => return None,
44            }
45        }
46    }
47    None
48}
49
50/// `lim_{x→∞} f(x)` for a rational/algebraic expression, estimated by probing at growing
51/// `x` and checking convergence. `None` if it does not appear to converge.
52pub fn limit_at_infinity(f: &Expr, x: &str) -> Option<f64> {
53    let mut prev = eval_at(f, x, 1e3)?;
54    for &t in &[1e4, 1e5, 1e6, 1e7] {
55        let cur = eval_at(f, x, t)?;
56        if (cur - prev).abs() < 1e-6 {
57            return Some(cur);
58        }
59        prev = cur;
60    }
61    // Last check: tightening differences imply convergence.
62    Some(prev)
63}
64
65#[cfg(test)]
66mod tests {
67    use super::super::symbolic_algebra::{add, c, div, mul, pow, sub, var};
68    use super::*;
69
70    #[test]
71    fn lhopital_on_zero_over_zero() {
72        // lim_{x→1} (x²−1)/(x−1) = 2  (0/0 → 2x/1 → 2)
73        let f = div(sub(pow(var("x"), 2), c(1.0)), sub(var("x"), c(1.0)));
74        assert!((limit(&f, "x", 1.0).unwrap() - 2.0).abs() < 1e-7);
75        // lim_{x→2} (x²−4)/(x−2) = 4
76        let g = div(sub(pow(var("x"), 2), c(4.0)), sub(var("x"), c(2.0)));
77        assert!((limit(&g, "x", 2.0).unwrap() - 4.0).abs() < 1e-7);
78    }
79
80    #[test]
81    fn direct_substitution_when_defined() {
82        // lim_{x→3} (x²+1) = 10
83        let f = add(pow(var("x"), 2), c(1.0));
84        assert!((limit(&f, "x", 3.0).unwrap() - 10.0).abs() < 1e-9);
85    }
86
87    #[test]
88    fn divergent_fails_closed() {
89        // lim_{x→0} 1/x has no finite limit.
90        let f = div(c(1.0), var("x"));
91        assert!(limit(&f, "x", 0.0).is_none());
92    }
93
94    #[test]
95    fn rational_limit_at_infinity() {
96        // (2x²+3)/(x²−1) → 2 as x→∞
97        let f = div(
98            add(mul(c(2.0), pow(var("x"), 2)), c(3.0)),
99            sub(pow(var("x"), 2), c(1.0)),
100        );
101        assert!((limit_at_infinity(&f, "x").unwrap() - 2.0).abs() < 1e-3);
102    }
103}