algebraeon_rings/matrix/
polynomial.rs

1use super::*;
2
3impl<FS: FieldSignature, FSB: BorrowedStructure<FS>> MatrixStructure<FS, FSB> {
4    pub fn presentation_matrix(
5        &self,
6        m: Matrix<FS::Set>,
7    ) -> Result<Matrix<Polynomial<FS::Set>>, MatOppErr> {
8        let n = m.rows();
9        if n == m.cols() {
10            let poly_ring = self.ring().polynomials();
11            let poly_mat_struct = MatrixStructure::new(poly_ring.clone());
12            Ok(poly_mat_struct
13                .add(
14                    &m.apply_map(|x| Polynomial::constant(x.clone())),
15                    &poly_mat_struct.neg(
16                        poly_mat_struct.diag(&(0..n).map(|_i| poly_ring.var()).collect::<Vec<_>>()),
17                    ),
18                )
19                .unwrap())
20        } else {
21            Err(MatOppErr::NotSquare)
22        }
23    }
24
25    pub fn minimal_polynomial(&self, m: Matrix<FS::Set>) -> Result<Polynomial<FS::Set>, MatOppErr> {
26        match self.presentation_matrix(m) {
27            Ok(pres_mat) => {
28                let poly_ring = self.ring().polynomials();
29                let poly_mat_struct = MatrixStructure::new(poly_ring.clone());
30                let (_u, s, _v, k) = poly_mat_struct.smith_algorithm(pres_mat);
31                debug_assert!(k > 0); //cant be all zero because we are taking SNF of a non-zero matrix
32                Ok(s.at(k - 1, k - 1).unwrap().clone())
33            }
34            Err(MatOppErr::NotSquare) => Err(MatOppErr::NotSquare),
35            Err(_) => panic!(),
36        }
37    }
38
39    pub fn characteristic_polynomial(
40        &self,
41        m: Matrix<FS::Set>,
42    ) -> Result<Polynomial<FS::Set>, MatOppErr> {
43        match self.presentation_matrix(m) {
44            Ok(pres_mat) => {
45                let poly_ring = self.ring().polynomials();
46                let poly_mat_struct = MatrixStructure::new(poly_ring.clone());
47                let (_u, s, _v, k) = poly_mat_struct.smith_algorithm(pres_mat);
48                debug_assert!(k > 0); //cant be all zero because we are taking SNF of a non-zero matrix
49                let mut char_poly = poly_ring.one();
50                for i in 0..k {
51                    poly_ring.mul_mut(&mut char_poly, s.at(i, i).unwrap());
52                }
53                Ok(char_poly)
54            }
55            Err(MatOppErr::NotSquare) => Err(MatOppErr::NotSquare),
56            Err(_) => panic!(),
57        }
58    }
59}
60
61impl<F: MetaType> Matrix<F>
62where
63    F::Signature: FieldSignature,
64{
65    pub fn presentation_matrix(&self) -> Result<Matrix<Polynomial<F>>, MatOppErr> {
66        Self::structure().presentation_matrix(self.clone())
67    }
68
69    pub fn minimal_polynomial(&self) -> Result<Polynomial<F>, MatOppErr> {
70        Self::structure().minimal_polynomial(self.clone())
71    }
72
73    pub fn characteristic_polynomial(&self) -> Result<Polynomial<F>, MatOppErr> {
74        Self::structure().characteristic_polynomial(self.clone())
75    }
76}
77
78#[cfg(test)]
79mod tests {
80    use super::*;
81
82    #[test]
83    fn min_and_char_polys() {
84        {
85            let a = Matrix::<Integer>::from_rows(vec![
86                vec![Integer::from(0), Integer::from(4), Integer::from(4)],
87                vec![Integer::from(1), Integer::from(4), Integer::from(16)],
88            ]);
89            let (_u, s, _v, _k) = a.smith_algorithm();
90
91            assert_eq!(
92                s,
93                Matrix::from_rows(vec![
94                    vec![Integer::from(1), Integer::from(0), Integer::from(0),],
95                    vec![Integer::from(0), Integer::from(4), Integer::from(0)],
96                ])
97            );
98        }
99
100        {
101            let a = Matrix::<Rational>::from_rows(vec![
102                vec![Rational::from(0), Rational::from(0), Rational::from(0)],
103                vec![Rational::from(0), Rational::from(0), Rational::from(1)],
104                vec![Rational::from(0), Rational::from(0), Rational::from(0)],
105            ]);
106            let min_p = a.clone().minimal_polynomial().unwrap();
107            let char_p = a.clone().characteristic_polynomial().unwrap();
108            assert_eq!(
109                &min_p,
110                &Polynomial::from_coeffs(vec![
111                    Rational::from(0),
112                    Rational::from(0),
113                    Rational::from(1)
114                ])
115            );
116            assert_eq!(
117                &char_p,
118                &Polynomial::from_coeffs(vec![
119                    Rational::from(0),
120                    Rational::from(0),
121                    Rational::from(0),
122                    Rational::from(1)
123                ])
124            );
125        }
126
127        {
128            let a = Matrix::<Rational>::from_rows(vec![
129                vec![
130                    Rational::from(4),
131                    Rational::from(0),
132                    Rational::from(0),
133                    Rational::from(0),
134                ],
135                vec![
136                    Rational::from(0),
137                    Rational::from(4),
138                    Rational::from(0),
139                    Rational::from(0),
140                ],
141                vec![
142                    Rational::from(0),
143                    Rational::from(1),
144                    Rational::from(4),
145                    Rational::from(0),
146                ],
147                vec![
148                    Rational::from(0),
149                    Rational::from(0),
150                    Rational::from(1),
151                    Rational::from(4),
152                ],
153            ]);
154            let min_p = a.clone().minimal_polynomial().unwrap();
155            let char_p = a.clone().characteristic_polynomial().unwrap();
156            assert_eq!(
157                &min_p,
158                &Polynomial::from_coeffs(vec![Rational::from(-4), Rational::from(1),])
159                    .nat_pow(&Natural::from(3u8))
160            );
161            assert_eq!(
162                &char_p,
163                &Polynomial::from_coeffs(vec![Rational::from(-4), Rational::from(1),])
164                    .nat_pow(&Natural::from(4u8))
165            );
166        }
167    }
168}
algebraeon_rings/matrix/polynomial.rs

algebraeon_rings/matrix/
polynomial.rs
