math_rs/matrix/float/
f32.rs

1use crate::result::{MathError, Result};
2
3use crate::matrix::traits::{Identity, Matrix, Parseable};
4
5#[derive(Debug, Clone)]
6pub struct MatrixF32 {
7    pub content: Vec<Vec<f32>>,
8    tolerance: f32,
9}
10
11#[macro_export]
12macro_rules! matrix_f32 {
13    ($expression:tt, $tol:expr) => {
14        MatrixF32::parse($expression, $tol)
15    };
16}
17pub use matrix_f32;
18
19impl Matrix for MatrixF32 {
20    type T = f32;
21    fn columns(&self) -> usize {
22        self.content
23            .iter()
24            .next()
25            .map(|row| row.len())
26            .expect("There is no row") // TODO: Arrange this in a better way
27    }
28
29    fn rows(&self) -> usize {
30        self.content.len()
31    }
32
33    fn is_square(&self) -> bool {
34        self.columns() == self.rows()
35    }
36
37    fn is_symmetric(&self) -> bool {
38        if !self.is_square() {
39            return false;
40        }
41        for i in 0..self.rows() {
42            for j in 0..self.columns() {
43                if self.get(i, j).unwrap() != self.get(j, i).unwrap() {
44                    return false;
45                }
46            }
47        }
48        true
49    }
50
51    fn get(&self, row: usize, column: usize) -> Result<&Self::T> {
52        if row > self.rows() || column > self.columns() {
53            return Err(MathError::MatrixError(format!("Cannot get element from position ({row},{column}) if the matrix has {}x{} dimensions!", self.rows(), self.columns())));
54        }
55        let Some(matrix_row) = self.content.get(row) else {
56            return Err(MathError::MatrixError(format!("Cannot get row {row} if the matrix has {}x{} dimensions!", self.rows(), self.columns())));
57        };
58        matrix_row.get(column).ok_or(MathError::MatrixError(format!(
59            "Cannot get element from position ({row},{column}) if the matrix has {}x{} dimensions!",
60            self.rows(),
61            self.columns()
62        )))
63    }
64
65    fn get_mut(&mut self, row: usize, column: usize) -> Result<&mut Self::T> {
66        if row > self.rows() || column > self.columns() {
67            return Err(MathError::MatrixError(format!("Cannot get element from position ({row},{column}) if the matrix has {}x{} dimensions!", self.rows(), self.columns())));
68        }
69        let rows = self.rows();
70        let columns = self.columns();
71        let Some(matrix_row) = self.content.get_mut(row) else {
72            return Err(MathError::MatrixError(format!("Cannot get row {row} if the matrix has {}x{} dimensions!", rows, columns)));
73        };
74        matrix_row
75            .get_mut(column)
76            .ok_or(MathError::MatrixError(format!(
77            "Cannot get element from position ({row},{column}) if the matrix has {}x{} dimensions!",
78            rows,
79            columns
80        )))
81    }
82
83    fn set(&mut self, row: usize, column: usize, value: Self::T) -> Result<()> {
84        *self.get_mut(row, column)? = value;
85        Ok(())
86    }
87
88    fn swap_rows(&mut self, row1: usize, row2: usize) -> Result<()> {
89        if row1 > self.rows() || row2 > self.rows() {
90            return Err(MathError::MatrixError(format!(
91                "Cannot swap rows {row1} and {row2} if the matrix has {} rows!",
92                self.rows()
93            )));
94        }
95        self.content.swap(row1, row2);
96        Ok(())
97    }
98
99    fn transpose(&self) -> Self {
100        let mut new_matrix = self.clone(); // Clone as we want a new matrix
101        for i in 0..self.rows() {
102            for j in 0..self.columns() {
103                new_matrix
104                    .set(j, i, self.get(i, j).unwrap().clone())
105                    .unwrap();
106            }
107        }
108        new_matrix
109    }
110
111    fn gaussian_triangulation(&self) -> Result<Self> {
112        let mut new_matrix = self.clone();
113        for i in 0..self.rows() {
114            let mut max_row = i;
115            for j in i + 1..self.rows() {
116                if new_matrix.get(j, i).unwrap().abs() > new_matrix.get(max_row, i).unwrap().abs() {
117                    max_row = j;
118                }
119            }
120            new_matrix.swap_rows(i, max_row)?;
121            if i < self.rows() - 1 && new_matrix.get(i, i).unwrap().abs() <= new_matrix.tolerance()
122            {
123                return Err(MathError::MatrixError(
124                    "Cannot perform gaussian elimination on a matrix with zero pivot".to_string(),
125                ));
126            }
127            for j in i + 1..self.rows() {
128                let factor = new_matrix.get(j, i).unwrap() / new_matrix.get(i, i).unwrap();
129                for k in i..self.columns() {
130                    let new_value =
131                        new_matrix.get(j, k).unwrap() - factor * new_matrix.get(i, k).unwrap();
132                    new_matrix.set(j, k, new_value).unwrap();
133                }
134            }
135        }
136        Ok(new_matrix)
137    }
138
139    fn lu_decomposition(&self) -> Result<(Self, Self)> {
140        if !self.is_square() {
141            return Err(MathError::MatrixError(
142                "Cannot perform LU decomposition on a non-square matrix".to_string(),
143            ));
144        }
145        let mut lower = Self::id(self.rows(), self.tolerance());
146        let mut upper = self.clone();
147        for i in 0..self.rows() {
148            for j in 0..self.columns() {
149                if j < i {
150                    let factor = upper.get(i, j).unwrap() / upper.get(j, j).unwrap();
151                    lower.set(i, j, factor).unwrap();
152                    for k in j..self.columns() {
153                        let new_value =
154                            upper.get(i, k).unwrap() - factor * upper.get(j, k).unwrap();
155                        upper.set(i, k, new_value).unwrap();
156                    }
157                }
158            }
159        }
160        Ok((lower, upper))
161    }
162
163    fn determinant_using_gauss(&self) -> Option<Self::T> {
164        let gaussian_elimination_result = self.gaussian_triangulation().ok()?;
165        let mut mult = Self::T::id(0, 0.0);
166        for i in 0..gaussian_elimination_result.rows() {
167            for j in 0..gaussian_elimination_result.columns() {
168                if i == j {
169                    let value = gaussian_elimination_result.get(i, j).unwrap();
170                    mult = mult * *value;
171                }
172            }
173        }
174        Some(mult)
175    }
176
177    fn determinant_using_lu(&self) -> Option<Self::T> {
178        let (_, upper) = self.lu_decomposition().ok()?;
179        let mut mult = Self::T::id(0, 0.0);
180        for i in 0..self.rows() {
181            for j in 0..self.columns() {
182                if i == j {
183                    let upper_value = upper.get(i, j).unwrap();
184                    mult = mult * *upper_value;
185                }
186            }
187        }
188        Some(mult)
189    }
190
191    fn cholesky_decomposition(&self) -> Result<Self> {
192        todo!("To be done");
193    }
194}
195
196impl MatrixF32 {
197    pub fn new(content: Vec<Vec<f32>>, tolerance: f32) -> Result<Self> {
198        let mut content_iter = content.iter();
199        if let Some(length) = content_iter.next().map(|row| row.len()) {
200            while let Some(next_length) = content_iter.next().map(|row| row.len()) {
201                if length != next_length {
202                    return Err(MathError::MatrixError(
203                        "Cannot build matrix from rows with different lenght".to_string(),
204                    ));
205                }
206            }
207        }
208        Ok(Self { content, tolerance })
209    }
210
211    pub fn tolerance(&self) -> f32 {
212        self.tolerance
213    }
214}
215
216impl TryFrom<&str> for MatrixF32 {
217    /// Performs the conversion from a string to the matrix, with default tolerance 1e-4
218    fn try_from(value: &str) -> Result<Self> {
219        MatrixF32::parse(value, 1e-4)
220    }
221
222    type Error = MathError;
223}
224
225#[cfg(test)]
226mod test {
227    use crate::matrix::traits::{Matrix, Parseable, CheckedAdd};
228
229    use super::{matrix_f32, MatrixF32};
230    use pretty_assertions;
231
232    const TOLERANCE: f32 = 1e-10;
233
234    #[test]
235    fn get_1() {
236        let matrix = MatrixF32::new(vec![vec![1f32, 2f32], vec![3f32, 4f32]], TOLERANCE).unwrap();
237        let item = matrix.get(0, 1).unwrap();
238        pretty_assertions::assert_eq!(item.clone(), 2f32)
239    }
240
241    #[test]
242    fn get_2() {
243        let matrix = MatrixF32::new(
244            vec![
245                vec![1.1, 2.2, 3.3],
246                vec![4.4, 5.5, 6.6],
247                vec![7.7, 8.8, 9.9],
248            ],
249            TOLERANCE,
250        )
251        .unwrap();
252        vec![((0, 0), 1.1), ((0, 2), 3.3), ((1, 1), 5.5), ((2, 0), 7.7)]
253            .into_iter()
254            .for_each(|item| {
255                pretty_assertions::assert_eq!(
256                    matrix.get(item.0 .0, item.0 .1).unwrap().clone(),
257                    item.1
258                )
259            })
260    }
261
262    #[test]
263    fn equal_matrices() {
264        let matrix_a = MatrixF32::new(vec![vec![1f32, 1f32], vec![2f32, 2f32]], TOLERANCE)
265            .expect("Should've been able to built this matrix");
266        let matrix_b = MatrixF32::new(vec![vec![1f32, 1f32], vec![2f32, 2f32]], TOLERANCE)
267            .expect("Should've been able to built this matrix");
268        pretty_assertions::assert_eq!(matrix_a, matrix_b)
269    }
270
271    #[test]
272    fn create_matrix_from_macro_1() {
273        let matrix = matrix_f32!("{{1,2},{3,4}}", TOLERANCE)
274            .expect("Should have been able to build matrix from macro");
275        println!("{matrix}");
276        let other = matrix_f32!("{{1,2},{3,4}}", TOLERANCE).expect("asdf");
277        pretty_assertions::assert_eq!(matrix, other)
278    }
279
280    #[test]
281    fn sum_with_macro() {
282        let result = (matrix_f32!("{{1,1},{1,1}}", TOLERANCE)
283            .expect("asdf")
284            .checked_add(&matrix_f32!("{{2,2},{2,2}}", TOLERANCE).expect("asdf")))
285        .expect("asdf");
286        println!("{result}")
287    }
288
289    #[test]
290    fn gaussian_elimination_1() {
291        let matrix = matrix_f32!("{{1,2,3},{4,5,6},{7,8,9}}", TOLERANCE).expect("asdf");
292        let gauss = matrix.gaussian_triangulation().expect("asdf");
293        pretty_assertions::assert_eq!(
294            gauss,
295            matrix_f32!(
296                "{{+7.0000000000000, +8.0000000000000, +9.0000000000000},
297                {+0.0000000000000, +0.8571428060532, +1.7142856121063 },
298                {+0.0000000000000, +0.0000000000000, +0.0000000000000}}",
299                TOLERANCE
300            )
301            .expect("asdf")
302        );
303    }
304
305    #[test]
306    fn gaussian_elimination_2() {
307        let matrix =
308            matrix_f32!("{{1,2,1,2,1},{-1,2,-3,2,1},{0,1,-3,2,1}}", TOLERANCE).expect("asdf");
309        let gauss = matrix.gaussian_triangulation().expect("asdf");
310        pretty_assertions::assert_eq!(
311            gauss,
312            matrix_f32!("{{1,2,1,2,1},{0,4,-2,4,2},{0,0,-2.5,1,0.5}}", TOLERANCE).expect("asdf")
313        );
314    }
315
316    #[test]
317    fn lu_decomposition_1() {
318        let matrix = matrix_f32!("{{1,2,3},{4,5,6},{7,8,9}}", TOLERANCE).expect("asdf");
319        let (l, u) = matrix.lu_decomposition().expect("asdf");
320        pretty_assertions::assert_eq!(
321            l,
322            matrix_f32!("{{1,0,0},{4,1,0},{7,2,1}}", TOLERANCE).expect("asdf")
323        );
324        pretty_assertions::assert_eq!(
325            u,
326            matrix_f32!("{{1,2,3},{0,-3,-6},{0,0,0}}", TOLERANCE).expect("asdf")
327        );
328    }
329
330    #[test]
331    fn determinant_1() {
332        let matrix = matrix_f32!("{{1,2,3},{4,5,6},{7,8,9}}", TOLERANCE).expect("asdf");
333        let determinant = matrix.determinant_using_gauss().expect("asdf");
334        pretty_assertions::assert_eq!(determinant, 0f32);
335    }
336
337    #[test]
338    fn determinant_using_lu_1() {
339        let matrix = matrix_f32!("{{1,2,3},{4,5,6},{7,8,1}}", TOLERANCE).expect("asdf");
340        let (lower, upper) = matrix.lu_decomposition().expect("asdf");
341        println!("{lower}");
342        println!("{upper}");
343        let determinant = matrix.determinant_using_lu().expect("asdf");
344        pretty_assertions::assert_eq!(determinant, 24f32);
345    }
346}