From 9d51f7c457671bfcbab9a1d62d416e1a83e6ad8a Mon Sep 17 00:00:00 2001 From: Rasmus Munk Larsen Date: Fri, 1 Apr 2016 16:48:38 -0700 Subject: Add rcond method to LDLT. --- test/cholesky.cpp | 17 +++++++++++++++++ 1 file changed, 17 insertions(+) (limited to 'test/cholesky.cpp') diff --git a/test/cholesky.cpp b/test/cholesky.cpp index 8a21cdbd5..148a0b388 100644 --- a/test/cholesky.cpp +++ b/test/cholesky.cpp @@ -160,6 +160,15 @@ template void cholesky(const MatrixType& m) matX = ldltlo.solve(matB); VERIFY_IS_APPROX(symm * matX, matB); + // Verify that the estimated condition number is within a factor of 10 of the + // truth. + const MatrixType symmLo_inverse = ldltlo.solve(MatrixType::Identity(rows,cols)); + RealScalar rcond = (RealScalar(1) / matrix_l1_norm(symmLo)) / + matrix_l1_norm(symmLo_inverse); + RealScalar rcond_est = ldltlo.rcond(); + VERIFY(rcond_est > rcond / 10 && rcond_est < rcond * 10); + + LDLT ldltup(symmUp); VERIFY_IS_APPROX(symm, ldltup.reconstructedMatrix()); vecX = ldltup.solve(vecB); @@ -167,6 +176,14 @@ template void cholesky(const MatrixType& m) matX = ldltup.solve(matB); VERIFY_IS_APPROX(symm * matX, matB); + // Verify that the estimated condition number is within a factor of 10 of the + // truth. + const MatrixType symmUp_inverse = ldltup.solve(MatrixType::Identity(rows,cols)); + rcond = (RealScalar(1) / matrix_l1_norm(symmUp)) / + matrix_l1_norm(symmUp_inverse); + rcond_est = ldltup.rcond(); + VERIFY(rcond_est > rcond / 10 && rcond_est < rcond * 10); + VERIFY_IS_APPROX(MatrixType(ldltlo.matrixL().transpose().conjugate()), MatrixType(ldltlo.matrixU())); VERIFY_IS_APPROX(MatrixType(ldltlo.matrixU().transpose().conjugate()), MatrixType(ldltlo.matrixL())); VERIFY_IS_APPROX(MatrixType(ldltup.matrixL().transpose().conjugate()), MatrixType(ldltup.matrixU())); -- cgit v1.2.3