mathhook_core/algebra/solvers/polynomial/
solver.rs

1//! Core polynomial solver logic for cubic and quartic equations
2
3use crate::algebra::solvers::{EquationSolver, SolverResult};
4use crate::core::constants::EPSILON;
5use crate::core::{Expression, Number, Symbol};
6use crate::simplify::Simplify;
7use crate::symbol;
8
9/// Polynomial equation solver
10#[derive(Debug, Clone)]
11pub struct PolynomialSolver;
12
13impl Default for PolynomialSolver {
14    fn default() -> Self {
15        Self::new()
16    }
17}
18
19impl PolynomialSolver {
20    pub fn new() -> Self {
21        Self
22    }
23
24    /// Find the degree of polynomial in given variable
25    /// Static helper function for recursive degree calculation.
26    /// Does not require instance state, only performs expression traversal.
27    pub fn find_polynomial_degree(expr: &Expression, variable: &Symbol) -> u32 {
28        match expr {
29            Expression::Pow(base, exp) if **base == Expression::symbol(variable.clone()) => {
30                match exp.as_ref() {
31                    Expression::Number(Number::Integer(n)) => *n as u32,
32                    _ => 1,
33                }
34            }
35            Expression::Mul(factors) => factors
36                .iter()
37                .map(|f| Self::find_polynomial_degree(f, variable))
38                .max()
39                .unwrap_or(0),
40            Expression::Add(terms) => terms
41                .iter()
42                .map(|t| Self::find_polynomial_degree(t, variable))
43                .max()
44                .unwrap_or(0),
45            _ if *expr == Expression::symbol(variable.clone()) => 1,
46            _ => 0,
47        }
48    }
49
50    /// Solve cubic equation (simplified implementation)
51    pub fn solve_cubic(&self, equation: &Expression, variable: &Symbol) -> SolverResult {
52        if let Expression::Add(terms) = equation {
53            if terms.len() == 2 {
54                let (power_term, constant_term) = match (&terms[0], &terms[1]) {
55                    (Expression::Number(Number::Integer(c)), p @ Expression::Pow(..)) => (p, c),
56                    (p @ Expression::Pow(..), Expression::Number(Number::Integer(c))) => (p, c),
57                    _ => {
58                        return self.solve_cubic_rational_root_theorem(equation, variable);
59                    }
60                };
61
62                if let Expression::Pow(base, exp) = power_term {
63                    if **base == Expression::symbol(variable.clone())
64                        && **exp == Expression::integer(3)
65                    {
66                        let cube_root_value = (-constant_term as f64).cbrt();
67
68                        if cube_root_value.fract().abs() < EPSILON {
69                            let real_root = Expression::integer(cube_root_value as i64);
70                            return SolverResult::Partial(vec![real_root]);
71                        }
72                    }
73                }
74            }
75        }
76
77        self.solve_cubic_rational_root_theorem(equation, variable)
78    }
79
80    /// Try to solve cubic using rational root theorem
81    pub fn solve_cubic_rational_root_theorem(
82        &self,
83        equation: &Expression,
84        variable: &Symbol,
85    ) -> SolverResult {
86        let potential_roots = vec![-3, -2, -1, 0, 1, 2, 3];
87        let mut found_roots = Vec::new();
88
89        for &root in &potential_roots {
90            let test_value = Expression::integer(root);
91            if self
92                .evaluate_polynomial_at(equation, variable, &test_value)
93                .is_zero()
94            {
95                found_roots.push(Expression::integer(root));
96            }
97        }
98
99        if found_roots.len() >= 3 {
100            SolverResult::Multiple(found_roots)
101        } else if !found_roots.is_empty() {
102            SolverResult::Partial(found_roots)
103        } else {
104            SolverResult::NoSolution
105        }
106    }
107
108    /// Solve quartic equation (simplified implementation)
109    pub fn solve_quartic(&self, equation: &Expression, variable: &Symbol) -> SolverResult {
110        if let Expression::Add(terms) = equation {
111            if terms.len() == 2 {
112                let (power_term, constant_term) = match (&terms[0], &terms[1]) {
113                    (Expression::Number(Number::Integer(c)), p @ Expression::Pow(..)) => (p, c),
114                    (p @ Expression::Pow(..), Expression::Number(Number::Integer(c))) => (p, c),
115                    _ => {
116                        return self.solve_quartic_rational_root_theorem(equation, variable);
117                    }
118                };
119
120                if let Expression::Pow(base, exp) = power_term {
121                    if **base == Expression::symbol(variable.clone())
122                        && **exp == Expression::integer(4)
123                    {
124                        let fourth_root_value = (-constant_term as f64).powf(0.25);
125
126                        if fourth_root_value.is_finite() {
127                            let real_root = fourth_root_value.abs();
128
129                            let real_root_expr = if real_root.fract().abs() < EPSILON {
130                                Expression::integer(real_root as i64)
131                            } else {
132                                Expression::Number(Number::float(real_root))
133                            };
134
135                            let neg_real_root_expr = if real_root.fract().abs() < EPSILON {
136                                Expression::integer(-(real_root as i64))
137                            } else {
138                                Expression::Number(Number::float(-real_root))
139                            };
140
141                            return SolverResult::Partial(vec![real_root_expr, neg_real_root_expr]);
142                        }
143                    }
144                }
145            }
146        }
147
148        self.solve_quartic_rational_root_theorem(equation, variable)
149    }
150
151    /// Try to solve quartic using rational root theorem
152    pub fn solve_quartic_rational_root_theorem(
153        &self,
154        equation: &Expression,
155        variable: &Symbol,
156    ) -> SolverResult {
157        let potential_roots = vec![-3, -2, -1, 0, 1, 2, 3];
158        let mut found_roots = Vec::new();
159
160        for &root in &potential_roots {
161            let test_value = Expression::integer(root);
162            if self
163                .evaluate_polynomial_at(equation, variable, &test_value)
164                .is_zero()
165            {
166                found_roots.push(Expression::integer(root));
167            }
168        }
169
170        if found_roots.len() >= 4 {
171            SolverResult::Multiple(found_roots)
172        } else if !found_roots.is_empty() {
173            SolverResult::Partial(found_roots)
174        } else {
175            SolverResult::NoSolution
176        }
177    }
178
179    /// Evaluate polynomial at specific value
180    pub fn evaluate_polynomial_at(
181        &self,
182        poly: &Expression,
183        variable: &Symbol,
184        value: &Expression,
185    ) -> Expression {
186        Self::substitute_variable(poly, variable, value).simplify()
187    }
188
189    /// Substitute variable with value in expression
190    ///
191    /// Static helper function for recursive substitution.
192    /// Does not require instance state, only performs expression traversal and replacement.
193    pub fn substitute_variable(
194        expr: &Expression,
195        variable: &Symbol,
196        value: &Expression,
197    ) -> Expression {
198        match expr {
199            _ if *expr == Expression::symbol(variable.clone()) => value.clone(),
200            Expression::Add(terms) => {
201                let new_terms: Vec<Expression> = terms
202                    .iter()
203                    .map(|t| Self::substitute_variable(t, variable, value))
204                    .collect();
205                Expression::add(new_terms).simplify()
206            }
207            Expression::Mul(factors) => {
208                let new_factors: Vec<Expression> = factors
209                    .iter()
210                    .map(|f| Self::substitute_variable(f, variable, value))
211                    .collect();
212                Expression::mul(new_factors).simplify()
213            }
214            Expression::Pow(base, exp) => {
215                let new_base = Self::substitute_variable(base, variable, value);
216                let new_exp = Self::substitute_variable(exp, variable, value);
217                Expression::pow(new_base, new_exp).simplify()
218            }
219            Expression::Function { name, args } => {
220                let new_args: Vec<Expression> = args
221                    .iter()
222                    .map(|a| Self::substitute_variable(a, variable, value))
223                    .collect();
224                Expression::function(name, new_args).simplify()
225            }
226            _ => expr.clone(),
227        }
228    }
229
230    /// Extract constant term and leading coefficient from polynomial
231    pub fn extract_constant_and_leading(
232        &self,
233        poly: &Expression,
234        _variable: &Symbol,
235    ) -> (i64, i64) {
236        match poly {
237            Expression::Add(terms) => {
238                let mut constant = 0i64;
239                let mut leading = 1i64;
240
241                for term in terms.iter() {
242                    if let Expression::Number(Number::Integer(n)) = term {
243                        constant = *n;
244                    } else if let Expression::Pow(_, exp) = term {
245                        if let Expression::Number(Number::Integer(_)) = exp.as_ref() {
246                            leading = 1;
247                        }
248                    }
249                }
250
251                (constant, leading)
252            }
253            _ => (0, 1),
254        }
255    }
256
257    /// Get all divisors of a number
258    pub fn get_divisors(&self, n: i64) -> Vec<i64> {
259        if n == 0 {
260            return vec![1];
261        }
262
263        let n = n.abs();
264        let mut divisors = Vec::new();
265
266        for i in 1..=n {
267            if i * i > n {
268                break;
269            }
270            if n % i == 0 {
271                divisors.push(i);
272                if i != n / i {
273                    divisors.push(n / i);
274                }
275            }
276        }
277
278        divisors.sort();
279        divisors
280    }
281}
282
283impl EquationSolver for PolynomialSolver {
284    #[inline(always)]
285    fn solve(&self, equation: &Expression, variable: &Symbol) -> SolverResult {
286        let degree = Self::find_polynomial_degree(equation, variable);
287
288        match degree {
289            3 => self.solve_cubic(equation, variable),
290            4 => self.solve_quartic(equation, variable),
291            _ => SolverResult::NoSolution,
292        }
293    }
294
295    fn solve_with_explanation(
296        &self,
297        equation: &Expression,
298        variable: &Symbol,
299    ) -> (
300        SolverResult,
301        crate::educational::step_by_step::StepByStepExplanation,
302    ) {
303        super::educational::solve_with_explanation(self, equation, variable)
304    }
305
306    fn can_solve(&self, equation: &Expression) -> bool {
307        let degree = Self::find_polynomial_degree(equation, &symbol!(x));
308        (3..=4).contains(&degree)
309    }
310}