mathhook_core/calculus/integrals/rational/
linear.rs

1//! Integration of linear factors in partial fraction decomposition
2//!
3//! Implements Heaviside's method for computing coefficients of repeated linear factors.
4
5use crate::calculus::derivatives::Derivative;
6use crate::core::{Expression, Symbol};
7use crate::simplify::Simplify;
8
9use super::helpers::{factorial, substitute_variable};
10
11/// Integrate linear partial fraction using Heaviside's method
12///
13/// For simple pole (power=1): Uses cover-up method
14/// For repeated pole (power>1): Uses Heaviside's method with derivatives
15///
16/// # Mathematical Basis
17///
18/// For `P(x)/Q(x)` with `Q(x) = (x-r)^n · R(x)`, we want:
19/// ```text
20/// P(x)/Q(x) = A₁/(x-r) + A₂/(x-r)² + ... + Aₙ/(x-r)ⁿ + (other terms)
21/// ```
22///
23/// Heaviside's method: Define `g(x) = (x-r)ⁿ · P(x)/Q(x) = P(x)/R(x)`
24/// Then: `Aₖ = (1/(n-k)!) · [d^(n-k)/dx^(n-k) g(x)]|ₓ₌ᵣ`
25///
26/// # Integration Formulas
27///
28/// - `∫A/(x-r) dx = A·ln|x-r| + C`
29/// - `∫A/(x-r)ⁿ dx = -A/((n-1)(x-r)^(n-1)) + C` for `n > 1`
30///
31/// # Arguments
32///
33/// * `numerator` - Numerator polynomial `P(x)`
34/// * `denominator` - Full denominator `Q(x)`
35/// * `root` - The root `r` of the linear factor
36/// * `power` - Multiplicity `n` of the root
37/// * `var` - Integration variable
38///
39/// # Returns
40///
41/// Integrated expression or `None` if integration fails
42///
43/// # Examples
44///
45/// ```
46/// use mathhook_core::{Expression, symbol};
47/// use mathhook_core::calculus::integrals::rational::linear::integrate_linear_factor;
48/// use mathhook_core::simplify::Simplify;
49///
50/// let x = symbol!(x);
51///
52/// let numerator = Expression::integer(1);
53/// let denominator = Expression::pow(
54///     Expression::add(vec![
55///         Expression::symbol(x.clone()),
56///         Expression::integer(-2),
57///     ]),
58///     Expression::integer(2),
59/// );
60/// let root = Expression::integer(2);
61///
62/// let result = integrate_linear_factor(&numerator, &denominator, &root, 2, &x);
63/// assert!(result.is_some());
64/// ```
65pub fn integrate_linear_factor(
66    numerator: &Expression,
67    denominator: &Expression,
68    root: &Expression,
69    power: i64,
70    var: &Symbol,
71) -> Option<Expression> {
72    let mut result = Expression::integer(0);
73
74    for k in 1..=power {
75        let coeff = if power == 1 {
76            substitute_variable(numerator, var, root).simplify()
77        } else {
78            compute_heaviside_coefficient(numerator, denominator, root, power, k, var)?
79        };
80
81        let x_minus_r = Expression::add(vec![
82            Expression::symbol(var.clone()),
83            Expression::mul(vec![Expression::integer(-1), root.clone()]),
84        ]);
85
86        let term = if k == 1 {
87            Expression::mul(vec![
88                coeff,
89                Expression::function("ln", vec![Expression::function("abs", vec![x_minus_r])]),
90            ])
91        } else {
92            Expression::mul(vec![
93                Expression::integer(-1),
94                coeff,
95                Expression::rational(1, k - 1),
96                Expression::pow(x_minus_r, Expression::integer(-(k - 1))),
97            ])
98        };
99
100        result = Expression::add(vec![result, term]);
101    }
102
103    Some(result)
104}
105
106/// Compute Heaviside coefficient for repeated linear factors
107///
108/// For `(x-r)^n` in denominator, coefficient k is:
109/// ```text
110/// Aₖ = (1/(n-k)!) · [d^(n-k)/dx^(n-k) of g(x)]|ₓ₌ᵣ
111/// ```
112/// where `g(x) = (x-r)^n · P(x)/Q(x)`
113///
114/// # Arguments
115///
116/// * `numerator` - Numerator polynomial `P(x)`
117/// * `denominator` - Full denominator `Q(x)`
118/// * `root` - The root `r`
119/// * `total_power` - Multiplicity `n`
120/// * `k` - Coefficient index (1 to n)
121/// * `var` - Integration variable
122///
123/// # Returns
124///
125/// The coefficient `Aₖ` or `None` if computation fails
126///
127/// # Examples
128///
129/// ```
130/// use mathhook_core::{Expression, symbol};
131/// use mathhook_core::calculus::integrals::rational::linear::compute_heaviside_coefficient;
132///
133/// let x = symbol!(x);
134/// let numerator = Expression::integer(1);
135/// let denominator = Expression::pow(
136///     Expression::add(vec![
137///         Expression::symbol(x.clone()),
138///         Expression::integer(-1),
139///     ]),
140///     Expression::integer(3),
141/// );
142/// let root = Expression::integer(1);
143///
144/// let coeff = compute_heaviside_coefficient(&numerator, &denominator, &root, 3, 1, &x);
145/// assert!(coeff.is_some());
146/// ```
147pub fn compute_heaviside_coefficient(
148    numerator: &Expression,
149    denominator: &Expression,
150    root: &Expression,
151    total_power: i64,
152    k: i64,
153    var: &Symbol,
154) -> Option<Expression> {
155    let x_minus_r = Expression::add(vec![
156        Expression::symbol(var.clone()),
157        Expression::mul(vec![Expression::integer(-1), root.clone()]),
158    ]);
159    let x_minus_r_pow_n = Expression::pow(x_minus_r, Expression::integer(total_power));
160
161    let g = Expression::mul(vec![
162        x_minus_r_pow_n,
163        Expression::mul(vec![
164            numerator.clone(),
165            Expression::pow(denominator.clone(), Expression::integer(-1)),
166        ]),
167    ])
168    .simplify();
169
170    let derivative_order = total_power - k;
171
172    let derivative_result = if derivative_order == 0 {
173        g
174    } else {
175        g.nth_derivative(var.clone(), derivative_order as u32)
176    };
177
178    let evaluated = substitute_variable(&derivative_result, var, root).simplify();
179
180    let fact = factorial(derivative_order);
181    Some(Expression::mul(vec![
182        Expression::rational(1, fact),
183        evaluated,
184    ]))
185}
mathhook_core/calculus/integrals/rational/linear.rs

mathhook_core/calculus/integrals/rational/
linear.rs
