1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126

use crate::Matrix;
use std::error::Error;

impl<T> Matrix<T>
where T: std::ops::Mul<Output=T>
+ std::ops::Sub<Output=T>
+ Clone
+ std::cmp::PartialOrd
+ std::ops::Div<Output=T>
+ std::ops::Add<Output=T> {
    
    pub fn inverse(&self, zero: T, one: T) -> Result<Matrix<T>, Box<dyn Error>> {

        let (rows, columns) = self.dimensions()?;

        if rows == columns {

            let determinant = self.determinant()?;

            if determinant != zero {

                let mut a = self.clone();

                let identity = Matrix::identity(rows, zero.clone(), one.clone());

                let mut inverse = identity.clone();

                for c in 0..columns {

                    if self.0[c][c] != one {
                    
                        if self.0[c][c] != zero {

                            inverse.0[c] = inverse.0[c]
                                .iter()
                                .map(|x| x.clone() / a.0[c][c].clone())
                                .collect();

                            a.0[c] = a.0[c]
                                .iter()
                                .map(|x| x.clone() / a.0[c][c].clone())
                                .collect();

                        } else {

                            inverse.0[c] = inverse.0[c]
                                .iter()
                                .map(|x| x.clone() + one.clone())
                                .collect();

                            a.0[c] = a.0[c]
                                .iter()
                                .map(|x| x.clone() + one.clone())
                                .collect();

                        }
                    
                    }

                    for r in 0..rows {

                        if c != r {

                            inverse.0[r] = inverse.0[r]
                                .iter()
                                .enumerate()
                                .map(|(i,x)| x.clone() - (a.0[r][c].clone() * inverse.0[c][i].clone()))
                                .collect();
                            
                            a.0[r] = a.0[r]
                                .iter()
                                .enumerate()
                                .map(|(i,x)| x.clone() - (a.0[r][c].clone() * a.0[c][i].clone()))
                                .collect();

                        }

                    }

                }

                if (self.clone() * inverse.clone())? == identity {
                    
                    Ok(inverse)

                } else {
                    
                    Err("Cannot find inverse!")?

                }
                

            } else {

                Err("Matrix is not invertible!")?

            }

        } else {
            
            Err("Non Square matrix!")?

        }
        
        

    }

}

#[cfg(test)]
mod tests {
    
    use super::*;
    #[test]
    fn test_matrix_inverse() {

        let a = Matrix(vec![vec![1,2], vec![1,3]]);

        let a_inverse = Matrix(vec![vec![3,-2], vec![-1,1]]);

        assert_eq!(a.inverse(0,1).unwrap(), a_inverse);
         
    }

}