use super::Matrix;
use crate::num::MinusOne;
use crate::num::One;
use std::fmt::Debug;
use std::ops::{Add, Div, Mul, Sub, SubAssign};
impl<K> Matrix<K>
where
K: Debug
+ Copy
+ One
+ MinusOne
+ Default
+ PartialEq
+ Add<Output = K>
+ Sub<Output = K>
+ SubAssign<K>
+ Mul<Output = K>
+ Div<Output = K>,
{
fn find_pivot(&self, start_row: usize, col: usize) -> Option<usize> {
for row in start_row..self.rows() {
if self.elements[row][col] != K::default() {
return Some(row);
}
}
None
}
fn div_row(&mut self, row: usize, scalar: K) {
self.elements[row]
.iter_mut()
.for_each(|element| *element = *element / scalar);
}
pub(crate) fn gaussian_elimination(
&self,
mut det: Option<&mut K>,
mut identity_mat: Option<&mut Matrix<K>>,
) -> Matrix<K> {
if self.elements.is_empty() || self.elements[0].is_empty() {
return self.clone();
}
let mut result = self.clone();
let mut pivot_row = 0;
let mut sign = K::one();
for col in 0..result.cols() {
if let Some(pivot) = result.find_pivot(pivot_row, col) {
if pivot != pivot_row {
result.elements.swap(pivot_row, pivot);
if let Some(ref mut identity) = identity_mat {
identity.elements.swap(pivot_row, pivot);
}
sign = sign * K::minus_one();
}
let pivot_value = result.elements[pivot_row][col];
if let Some(ref mut determinant) = det {
**determinant = **determinant * pivot_value;
}
result.div_row(pivot_row, pivot_value);
if let Some(ref mut identity) = identity_mat {
identity.div_row(pivot_row, pivot_value);
}
for row in 0..result.rows() {
if row == pivot_row {
continue;
}
let factor = result.elements[row][col];
if factor == K::default() {
continue;
}
for column in col..result.cols() {
let pivot_value = result.elements[pivot_row][column];
result.elements[row][column] -= factor * pivot_value;
}
if let Some(ref mut identity) = identity_mat {
for column in 0..identity.cols() {
let pivot_identity_value = identity.elements[pivot_row][column];
identity.elements[row][column] -= factor * pivot_identity_value;
}
}
}
pivot_row += 1;
} else {
if let Some(ref mut determinant) = det {
**determinant = K::default();
}
}
}
if let Some(ref mut determinant) = det {
**determinant = **determinant * sign;
}
result
}
}