algebraeon_rings/finite_fields/
polynomial.rs

1use crate::{matrix::*, polynomial::*, structure::*};
2use algebraeon_nzq::*;
3use algebraeon_sets::structure::*;
4use itertools::Itertools;
5
6// Useful: https://en.wikipedia.org/wiki/Factorization_of_polynomials_over_finite_fields
7/*
8Factorization over finite fields typically happens in the following steps:
91. Monic factorization: Remove a unit so that the polynomial to be factored is monic.
102. Squarefree factorization: Factor into squarefree factors.
113. Distinct degree factorization: Factor a squarefree polynomial into polynomials such that all irreducible factors of each polynomial have equal degree.
124. Equal degree factorization: Factor a squarefree polynomial each of whose irreducible factors all have the same known degree d.
13
141 is trivial.
15There is a standard algorithm for 2.
16Berlekamps algorithm does 3 and 4 at the same time.
17There is a standard algorithm for 3.
18Cantor–Zassenhaus algorithm does 4.
19*/
20
21/// Store a monic factorization
22#[derive(Debug, Clone)]
23pub struct MonicFactored<FS: FieldSignature, FSB: BorrowedStructure<FS>> {
24    poly_ring: PolynomialStructure<FS, FSB>,
25    unit: FS::Set,              // a unit
26    monic: Polynomial<FS::Set>, // a monic polynomial
27}
28
29/// Store a squarefree factorization
30#[derive(Debug, Clone)]
31pub struct SquarefreeFactored<FS: FiniteFieldSignature, FSB: BorrowedStructure<FS>> {
32    poly_ring: PolynomialStructure<FS, FSB>,
33    unit: FS::Set,                                           // a unit
34    squarefree_factors: Vec<(Polynomial<FS::Set>, Natural)>, // squarefree monic polynomials and their multiplicities
35}
36
37#[derive(Debug, Clone)]
38struct DistinctDegreeFactor<FS: FiniteFieldSignature> {
39    irreducible_factor_degree: usize, // the degree of the irreducible factors
40    polynomial: Polynomial<FS::Set>, // a monic squarefree polynomial equal to a product of irreducibles all of the same degree
41}
42/// Store a distinct degree factorization
43#[derive(Debug, Clone)]
44pub struct DistinctDegreeFactored<FS: FiniteFieldSignature, FSB: BorrowedStructure<FS>> {
45    poly_ring: PolynomialStructure<FS, FSB>,
46    unit: FS::Set, // a unit
47    distinct_degree_factors: Vec<(DistinctDegreeFactor<FS>, Natural)>,
48}
49
50impl<FS: FieldSignature, FSB: BorrowedStructure<FS>> MonicFactored<FS, FSB> {
51    // pub fn factor(poly_ring: PolynomialStructure<FS>, poly: Polynomial<FS::Set>) -> Self {
52    //     let (unit, monic) = poly_ring.factor_fav_assoc(&poly);
53    //     let unit = poly_ring.as_constant(&unit).unwrap();
54    //     debug_assert!(poly_ring.is_monic(&monic));
55    //     Self {
56    //         poly_ring,
57    //         unit,
58    //         monic,
59    //     }
60    // }
61
62    pub fn new_monic_unchecked(
63        poly_ring: PolynomialStructure<FS, FSB>,
64        monic: Polynomial<FS::Set>,
65    ) -> Self {
66        debug_assert!(poly_ring.is_monic(&monic));
67        let unit = poly_ring.coeff_ring().one();
68        Self {
69            poly_ring,
70            unit,
71            monic,
72        }
73    }
74
75    pub fn scalar_part(&self) -> &FS::Set {
76        &self.unit
77    }
78
79    pub fn monic_part(&self) -> &Polynomial<FS::Set> {
80        &self.monic
81    }
82
83    pub fn into_polynomial(self) -> Polynomial<FS::Set> {
84        self.poly_ring
85            .mul(&Polynomial::constant(self.unit), &self.monic)
86    }
87}
88
89impl<FS: FiniteFieldSignature, FSB: BorrowedStructure<FS>> SquarefreeFactored<FS, FSB> {
90    pub fn unit_unchecked(poly_ring: PolynomialStructure<FS, FSB>, unit: FS::Set) -> Self {
91        debug_assert!(poly_ring.coeff_ring().is_unit(&unit));
92        Self {
93            poly_ring,
94            unit,
95            squarefree_factors: vec![],
96        }
97    }
98
99    pub fn new_squarefree_poly_unchecked(
100        poly_ring: PolynomialStructure<FS, FSB>,
101        poly: Polynomial<FS::Set>,
102    ) -> Self {
103        debug_assert!(!poly_ring.is_zero(&poly));
104        let deg = poly_ring.degree(&poly).unwrap();
105        if deg == 0 {
106            let unit = poly_ring.as_constant(&poly).unwrap();
107            Self {
108                poly_ring,
109                unit,
110                squarefree_factors: vec![],
111            }
112        } else {
113            let unit = poly_ring.coeff_ring().one();
114            Self {
115                poly_ring,
116                unit,
117                squarefree_factors: vec![(poly, Natural::ONE)],
118            }
119        }
120    }
121
122    pub fn pow(mut self, power: &Natural) -> Self {
123        for (_, poly_power) in &mut self.squarefree_factors {
124            *poly_power *= power;
125        }
126        self
127    }
128
129    pub fn into_monic_factored(self) -> MonicFactored<FS, FSB> {
130        let monic = self.poly_ring.product(
131            self.squarefree_factors
132                .iter()
133                .map(|(f, k)| self.poly_ring.nat_pow(f, k))
134                .collect(),
135        );
136        MonicFactored {
137            poly_ring: self.poly_ring,
138            unit: self.unit,
139            monic,
140        }
141    }
142
143    fn mul(&mut self, other: &Self) {
144        assert_eq!(self.poly_ring, other.poly_ring);
145        let poly_ring = &self.poly_ring;
146        self.unit = poly_ring.coeff_ring().mul(&self.unit, &other.unit);
147        'OTHER_LOOP: for (poly, power) in &other.squarefree_factors {
148            for (self_poly, self_power) in &mut self.squarefree_factors {
149                if poly_ring.equal(poly, self_poly) {
150                    *self_power += power;
151                    continue 'OTHER_LOOP;
152                }
153            }
154            self.squarefree_factors.push((poly.clone(), power.clone()));
155        }
156    }
157}
158
159impl<FS: FiniteFieldSignature, FSB: BorrowedStructure<FS>> DistinctDegreeFactored<FS, FSB> {
160    pub fn into_squarefree_factored(self) -> SquarefreeFactored<FS, FSB> {
161        let squarefree_factors = self
162            .distinct_degree_factors
163            .into_iter()
164            .map(|(ddf, k)| (ddf.polynomial, k))
165            .collect();
166        SquarefreeFactored {
167            poly_ring: self.poly_ring,
168            unit: self.unit,
169            squarefree_factors,
170        }
171    }
172}
173
174impl<
175    FS: FiniteFieldSignature,
176    FSB: BorrowedStructure<FS>,
177    FSPB: BorrowedStructure<PolynomialStructure<FS, FSB>>,
178    NB: BorrowedStructure<NaturalCanonicalStructure>,
179> FactoringStructure<PolynomialStructure<FS, FSB>, FSPB, NaturalCanonicalStructure, NB>
180where
181    PolynomialStructure<FS, FSB>: SetSignature<Set = Polynomial<FS::Set>>
182        + FactoringMonoidSignature<FactoredExponent = NaturalCanonicalStructure>,
183{
184    pub fn into_distinct_degree_factored(
185        &self,
186        a: NonZeroFactored<Polynomial<FS::Set>, Natural>,
187    ) -> DistinctDegreeFactored<FS, FSB> {
188        let poly_ring = self.objects().clone();
189        let (unit, factors) = a.into_unit_and_powers();
190        let unit = poly_ring.as_constant(&unit).unwrap();
191        let distinct_degree_factors = factors
192            .into_iter()
193            .map(|(f, k)| {
194                (
195                    DistinctDegreeFactor {
196                        irreducible_factor_degree: poly_ring.degree(&f).unwrap(),
197                        polynomial: f,
198                    },
199                    k,
200                )
201            })
202            .collect();
203        DistinctDegreeFactored {
204            poly_ring,
205            unit,
206            distinct_degree_factors,
207        }
208    }
209}
210
211impl<FS: FieldSignature, FSB: BorrowedStructure<FS>> PolynomialStructure<FS, FSB>
212where
213    PolynomialStructure<FS, FSB>: SetSignature<Set = Polynomial<FS::Set>>,
214{
215    /// monic factorization
216    pub fn factorize_monic(&self, poly: &Polynomial<FS::Set>) -> Option<MonicFactored<FS, FSB>> {
217        if self.is_zero(poly) {
218            None
219        } else {
220            let (unit, monic) = self.factor_fav_assoc(poly);
221            let unit = self.as_constant(&unit).unwrap();
222            Some(MonicFactored {
223                poly_ring: self.clone(),
224                unit,
225                monic,
226            })
227        }
228    }
229}
230
231impl<F: MetaType> Polynomial<F>
232where
233    F::Signature: FiniteFieldSignature,
234    PolynomialStructure<F::Signature, F::Signature>: SetSignature<Set = Polynomial<F>>,
235    Self: MetaType<Signature = PolynomialStructure<F::Signature, F::Signature>>,
236{
237    pub fn factorize_monic(&self) -> Option<MonicFactored<F::Signature, F::Signature>> {
238        Self::structure().factorize_monic(self)
239    }
240}
241
242impl<FS: FiniteFieldSignature, FSB: BorrowedStructure<FS>> MonicFactored<FS, FSB>
243where
244    PolynomialStructure<FS, FSB>: SetSignature<Set = Polynomial<FS::Set>>,
245{
246    /// squarefree factorization
247    pub fn factorize_squarefree(&self) -> SquarefreeFactored<FS, FSB> {
248        let mut factors =
249            SquarefreeFactored::unit_unchecked(self.poly_ring.clone(), self.unit.clone());
250
251        //let w be the squarefree product of all factors of f of multiplicity not divisible by p
252        let mut c = self
253            .poly_ring
254            .gcd(&self.monic, &self.poly_ring.derivative(self.monic.clone()));
255        let mut w = self.poly_ring.try_divide(&self.monic, &c).unwrap();
256
257        //find all the factors in w
258        let mut i = Natural::ONE;
259        while !self.poly_ring.equal(&w, &self.poly_ring.one()) {
260            let y = self.poly_ring.gcd(&w, &c);
261            factors.mul(
262                &SquarefreeFactored::new_squarefree_poly_unchecked(
263                    self.poly_ring.clone(),
264                    self.poly_ring.try_divide(&w, &y).unwrap(),
265                )
266                .pow(&i),
267            );
268            let c_over_y = self.poly_ring.try_divide(&c, &y).unwrap();
269            (w, c) = (y, c_over_y);
270            i += Natural::ONE;
271        }
272
273        // c is now the product, with multiplicity, of the remaining factors of f
274        if !self.poly_ring.equal(&c, &self.poly_ring.one()) {
275            //c = c^{1/p}
276            let p = self.poly_ring.coeff_ring().characteristic_and_power().0;
277            let mut reduced_c_coeffs = vec![];
278            for (k, coeff) in self.poly_ring.into_coeffs(c).into_iter().enumerate() {
279                if Natural::from(k) % &p == Natural::ZERO {
280                    reduced_c_coeffs.push(coeff.clone());
281                } else {
282                    debug_assert!(self.poly_ring.coeff_ring().is_zero(&coeff));
283                }
284            }
285            let reduced_c: Polynomial<<FS as SetSignature>::Set> =
286                Polynomial::from_coeffs(reduced_c_coeffs);
287            factors.mul(
288                &MonicFactored::new_monic_unchecked(self.poly_ring.clone(), reduced_c)
289                    .factorize_squarefree()
290                    .pow(&p),
291            );
292        }
293        factors
294    }
295}
296
297impl<FS: FiniteFieldSignature, FSB: BorrowedStructure<FS>> PolynomialStructure<FS, FSB>
298where
299    PolynomialStructure<FS, FSB>: SetSignature<Set = Polynomial<FS::Set>>,
300{
301    fn find_factor_by_berlekamps_algorithm(
302        &self,
303        f: Polynomial<FS::Set>,
304    ) -> FindFactorResult<Polynomial<FS::Set>> {
305        debug_assert!(self.is_squarefree(&f));
306
307        let f_deg = self.degree(&f).unwrap();
308        let all_elems = self.coeff_ring().list_all_elements();
309        let q = all_elems.len();
310
311        let row_polys = (0..f_deg)
312            .map(|i| {
313                self.rem(
314                    &self.add(&self.var_pow(i * q), &self.neg(&self.var_pow(i))),
315                    &f,
316                )
317            })
318            .collect::<Vec<_>>();
319        let mat = Matrix::construct(f_deg, f_deg, |i, j| {
320            self.coeff(&row_polys[i], j).into_owned()
321        });
322        // mat.pprint();
323        //the column kernel gives a basis of berlekamp subalgebra - all polynomials g such that g^q=g
324        let mat_struct = MatrixStructure::<FS, _>::new(self.coeff_ring());
325        let ker = mat_struct.row_kernel(mat);
326        // ker.pprint();
327        let ker_rank = ker.rank();
328        let ker_basis = ker
329            .basis()
330            .into_iter()
331            .map(|col| Polynomial::from_coeffs((0..f_deg).map(|c| col[c].clone()).collect()))
332            .collect::<Vec<_>>();
333
334        for berlekamp_subspace_coeffs in (0..ker_rank)
335            .map(|_i| all_elems.clone().into_iter())
336            .multi_cartesian_product()
337        {
338            let h = self.sum(
339                (0..ker_rank)
340                    .map(|i| {
341                        self.mul(
342                            &Polynomial::constant(berlekamp_subspace_coeffs[i].clone()),
343                            &ker_basis[i],
344                        )
345                    })
346                    .collect(),
347            );
348            //g is a possible non-trivial factor
349            let g = self.gcd(&h, &f);
350            // println!("g = {}", g);
351            let g_deg = self.degree(&g).unwrap();
352            if g_deg != 0 && g_deg != f_deg {
353                let f_over_g = self.try_divide(&f, &g).unwrap();
354                return FindFactorResult::Composite(g, f_over_g);
355            }
356        }
357        FindFactorResult::Irreducible
358    }
359}
360
361impl<FS: FiniteFieldSignature, FSB: BorrowedStructure<FS>> SquarefreeFactored<FS, FSB>
362where
363    PolynomialStructure<FS, FSB>: SetSignature<Set = Polynomial<FS::Set>>,
364{
365    /// use Berlekamps algorithm for a full factorization from a squarefree
366    pub fn factorize_berlekamps(&self) -> Factored<Polynomial<FS::Set>, Natural>
367    where
368        PolynomialStructure<FS, FSB>:
369            FactoringMonoidSignature<FactoredExponent = NaturalCanonicalStructure>,
370    {
371        let mut factors = self
372            .poly_ring
373            .factorizations()
374            .new_unit_and_powers_unchecked(Polynomial::constant(self.unit.clone()), vec![]);
375        for (sqfree_poly, power) in &self.squarefree_factors {
376            self.poly_ring.factorizations().mul_mut(
377                &mut factors,
378                &self.poly_ring.factorizations().pow(
379                    &factorize_by_find_factor(&self.poly_ring, sqfree_poly.clone(), &|f| {
380                        self.poly_ring.find_factor_by_berlekamps_algorithm(f)
381                    }),
382                    power,
383                ),
384            );
385        }
386        factors
387    }
388}
389
390impl<FS: FiniteFieldSignature, FSB: BorrowedStructure<FS>> SquarefreeFactored<FS, FSB>
391where
392    PolynomialStructure<FS, FSB>: SetSignature<Set = Polynomial<FS::Set>>,
393{
394    /// distinct degree factorization
395    pub fn factorize_distinct_degree(&self) -> DistinctDegreeFactored<FS, FSB> {
396        // https://en.wikipedia.org/wiki/Factorization_of_polynomials_over_finite_fields#Distinct-degree_factorization
397        let (p, k) = self.poly_ring.coeff_ring().characteristic_and_power();
398        let q = p.nat_pow(&k);
399        let mut distinct_degree_factors = vec![];
400        for (poly, sqfree_poly_multiplicity) in &self.squarefree_factors {
401            // a key step in this algorithm is the computation of g = gcd(f, x^{q^i} - x)
402            // the naive approach isn't great since x^{q^i} gets quite big
403            // instead the GCD can be computed after first reducing x^{q^i} modulo poly, since f divides poly
404            // to quickly compute x^{q^i} we can compute x^{q^1}, x^{q^2}, ... in successive loops, each time raising the previous value to the power of q
405            // raising polynomials mod poly to the power of q is a linear map (freshmans dream over finite field)
406            // so the qth power can be obtained by pre-computing the qth power matrix for polynomials mod poly
407
408            let n = self.poly_ring.degree(poly).unwrap();
409            debug_assert!(n >= 1);
410
411            let mod_poly_ring = self.poly_ring.quotient_ring(poly.clone());
412            let mat_structure = MatrixStructure::new(self.poly_ring.coeff_ring().clone());
413            let xq = mod_poly_ring.nat_pow(&self.poly_ring.var(), &q);
414            let qth_power_matrix = Matrix::join_cols(
415                n,
416                (0..n)
417                    .map(|c| {
418                        //compute (x^c)^q mod poly as a length n column vector of coefficients
419                        mod_poly_ring.to_col(&mod_poly_ring.nat_pow(&self.poly_ring.var_pow(c), &q))
420                    })
421                    .collect(),
422            );
423
424            let mut i = 1;
425            let mut xqi = xq.clone(); // x^{q^i} mod poly
426            let mut f = poly.clone();
427            while self.poly_ring.degree(&f).unwrap() >= 2 * i {
428                debug_assert!(self.poly_ring.is_monic(&f));
429
430                let g = self
431                    .poly_ring
432                    .factorize_monic(
433                        &self.poly_ring.gcd_by_primitive_subresultant(
434                            f.clone(),
435                            self.poly_ring
436                                .add(&xqi, &self.poly_ring.neg(&self.poly_ring.var())),
437                        ),
438                    )
439                    .unwrap()
440                    .monic;
441                debug_assert!(self.poly_ring.is_monic(&g));
442
443                #[cfg(debug_assertions)]
444                {
445                    debug_assert!(
446                        mod_poly_ring.equal(
447                            &xqi,
448                            &self
449                                .poly_ring
450                                .var_pow(q.nat_pow(&i.into()).try_into().unwrap())
451                        )
452                    );
453                    let g_naive = self
454                        .poly_ring
455                        .factorize_monic(
456                            &self.poly_ring.gcd_by_primitive_subresultant(
457                                f.clone(),
458                                self.poly_ring.add(
459                                    &self
460                                        .poly_ring
461                                        .var_pow(q.nat_pow(&i.into()).try_into().unwrap()),
462                                    &self.poly_ring.neg(&self.poly_ring.var()),
463                                ),
464                            ),
465                        )
466                        .unwrap()
467                        .monic;
468                    // This GCD was the naive and expensive way to do it. Use it to verify the fast and fancy way.
469                    debug_assert!(self.poly_ring.is_monic(&g_naive));
470                    debug_assert!(self.poly_ring.equal(&g, &g_naive));
471                }
472
473                if !self.poly_ring.equal(&g, &self.poly_ring.one()) {
474                    f = self.poly_ring.try_divide(&f, &g).unwrap();
475                    distinct_degree_factors.push((
476                        DistinctDegreeFactor {
477                            irreducible_factor_degree: i,
478                            polynomial: g,
479                        },
480                        sqfree_poly_multiplicity.clone(),
481                    ));
482                }
483                i += 1;
484                xqi = mod_poly_ring.from_col(
485                    mat_structure
486                        .mul(&qth_power_matrix, &mod_poly_ring.to_col(&xqi))
487                        .unwrap(),
488                );
489            }
490            if !self.poly_ring.equal(&f, &self.poly_ring.one()) {
491                distinct_degree_factors.push((
492                    DistinctDegreeFactor {
493                        irreducible_factor_degree: self.poly_ring.degree(&f).unwrap(),
494                        polynomial: f,
495                    },
496                    sqfree_poly_multiplicity.clone(),
497                ));
498            }
499        }
500        DistinctDegreeFactored {
501            poly_ring: self.poly_ring.clone(),
502            unit: self.unit.clone(),
503            distinct_degree_factors,
504        }
505    }
506}
507
508impl<FS: FiniteFieldSignature, FSB: BorrowedStructure<FS>> DistinctDegreeFactored<FS, FSB>
509where
510    PolynomialStructure<FS, FSB>: SetSignature<Set = Polynomial<FS::Set>>,
511{
512    /// Cantor–Zassenhaus algorithm for equal degree factorization
513    pub fn factorize_cantor_zassenhaus(&self) -> Factored<Polynomial<FS::Set>, Natural>
514    where
515        PolynomialStructure<FS, FSB>:
516            FactoringMonoidSignature<FactoredExponent = NaturalCanonicalStructure>,
517    {
518        let poly_ring = &self.poly_ring;
519        let mut fs = poly_ring
520            .factorizations()
521            .new_unit_unchecked(Polynomial::constant(self.unit.clone()));
522        for (ddf, mult) in &self.distinct_degree_factors {
523            let d = ddf.irreducible_factor_degree;
524            let n = self.poly_ring.degree(&ddf.polynomial).unwrap();
525            debug_assert_eq!(n % d, 0);
526
527            let finite_field = self.poly_ring.coeff_ring();
528            let (p, k) = finite_field.characteristic_and_power();
529            let q = p.nat_pow(&k);
530            let mut prand_elements = finite_field.generate_random_elements(0);
531
532            let mut to_factor = vec![ddf.polynomial.clone()];
533            loop {
534                // println!(
535                //     "{} {:?}",
536                //     d,
537                //     to_factor
538                //         .iter()
539                //         .map(|u| self.poly_ring.degree(u).unwrap())
540                //         .collect::<Vec<_>>()
541                // );
542                // Any polynomial in to_factor of degree d is irreducible.
543                to_factor = to_factor
544                    .into_iter()
545                    .filter_map(|f| {
546                        if self.poly_ring.degree(&f).unwrap() == d {
547                            poly_ring.factorizations().mul_mut(
548                                &mut fs,
549                                &poly_ring.factorizations().pow(
550                                    &poly_ring.factorizations().new_irreducible_unchecked(f),
551                                    mult,
552                                ),
553                            );
554                            None
555                        } else {
556                            Some(f)
557                        }
558                    })
559                    .collect();
560
561                // Is there anything left to factor?
562                if to_factor.is_empty() {
563                    break;
564                }
565
566                // Try to factor what's left using Cantor-Zassenhaus
567                // let h be a random polynomial of degree <n
568                let h = Polynomial::<FS::Set>::from_coeffs(
569                    (0..n).map(|_| prand_elements.next().unwrap()).collect(),
570                );
571                let g = if p == Natural::TWO {
572                    // when char = 2 use h + h^2 + h^4 + ... + h^{2^{kd-1}}
573                    // https://math.stackexchange.com/questions/1636518/how-do-i-apply-the-cantor-zassenhaus-algorithm-to-mathbbf-2
574                    let mut sum = self.poly_ring.zero();
575                    let mut square_powers = h.clone();
576                    let mut square_pow = 0usize;
577                    while Natural::from(square_pow) < &k * Natural::from(d) {
578                        self.poly_ring.add_mut(&mut sum, &square_powers);
579                        square_powers = self.poly_ring.mul(&square_powers, &square_powers);
580                        square_pow += 1;
581                    }
582                    sum
583                } else {
584                    // when char != 2 use h^{(q^d-1)/2}-1 mod f
585                    let poly_mod_f = self.poly_ring.quotient_ring(ddf.polynomial.clone());
586                    let a = (q.nat_pow(&d.into()) - Natural::ONE) / Natural::TWO;
587                    poly_mod_f.add(
588                        &poly_mod_f.nat_pow(&h, &a),
589                        &poly_mod_f.neg(&poly_mod_f.one()),
590                    )
591                };
592                to_factor = to_factor
593                    .into_iter()
594                    .flat_map(|u| {
595                        let gcd = self
596                            .poly_ring
597                            .factorize_monic(&self.poly_ring.subresultant_gcd(u.clone(), g.clone()))
598                            .unwrap()
599                            .monic;
600                        let gcd_deg = self.poly_ring.degree(&gcd).unwrap();
601                        if gcd_deg == 0 || gcd_deg == self.poly_ring.degree(&u).unwrap() {
602                            vec![u]
603                        } else {
604                            vec![self.poly_ring.try_divide(&u, &gcd).unwrap(), gcd]
605                        }
606                    })
607                    .collect();
608            }
609        }
610        fs
611    }
612}
613
614#[cfg(test)]
615mod tests {
616    use super::*;
617    use crate::finite_fields::{modulo::*, quaternary_field::QuaternaryField};
618
619    #[test]
620    fn test_distinct_degree_and_cantor_zassenhaus_factorization_f2() {
621        let x = &Polynomial::<Modulo<2>>::var().into_ergonomic();
622        let p = (x.pow(3) + x.pow(2) + 2) * (x.pow(3) + 2 * x.pow(2) + 1);
623        let p = p.into_verbose();
624
625        assert!(
626            Polynomial::<Modulo<2>>::structure().factorizations().equal(
627                &p.factorize_monic()
628                    .unwrap()
629                    .factorize_squarefree()
630                    .factorize_berlekamps(),
631                &p.factorize_monic()
632                    .unwrap()
633                    .factorize_squarefree()
634                    .factorize_distinct_degree()
635                    .factorize_cantor_zassenhaus()
636            )
637        );
638    }
639
640    #[test]
641    fn test_distinct_degree_and_cantor_zassenhaus_factorization_f3() {
642        let x = &Polynomial::<Modulo<3>>::var().into_ergonomic();
643        let p = (x.pow(3) + x.pow(2) + 2) * (x.pow(3) + 2 * x.pow(2) + 1);
644        let p = p.into_verbose();
645
646        assert!(
647            Polynomial::<Modulo<3>>::structure().factorizations().equal(
648                &p.factorize_monic()
649                    .unwrap()
650                    .factorize_squarefree()
651                    .factorize_berlekamps(),
652                &p.factorize_monic()
653                    .unwrap()
654                    .factorize_squarefree()
655                    .factorize_distinct_degree()
656                    .factorize_cantor_zassenhaus()
657            )
658        );
659    }
660
661    #[test]
662    fn test_factorize_over_f2_example1() {
663        let x = &Polynomial::<Modulo<2>>::var().into_ergonomic();
664        let p = (1 + x.pow(4) + x.pow(5)).pow(12).into_verbose();
665        let ans = Polynomial::<Modulo<2>>::structure()
666            .factorizations()
667            .new_unit_and_powers_unchecked(
668                Polynomial::one(),
669                vec![
670                    ((1 + x + x.pow(2)).into_verbose(), Natural::from(12u32)),
671                    ((1 + x + x.pow(3)).into_verbose(), Natural::from(12u32)),
672                ],
673            );
674
675        let f = Polynomial::<Modulo<2>>::structure().factorize_by_trying_all_factors(p.clone());
676        println!("{} = {}", p, f);
677        assert!(
678            Polynomial::<Modulo<2>>::structure()
679                .factorizations()
680                .equal(&f, &ans)
681        );
682
683        let f = p
684            .factorize_monic()
685            .unwrap()
686            .factorize_squarefree()
687            .factorize_berlekamps();
688        println!("{} = {}", p, f);
689        assert!(
690            Polynomial::<Modulo<2>>::structure()
691                .factorizations()
692                .equal(&f, &ans)
693        );
694    }
695
696    #[test]
697    fn test_factorize_over_f2_example2() {
698        let x = &Polynomial::<Modulo<2>>::var().into_ergonomic();
699        let p =
700            ((1 + x.pow(4) + x.pow(7)).pow(6) * (1 + x.pow(6) + x.pow(7)).pow(4)).into_verbose();
701        let ans = Polynomial::<Modulo<2>>::structure()
702            .factorizations()
703            .new_unit_and_powers_unchecked(
704                Polynomial::one(),
705                vec![
706                    (
707                        (1 + x.pow(4) + x.pow(7)).into_verbose(),
708                        Natural::from(6u32),
709                    ),
710                    (
711                        (1 + x.pow(6) + x.pow(7)).into_verbose(),
712                        Natural::from(4u32),
713                    ),
714                ],
715            );
716
717        let f = Polynomial::<Modulo<2>>::structure().factorize_by_trying_all_factors(p.clone());
718        println!("{} = {}", p, f);
719        assert!(
720            Polynomial::<Modulo<2>>::structure()
721                .factorizations()
722                .equal(&f, &ans)
723        );
724
725        let f = p
726            .factorize_monic()
727            .unwrap()
728            .factorize_squarefree()
729            .factorize_berlekamps();
730        println!("{} = {}", p, f);
731        assert!(
732            Polynomial::<Modulo<2>>::structure()
733                .factorizations()
734                .equal(&f, &ans)
735        );
736    }
737
738    #[test]
739    fn test_factorize_over_f5_example1() {
740        let x = &Polynomial::<Modulo<5>>::var().into_ergonomic();
741        let p = (1 + x.pow(4)).pow(5).into_verbose();
742        let fs = Polynomial::<Modulo<5>>::structure().into_factorizations();
743        let ans = fs.new_unit_and_powers_unchecked(
744            Polynomial::one(),
745            vec![
746                ((2 + x.pow(2)).into_verbose(), Natural::from(5u8)),
747                ((3 + x.pow(2)).into_verbose(), Natural::from(5u8)),
748            ],
749        );
750
751        let f = Polynomial::<Modulo<5>>::structure().factorize_by_trying_all_factors(p.clone());
752        println!("{} = {}", p, f);
753        assert!(fs.equal(&f, &ans));
754
755        let f = p
756            .factorize_monic()
757            .unwrap()
758            .factorize_squarefree()
759            .factorize_berlekamps();
760        println!("{} = {}", p, f);
761        assert!(fs.equal(&f, &ans));
762    }
763
764    #[test]
765    fn test_factorize_over_f5_example2() {
766        let x = &Polynomial::<Modulo<5>>::var().into_ergonomic();
767        let p = (3 + 2 * x.pow(2) + x.pow(4) + x.pow(6)).into_verbose();
768        let ans = Polynomial::<Modulo<5>>::structure()
769            .factorizations()
770            .new_unit_and_powers_unchecked(
771                Polynomial::one(),
772                vec![((2 + x.pow(2)).into_verbose(), Natural::from(3u8))],
773            );
774
775        let f = Polynomial::<Modulo<5>>::structure().factorize_by_trying_all_factors(p.clone());
776        println!("{} = {}", p, f);
777        assert!(
778            Polynomial::<Modulo<5>>::structure()
779                .factorizations()
780                .equal(&f, &ans)
781        );
782
783        let f = p
784            .factorize_monic()
785            .unwrap()
786            .factorize_squarefree()
787            .factorize_berlekamps();
788        println!("{} = {}", p, f);
789        assert!(
790            Polynomial::<Modulo<5>>::structure()
791                .factorizations()
792                .equal(&f, &ans)
793        );
794    }
795
796    #[test]
797    fn test_factorize_over_f31_example1() {
798        let x = &Polynomial::<Modulo<31>>::var().into_ergonomic();
799        let p = (1 + x.pow(27) + 8 * x.pow(30)).into_verbose();
800        let ans = Polynomial::<Modulo<31>>::structure()
801            .factorizations()
802            .new_unit_and_powers_unchecked(
803                Polynomial::constant(Modulo::from_int(Integer::from(8))),
804                vec![
805                    ((12 + x.pow(3)).into_verbose(), Natural::from(1u32)),
806                    (
807                        (25 + 27 * x + 3 * x.pow(2) + 3 * x.pow(3) + 29 * x.pow(4) + x.pow(5))
808                            .into_verbose()
809                            .clone(),
810                        Natural::from(1u32),
811                    ),
812                    (
813                        (1 + 24 * x + 3 * x.pow(2) + 15 * x.pow(3) + 12 * x.pow(4) + x.pow(5))
814                            .into_verbose()
815                            .clone(),
816                        Natural::from(1u32),
817                    ),
818                    (
819                        (21 + 12 * x.pow(3) + 22 * x.pow(6) + 4 * x.pow(9) + x.pow(12))
820                            .into_verbose()
821                            .clone(),
822                        Natural::from(1u32),
823                    ),
824                    (
825                        (5 + 11 * x + 3 * x.pow(2) + 13 * x.pow(3) + 21 * x.pow(4) + x.pow(5))
826                            .into_verbose()
827                            .clone(),
828                        Natural::from(1u32),
829                    ),
830                ],
831            );
832
833        println!("{:?}", p.factorize_monic());
834
835        let f = p
836            .factorize_monic()
837            .unwrap()
838            .factorize_squarefree()
839            .factorize_berlekamps();
840        println!("{} = {}", p, f);
841        assert!(
842            Polynomial::<Modulo<31>>::structure()
843                .factorizations()
844                .equal(&f, &ans)
845        );
846    }
847
848    #[test]
849    fn test_factorize_over_f4_example1() {
850        let x = &Polynomial::<QuaternaryField>::var().into_ergonomic();
851
852        let a = x - Polynomial::constant(QuaternaryField::One).into_ergonomic();
853        let b = x - Polynomial::constant(QuaternaryField::Alpha).into_ergonomic();
854        let c = x - Polynomial::constant(QuaternaryField::Beta).into_ergonomic();
855        let p = (1 - x.pow(3)).pow(48).into_verbose();
856        let ans = Polynomial::<QuaternaryField>::structure()
857            .factorizations()
858            .new_unit_and_powers_unchecked(
859                Polynomial::one(),
860                vec![
861                    (a.into_verbose(), Natural::from(48u32)),
862                    (b.into_verbose(), Natural::from(48u32)),
863                    (c.into_verbose(), Natural::from(48u32)),
864                ],
865            );
866
867        let f = QuaternaryField::structure()
868            .polynomials()
869            .factorize_by_trying_all_factors(p.clone());
870        println!("{} = {}", p, f);
871        assert!(
872            Polynomial::<QuaternaryField>::structure()
873                .factorizations()
874                .equal(&f, &ans)
875        );
876
877        let f = p
878            .factorize_monic()
879            .unwrap()
880            .factorize_squarefree()
881            .factorize_berlekamps();
882        println!("{} = {}", p, f);
883        assert!(
884            Polynomial::<QuaternaryField>::structure()
885                .factorizations()
886                .equal(&f, &ans)
887        );
888    }
889}
algebraeon_rings/finite_fields/polynomial.rs

algebraeon_rings/finite_fields/
polynomial.rs
