Add determinant method to Cholesky

This commit is contained in:
Jim Turner 2021-04-07 22:43:27 -04:00
parent f8173bd668
commit 50fed194f4
2 changed files with 28 additions and 0 deletions

View File

@ -138,6 +138,16 @@ where
self.solve_mut(&mut res);
res
}
/// Computes the determinant of the decomposed matrix.
pub fn determinant(&self) -> N {
let dim = self.chol.nrows();
let mut prod_diag = N::one();
for i in 0..dim {
prod_diag *= unsafe { *self.chol.get_unchecked((i, i)) };
}
prod_diag * prod_diag
}
}
impl<N: ComplexField, D: Dim> Cholesky<N, D>

View File

@ -79,6 +79,24 @@ macro_rules! gen_tests(
prop_assert!(id1.is_identity(1.0e-7) && id2.is_identity(1.0e-7));
}
#[test]
fn cholesky_determinant(n in PROPTEST_MATRIX_DIM) {
let m = RandomSDP::new(Dynamic::new(n), || random::<$scalar>().0).unwrap();
let lu_det = m.clone().lu().determinant();
let chol_det = m.cholesky().unwrap().determinant();
prop_assert!(relative_eq!(lu_det, chol_det, epsilon = 1.0e-7));
}
#[test]
fn cholesky_determinant_static(_n in PROPTEST_MATRIX_DIM) {
let m = RandomSDP::new(U4, || random::<$scalar>().0).unwrap();
let lu_det = m.clone().lu().determinant();
let chol_det = m.cholesky().unwrap().determinant();
prop_assert!(relative_eq!(lu_det, chol_det, epsilon = 1.0e-7));
}
#[test]
fn cholesky_rank_one_update(_n in PROPTEST_MATRIX_DIM) {
let mut m = RandomSDP::new(U4, || random::<$scalar>().0).unwrap();