ocas_poly/factor/mod.rs
1//! Polynomial factorization algorithms.
2//!
3//! This module groups square-free factorization with complete factorization
4//! over finite fields ([`finite_field`]) and, for lifting back to the integers,
5//! Hensel lifting ([`hensel`]).
6//!
7//! The top-level entry point for factoring a univariate polynomial over
8//! $\mathbb{Z}$ is [`DenseUnivariatePolynomial::factor`](crate::DenseUnivariatePolynomial::factor),
9//! and over a finite field
10//! [`factor_over_finite_field`](finite_field::factor_over_finite_field).
11
12use ocas_domain::EuclideanDomain;
13use ocas_domain::{FiniteField, IntegerDomain};
14
15use crate::dense::DenseUnivariatePolynomial;
16
17pub mod finite_field;
18pub mod hensel;
19pub mod multivariate;
20
21/// Result of a square-free factorization: list of (factor, multiplicity) pairs.
22pub type SquareFreeFactors<D> = Vec<(DenseUnivariatePolynomial<D>, usize)>;
23
24/// Result of a complete factorization: list of (factor, multiplicity) pairs
25/// where each factor is irreducible (or, over the integers, primitive and
26/// irreducible over $\mathbb{Q}$).
27pub type Factors<D> = Vec<(DenseUnivariatePolynomial<D>, usize)>;
28
29impl<D: EuclideanDomain> DenseUnivariatePolynomial<D> {
30 /// Compute the square-free factorization of this polynomial.
31 ///
32 /// Returns a list of (factor, multiplicity) pairs.
33 /// For example, `(x+1)^2 * (x-1)` yields `[(x+1, 2), (x-1, 1)]`.
34 ///
35 /// # Example
36 ///
37 /// ```
38 /// use ocas_domain::{IntegerDomain, Integer};
39 /// use ocas_poly::DenseUnivariatePolynomial;
40 ///
41 /// let d = IntegerDomain;
42 /// // (x+1)^2*(x-1) = x^3 + x^2 - x - 1
43 /// let p = DenseUnivariatePolynomial::from_coeffs(d, vec![
44 /// Integer::from(-1), Integer::from(-1), Integer::from(1), Integer::from(1),
45 /// ]);
46 /// let factors = p.square_free_factorization();
47 /// assert_eq!(factors.len(), 2);
48 /// ```
49 pub fn square_free_factorization(&self) -> SquareFreeFactors<D> {
50 let mut factors = SquareFreeFactors::new();
51 if self.is_zero() {
52 return factors;
53 }
54
55 // Step 1: make polynomial primitive.
56 let f = self.primitive_part();
57 let f_deriv = f.derivative();
58
59 // g = gcd(f, f')
60 let mut g = f.gcd(&f_deriv);
61 if g.is_zero() {
62 return factors;
63 }
64
65 // w = f / g contains each distinct irreducible factor exactly once.
66 let mut w = match f.div_rem(&g) {
67 Some((q, _)) => q,
68 None => return factors,
69 };
70
71 let mut k = 1;
72 while !w.is_one() {
73 // h = gcd(w, g)
74 let h = w.gcd(&g);
75 // z = w / h is the factor with multiplicity k.
76 if let Some((z, _)) = w.div_rem(&h)
77 && !z.is_one()
78 && !z.is_zero()
79 {
80 factors.push((z, k));
81 }
82
83 // Prepare for next iteration.
84 if let Some((q, _)) = g.div_rem(&h) {
85 g = q;
86 } else {
87 break;
88 }
89 w = h;
90 k += 1;
91 }
92
93 factors
94 }
95
96 /// Check whether this polynomial is square-free.
97 ///
98 /// A polynomial is square-free if gcd(p, p') = 1.
99 pub fn is_square_free(&self) -> bool {
100 if self.degree().unwrap_or(0) <= 1 {
101 return true;
102 }
103 let deriv = self.derivative();
104 let g = self.gcd(&deriv);
105 g.degree() == Some(0)
106 }
107}
108
109// ── factor() for integer polynomials ──────────────────────────────
110
111impl DenseUnivariatePolynomial<IntegerDomain> {
112 /// Completely factor this primitive integer polynomial into monic
113 /// irreducible factors with multiplicities.
114 ///
115 /// The input must be primitive (coefficient content = 1). Use
116 /// [`primitive_part`](crate::DenseUnivariatePolynomial::primitive_part)
117 /// to prepare an arbitrary integer polynomial before factoring.
118 ///
119 /// # Example
120 ///
121 /// ```
122 /// use ocas_domain::{Integer, IntegerDomain};
123 /// use ocas_poly::DenseUnivariatePolynomial;
124 ///
125 /// let d = IntegerDomain;
126 /// // x^2 - 1 = (x-1)(x+1)
127 /// let p = DenseUnivariatePolynomial::from_coeffs(d, vec![
128 /// Integer::from(-1), Integer::from(0), Integer::from(1),
129 /// ]);
130 /// let factors = p.factor();
131 /// assert_eq!(factors.len(), 2);
132 /// ```
133 pub fn factor(&self) -> Factors<IntegerDomain> {
134 hensel::factor_primitive(self)
135 }
136}
137
138// ── factor() for finite-field polynomials ─────────────────────────
139
140impl DenseUnivariatePolynomial<FiniteField> {
141 /// Completely factor this univariate polynomial over $\mathbb{F}_p$
142 /// into monic irreducible factors with multiplicities.
143 ///
144 /// # Example
145 ///
146 /// ```
147 /// use num_bigint::BigInt;
148 /// use ocas_domain::{Domain, FiniteField};
149 /// use ocas_poly::DenseUnivariatePolynomial;
150 ///
151 /// let f = FiniteField::new(BigInt::from(5));
152 /// // x^2 - 1 over F_5
153 /// let p = DenseUnivariatePolynomial::from_coeffs(
154 /// f.clone(), vec![f.element(4), f.element(0), f.element(1)]);
155 /// let factors = p.factor();
156 /// assert!(!factors.is_empty());
157 /// ```
158 pub fn factor(&self) -> Factors<FiniteField> {
159 finite_field::factor_over_finite_field(self)
160 }
161}
162
163#[cfg(test)]
164mod tests {
165 use super::*;
166 use ocas_domain::{Integer, IntegerDomain};
167
168 fn i(n: i64) -> Integer {
169 Integer::from(n)
170 }
171
172 #[test]
173 fn square_free_x_plus_1_cubed() {
174 let d = IntegerDomain;
175 // (x+1)^3 = x^3 + 3x^2 + 3x + 1
176 let p = DenseUnivariatePolynomial::from_coeffs(d, vec![i(1), i(3), i(3), i(1)]);
177 let factors = p.square_free_factorization();
178 assert!(!factors.is_empty());
179 for (factor, mult) in &factors {
180 if factor.degree() == Some(1) {
181 assert_eq!(*mult, 3);
182 }
183 }
184 }
185
186 #[test]
187 fn is_square_free_linear() {
188 let d = IntegerDomain;
189 let p = DenseUnivariatePolynomial::from_coeffs(d, vec![i(1), i(1)]); // x+1
190 assert!(p.is_square_free());
191 }
192
193 #[test]
194 fn is_not_square_free_perfect_square() {
195 let d = IntegerDomain;
196 let p = DenseUnivariatePolynomial::from_coeffs(d, vec![i(1), i(2), i(1)]); // (x+1)^2
197 assert!(!p.is_square_free());
198 }
199
200 #[test]
201 fn square_free_x2_minus_1() {
202 let d = IntegerDomain;
203 let p = DenseUnivariatePolynomial::from_coeffs(d, vec![i(-1), i(0), i(1)]); // x^2-1
204 let factors = p.square_free_factorization();
205 assert!(p.is_square_free());
206 assert!(!factors.is_empty());
207 }
208}