diff options
author | Yi Kong <yikong@google.com> | 2022-02-25 16:32:14 +0800 |
---|---|---|
committer | Yi Kong <yikong@google.com> | 2022-02-25 15:08:55 +0000 |
commit | 2aab794c004027d008d6b0b64165bf1961d5d2bb (patch) | |
tree | 83bb8f19c67bcafdb2ca4a98414af1b17392ec36 /test/cholesky.cpp | |
parent | ca5aa72016f062fd0712bcb86370478de332bca3 (diff) | |
download | eigen-2aab794c004027d008d6b0b64165bf1961d5d2bb.tar.gz |
Upgrade eigen to 3.4.0
Steps:
* Removed common files between Android copy and the matching upstream copy
* Obtained latest upstream tarball (see README.version)
* Extracted over the directory
Bug: 148287349
Test: presubmit
Change-Id: Iee2744719075fdf000b315e973645923da766111
Diffstat (limited to 'test/cholesky.cpp')
-rw-r--r-- | test/cholesky.cpp | 87 |
1 files changed, 55 insertions, 32 deletions
diff --git a/test/cholesky.cpp b/test/cholesky.cpp index 8ad5ac639..0b1a7b45b 100644 --- a/test/cholesky.cpp +++ b/test/cholesky.cpp @@ -7,18 +7,16 @@ // Public License v. 2.0. If a copy of the MPL was not distributed // with this file, You can obtain one at http://mozilla.org/MPL/2.0/. -#ifndef EIGEN_NO_ASSERTION_CHECKING -#define EIGEN_NO_ASSERTION_CHECKING -#endif - #define TEST_ENABLE_TEMPORARY_TRACKING #include "main.h" #include <Eigen/Cholesky> #include <Eigen/QR> +#include "solverbase.h" template<typename MatrixType, int UpLo> typename MatrixType::RealScalar matrix_l1_norm(const MatrixType& m) { + if(m.cols()==0) return typename MatrixType::RealScalar(0); MatrixType symm = m.template selfadjointView<UpLo>(); return symm.cwiseAbs().colwise().sum().maxCoeff(); } @@ -57,7 +55,6 @@ template<typename MatrixType,template <typename,int> class CholType> void test_c template<typename MatrixType> void cholesky(const MatrixType& m) { - typedef typename MatrixType::Index Index; /* this test covers the following files: LLT.h LDLT.h */ @@ -81,15 +78,17 @@ template<typename MatrixType> void cholesky(const MatrixType& m) } { + STATIC_CHECK(( internal::is_same<typename LLT<MatrixType,Lower>::StorageIndex,int>::value )); + STATIC_CHECK(( internal::is_same<typename LLT<MatrixType,Upper>::StorageIndex,int>::value )); + SquareMatrixType symmUp = symm.template triangularView<Upper>(); SquareMatrixType symmLo = symm.template triangularView<Lower>(); LLT<SquareMatrixType,Lower> chollo(symmLo); VERIFY_IS_APPROX(symm, chollo.reconstructedMatrix()); - vecX = chollo.solve(vecB); - VERIFY_IS_APPROX(symm * vecX, vecB); - matX = chollo.solve(matB); - VERIFY_IS_APPROX(symm * matX, matB); + + check_solverbase<VectorType, VectorType>(symm, chollo, rows, rows, 1); + check_solverbase<MatrixType, MatrixType>(symm, chollo, rows, cols, rows); const MatrixType symmLo_inverse = chollo.solve(MatrixType::Identity(rows,cols)); RealScalar rcond = (RealScalar(1) / matrix_l1_norm<MatrixType, Lower>(symmLo)) / @@ -97,7 +96,7 @@ template<typename MatrixType> void cholesky(const MatrixType& m) RealScalar rcond_est = chollo.rcond(); // Verify that the estimated condition number is within a factor of 10 of the // truth. - VERIFY(rcond_est > rcond / 10 && rcond_est < rcond * 10); + VERIFY(rcond_est >= rcond / 10 && rcond_est <= rcond * 10); // test the upper mode LLT<SquareMatrixType,Upper> cholup(symmUp); @@ -113,12 +112,12 @@ template<typename MatrixType> void cholesky(const MatrixType& m) rcond = (RealScalar(1) / matrix_l1_norm<MatrixType, Upper>(symmUp)) / matrix_l1_norm<MatrixType, Upper>(symmUp_inverse); rcond_est = cholup.rcond(); - VERIFY(rcond_est > rcond / 10 && rcond_est < rcond * 10); + VERIFY(rcond_est >= rcond / 10 && rcond_est <= rcond * 10); MatrixType neg = -symmLo; chollo.compute(neg); - VERIFY(chollo.info()==NumericalIssue); + VERIFY(neg.size()==0 || chollo.info()==NumericalIssue); VERIFY_IS_APPROX(MatrixType(chollo.matrixL().transpose().conjugate()), MatrixType(chollo.matrixU())); VERIFY_IS_APPROX(MatrixType(chollo.matrixU().transpose().conjugate()), MatrixType(chollo.matrixL())); @@ -143,6 +142,9 @@ template<typename MatrixType> void cholesky(const MatrixType& m) // LDLT { + STATIC_CHECK(( internal::is_same<typename LDLT<MatrixType,Lower>::StorageIndex,int>::value )); + STATIC_CHECK(( internal::is_same<typename LDLT<MatrixType,Upper>::StorageIndex,int>::value )); + int sign = internal::random<int>()%2 ? 1 : -1; if(sign == -1) @@ -156,10 +158,9 @@ template<typename MatrixType> void cholesky(const MatrixType& m) LDLT<SquareMatrixType,Lower> ldltlo(symmLo); VERIFY(ldltlo.info()==Success); VERIFY_IS_APPROX(symm, ldltlo.reconstructedMatrix()); - vecX = ldltlo.solve(vecB); - VERIFY_IS_APPROX(symm * vecX, vecB); - matX = ldltlo.solve(matB); - VERIFY_IS_APPROX(symm * matX, matB); + + check_solverbase<VectorType, VectorType>(symm, ldltlo, rows, rows, 1); + check_solverbase<MatrixType, MatrixType>(symm, ldltlo, rows, cols, rows); const MatrixType symmLo_inverse = ldltlo.solve(MatrixType::Identity(rows,cols)); RealScalar rcond = (RealScalar(1) / matrix_l1_norm<MatrixType, Lower>(symmLo)) / @@ -167,7 +168,7 @@ template<typename MatrixType> void cholesky(const MatrixType& m) RealScalar rcond_est = ldltlo.rcond(); // Verify that the estimated condition number is within a factor of 10 of the // truth. - VERIFY(rcond_est > rcond / 10 && rcond_est < rcond * 10); + VERIFY(rcond_est >= rcond / 10 && rcond_est <= rcond * 10); LDLT<SquareMatrixType,Upper> ldltup(symmUp); @@ -184,7 +185,7 @@ template<typename MatrixType> void cholesky(const MatrixType& m) rcond = (RealScalar(1) / matrix_l1_norm<MatrixType, Upper>(symmUp)) / matrix_l1_norm<MatrixType, Upper>(symmUp_inverse); rcond_est = ldltup.rcond(); - VERIFY(rcond_est > rcond / 10 && rcond_est < rcond * 10); + 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())); @@ -289,8 +290,6 @@ template<typename MatrixType> void cholesky_cplx(const MatrixType& m) // test mixing real/scalar types - typedef typename MatrixType::Index Index; - Index rows = m.rows(); Index cols = m.cols(); @@ -315,10 +314,9 @@ template<typename MatrixType> void cholesky_cplx(const MatrixType& m) LLT<RealMatrixType,Lower> chollo(symmLo); VERIFY_IS_APPROX(symm, chollo.reconstructedMatrix()); - vecX = chollo.solve(vecB); - VERIFY_IS_APPROX(symm * vecX, vecB); -// matX = chollo.solve(matB); -// VERIFY_IS_APPROX(symm * matX, matB); + + check_solverbase<VectorType, VectorType>(symm, chollo, rows, rows, 1); + //check_solverbase<MatrixType, MatrixType>(symm, chollo, rows, cols, rows); } // LDLT @@ -335,10 +333,9 @@ template<typename MatrixType> void cholesky_cplx(const MatrixType& m) LDLT<RealMatrixType,Lower> ldltlo(symmLo); VERIFY(ldltlo.info()==Success); VERIFY_IS_APPROX(symm, ldltlo.reconstructedMatrix()); - vecX = ldltlo.solve(vecB); - VERIFY_IS_APPROX(symm * vecX, vecB); -// matX = ldltlo.solve(matB); -// VERIFY_IS_APPROX(symm * matX, matB); + + check_solverbase<VectorType, VectorType>(symm, ldltlo, rows, rows, 1); + //check_solverbase<MatrixType, MatrixType>(symm, ldltlo, rows, cols, rows); } } @@ -373,6 +370,7 @@ template<typename MatrixType> void cholesky_definiteness(const MatrixType& m) VERIFY(ldlt.info()==Success); VERIFY(!ldlt.isNegative()); VERIFY(!ldlt.isPositive()); + VERIFY_IS_APPROX(mat,ldlt.reconstructedMatrix()); } { mat << 1, 2, 2, 1; @@ -380,6 +378,7 @@ template<typename MatrixType> void cholesky_definiteness(const MatrixType& m) VERIFY(ldlt.info()==Success); VERIFY(!ldlt.isNegative()); VERIFY(!ldlt.isPositive()); + VERIFY_IS_APPROX(mat,ldlt.reconstructedMatrix()); } { mat << 0, 0, 0, 0; @@ -387,6 +386,7 @@ template<typename MatrixType> void cholesky_definiteness(const MatrixType& m) VERIFY(ldlt.info()==Success); VERIFY(ldlt.isNegative()); VERIFY(ldlt.isPositive()); + VERIFY_IS_APPROX(mat,ldlt.reconstructedMatrix()); } { mat << 0, 0, 0, 1; @@ -394,6 +394,7 @@ template<typename MatrixType> void cholesky_definiteness(const MatrixType& m) VERIFY(ldlt.info()==Success); VERIFY(!ldlt.isNegative()); VERIFY(ldlt.isPositive()); + VERIFY_IS_APPROX(mat,ldlt.reconstructedMatrix()); } { mat << -1, 0, 0, 0; @@ -401,6 +402,7 @@ template<typename MatrixType> void cholesky_definiteness(const MatrixType& m) VERIFY(ldlt.info()==Success); VERIFY(ldlt.isNegative()); VERIFY(!ldlt.isPositive()); + VERIFY_IS_APPROX(mat,ldlt.reconstructedMatrix()); } } @@ -452,6 +454,18 @@ void cholesky_faillure_cases() VERIFY(ldlt.info()==NumericalIssue); VERIFY_IS_NOT_APPROX(mat,ldlt.reconstructedMatrix()); } + + // bug 1479 + { + mat.resize(4,4); + mat << 1, 2, 0, 1, + 2, 4, 0, 2, + 0, 0, 0, 1, + 1, 2, 1, 1; + ldlt.compute(mat); + VERIFY(ldlt.info()==NumericalIssue); + VERIFY_IS_NOT_APPROX(mat,ldlt.reconstructedMatrix()); + } } template<typename MatrixType> void cholesky_verify_assert() @@ -462,19 +476,23 @@ template<typename MatrixType> void cholesky_verify_assert() VERIFY_RAISES_ASSERT(llt.matrixL()) VERIFY_RAISES_ASSERT(llt.matrixU()) VERIFY_RAISES_ASSERT(llt.solve(tmp)) - VERIFY_RAISES_ASSERT(llt.solveInPlace(&tmp)) + VERIFY_RAISES_ASSERT(llt.transpose().solve(tmp)) + VERIFY_RAISES_ASSERT(llt.adjoint().solve(tmp)) + VERIFY_RAISES_ASSERT(llt.solveInPlace(tmp)) LDLT<MatrixType> ldlt; VERIFY_RAISES_ASSERT(ldlt.matrixL()) - VERIFY_RAISES_ASSERT(ldlt.permutationP()) + VERIFY_RAISES_ASSERT(ldlt.transpositionsP()) VERIFY_RAISES_ASSERT(ldlt.vectorD()) VERIFY_RAISES_ASSERT(ldlt.isPositive()) VERIFY_RAISES_ASSERT(ldlt.isNegative()) VERIFY_RAISES_ASSERT(ldlt.solve(tmp)) - VERIFY_RAISES_ASSERT(ldlt.solveInPlace(&tmp)) + VERIFY_RAISES_ASSERT(ldlt.transpose().solve(tmp)) + VERIFY_RAISES_ASSERT(ldlt.adjoint().solve(tmp)) + VERIFY_RAISES_ASSERT(ldlt.solveInPlace(tmp)) } -void test_cholesky() +EIGEN_DECLARE_TEST(cholesky) { int s = 0; for(int i = 0; i < g_repeat; i++) { @@ -493,6 +511,11 @@ void test_cholesky() CALL_SUBTEST_6( cholesky_cplx(MatrixXcd(s,s)) ); TEST_SET_BUT_UNUSED_VARIABLE(s) } + // empty matrix, regression test for Bug 785: + CALL_SUBTEST_2( cholesky(MatrixXd(0,0)) ); + + // This does not work yet: + // CALL_SUBTEST_2( cholesky(Matrix<double,0,0>()) ); CALL_SUBTEST_4( cholesky_verify_assert<Matrix3f>() ); CALL_SUBTEST_7( cholesky_verify_assert<Matrix3d>() ); |