cholesky.cpp (18340B)
1 // This file is part of Eigen, a lightweight C++ template library 2 // for linear algebra. 3 // 4 // Copyright (C) 2008 Gael Guennebaud <gael.guennebaud@inria.fr> 5 // 6 // This Source Code Form is subject to the terms of the Mozilla 7 // Public License v. 2.0. If a copy of the MPL was not distributed 8 // with this file, You can obtain one at http://mozilla.org/MPL/2.0/. 9 10 #define TEST_ENABLE_TEMPORARY_TRACKING 11 12 #include "main.h" 13 #include <Eigen/Cholesky> 14 #include <Eigen/QR> 15 #include "solverbase.h" 16 17 template<typename MatrixType, int UpLo> 18 typename MatrixType::RealScalar matrix_l1_norm(const MatrixType& m) { 19 if(m.cols()==0) return typename MatrixType::RealScalar(0); 20 MatrixType symm = m.template selfadjointView<UpLo>(); 21 return symm.cwiseAbs().colwise().sum().maxCoeff(); 22 } 23 24 template<typename MatrixType,template <typename,int> class CholType> void test_chol_update(const MatrixType& symm) 25 { 26 typedef typename MatrixType::Scalar Scalar; 27 typedef typename MatrixType::RealScalar RealScalar; 28 typedef Matrix<Scalar, MatrixType::RowsAtCompileTime, 1> VectorType; 29 30 MatrixType symmLo = symm.template triangularView<Lower>(); 31 MatrixType symmUp = symm.template triangularView<Upper>(); 32 MatrixType symmCpy = symm; 33 34 CholType<MatrixType,Lower> chollo(symmLo); 35 CholType<MatrixType,Upper> cholup(symmUp); 36 37 for (int k=0; k<10; ++k) 38 { 39 VectorType vec = VectorType::Random(symm.rows()); 40 RealScalar sigma = internal::random<RealScalar>(); 41 symmCpy += sigma * vec * vec.adjoint(); 42 43 // we are doing some downdates, so it might be the case that the matrix is not SPD anymore 44 CholType<MatrixType,Lower> chol(symmCpy); 45 if(chol.info()!=Success) 46 break; 47 48 chollo.rankUpdate(vec, sigma); 49 VERIFY_IS_APPROX(symmCpy, chollo.reconstructedMatrix()); 50 51 cholup.rankUpdate(vec, sigma); 52 VERIFY_IS_APPROX(symmCpy, cholup.reconstructedMatrix()); 53 } 54 } 55 56 template<typename MatrixType> void cholesky(const MatrixType& m) 57 { 58 /* this test covers the following files: 59 LLT.h LDLT.h 60 */ 61 Index rows = m.rows(); 62 Index cols = m.cols(); 63 64 typedef typename MatrixType::Scalar Scalar; 65 typedef typename NumTraits<Scalar>::Real RealScalar; 66 typedef Matrix<Scalar, MatrixType::RowsAtCompileTime, MatrixType::RowsAtCompileTime> SquareMatrixType; 67 typedef Matrix<Scalar, MatrixType::RowsAtCompileTime, 1> VectorType; 68 69 MatrixType a0 = MatrixType::Random(rows,cols); 70 VectorType vecB = VectorType::Random(rows), vecX(rows); 71 MatrixType matB = MatrixType::Random(rows,cols), matX(rows,cols); 72 SquareMatrixType symm = a0 * a0.adjoint(); 73 // let's make sure the matrix is not singular or near singular 74 for (int k=0; k<3; ++k) 75 { 76 MatrixType a1 = MatrixType::Random(rows,cols); 77 symm += a1 * a1.adjoint(); 78 } 79 80 { 81 STATIC_CHECK(( internal::is_same<typename LLT<MatrixType,Lower>::StorageIndex,int>::value )); 82 STATIC_CHECK(( internal::is_same<typename LLT<MatrixType,Upper>::StorageIndex,int>::value )); 83 84 SquareMatrixType symmUp = symm.template triangularView<Upper>(); 85 SquareMatrixType symmLo = symm.template triangularView<Lower>(); 86 87 LLT<SquareMatrixType,Lower> chollo(symmLo); 88 VERIFY_IS_APPROX(symm, chollo.reconstructedMatrix()); 89 90 check_solverbase<VectorType, VectorType>(symm, chollo, rows, rows, 1); 91 check_solverbase<MatrixType, MatrixType>(symm, chollo, rows, cols, rows); 92 93 const MatrixType symmLo_inverse = chollo.solve(MatrixType::Identity(rows,cols)); 94 RealScalar rcond = (RealScalar(1) / matrix_l1_norm<MatrixType, Lower>(symmLo)) / 95 matrix_l1_norm<MatrixType, Lower>(symmLo_inverse); 96 RealScalar rcond_est = chollo.rcond(); 97 // Verify that the estimated condition number is within a factor of 10 of the 98 // truth. 99 VERIFY(rcond_est >= rcond / 10 && rcond_est <= rcond * 10); 100 101 // test the upper mode 102 LLT<SquareMatrixType,Upper> cholup(symmUp); 103 VERIFY_IS_APPROX(symm, cholup.reconstructedMatrix()); 104 vecX = cholup.solve(vecB); 105 VERIFY_IS_APPROX(symm * vecX, vecB); 106 matX = cholup.solve(matB); 107 VERIFY_IS_APPROX(symm * matX, matB); 108 109 // Verify that the estimated condition number is within a factor of 10 of the 110 // truth. 111 const MatrixType symmUp_inverse = cholup.solve(MatrixType::Identity(rows,cols)); 112 rcond = (RealScalar(1) / matrix_l1_norm<MatrixType, Upper>(symmUp)) / 113 matrix_l1_norm<MatrixType, Upper>(symmUp_inverse); 114 rcond_est = cholup.rcond(); 115 VERIFY(rcond_est >= rcond / 10 && rcond_est <= rcond * 10); 116 117 118 MatrixType neg = -symmLo; 119 chollo.compute(neg); 120 VERIFY(neg.size()==0 || chollo.info()==NumericalIssue); 121 122 VERIFY_IS_APPROX(MatrixType(chollo.matrixL().transpose().conjugate()), MatrixType(chollo.matrixU())); 123 VERIFY_IS_APPROX(MatrixType(chollo.matrixU().transpose().conjugate()), MatrixType(chollo.matrixL())); 124 VERIFY_IS_APPROX(MatrixType(cholup.matrixL().transpose().conjugate()), MatrixType(cholup.matrixU())); 125 VERIFY_IS_APPROX(MatrixType(cholup.matrixU().transpose().conjugate()), MatrixType(cholup.matrixL())); 126 127 // test some special use cases of SelfCwiseBinaryOp: 128 MatrixType m1 = MatrixType::Random(rows,cols), m2(rows,cols); 129 m2 = m1; 130 m2 += symmLo.template selfadjointView<Lower>().llt().solve(matB); 131 VERIFY_IS_APPROX(m2, m1 + symmLo.template selfadjointView<Lower>().llt().solve(matB)); 132 m2 = m1; 133 m2 -= symmLo.template selfadjointView<Lower>().llt().solve(matB); 134 VERIFY_IS_APPROX(m2, m1 - symmLo.template selfadjointView<Lower>().llt().solve(matB)); 135 m2 = m1; 136 m2.noalias() += symmLo.template selfadjointView<Lower>().llt().solve(matB); 137 VERIFY_IS_APPROX(m2, m1 + symmLo.template selfadjointView<Lower>().llt().solve(matB)); 138 m2 = m1; 139 m2.noalias() -= symmLo.template selfadjointView<Lower>().llt().solve(matB); 140 VERIFY_IS_APPROX(m2, m1 - symmLo.template selfadjointView<Lower>().llt().solve(matB)); 141 } 142 143 // LDLT 144 { 145 STATIC_CHECK(( internal::is_same<typename LDLT<MatrixType,Lower>::StorageIndex,int>::value )); 146 STATIC_CHECK(( internal::is_same<typename LDLT<MatrixType,Upper>::StorageIndex,int>::value )); 147 148 int sign = internal::random<int>()%2 ? 1 : -1; 149 150 if(sign == -1) 151 { 152 symm = -symm; // test a negative matrix 153 } 154 155 SquareMatrixType symmUp = symm.template triangularView<Upper>(); 156 SquareMatrixType symmLo = symm.template triangularView<Lower>(); 157 158 LDLT<SquareMatrixType,Lower> ldltlo(symmLo); 159 VERIFY(ldltlo.info()==Success); 160 VERIFY_IS_APPROX(symm, ldltlo.reconstructedMatrix()); 161 162 check_solverbase<VectorType, VectorType>(symm, ldltlo, rows, rows, 1); 163 check_solverbase<MatrixType, MatrixType>(symm, ldltlo, rows, cols, rows); 164 165 const MatrixType symmLo_inverse = ldltlo.solve(MatrixType::Identity(rows,cols)); 166 RealScalar rcond = (RealScalar(1) / matrix_l1_norm<MatrixType, Lower>(symmLo)) / 167 matrix_l1_norm<MatrixType, Lower>(symmLo_inverse); 168 RealScalar rcond_est = ldltlo.rcond(); 169 // Verify that the estimated condition number is within a factor of 10 of the 170 // truth. 171 VERIFY(rcond_est >= rcond / 10 && rcond_est <= rcond * 10); 172 173 174 LDLT<SquareMatrixType,Upper> ldltup(symmUp); 175 VERIFY(ldltup.info()==Success); 176 VERIFY_IS_APPROX(symm, ldltup.reconstructedMatrix()); 177 vecX = ldltup.solve(vecB); 178 VERIFY_IS_APPROX(symm * vecX, vecB); 179 matX = ldltup.solve(matB); 180 VERIFY_IS_APPROX(symm * matX, matB); 181 182 // Verify that the estimated condition number is within a factor of 10 of the 183 // truth. 184 const MatrixType symmUp_inverse = ldltup.solve(MatrixType::Identity(rows,cols)); 185 rcond = (RealScalar(1) / matrix_l1_norm<MatrixType, Upper>(symmUp)) / 186 matrix_l1_norm<MatrixType, Upper>(symmUp_inverse); 187 rcond_est = ldltup.rcond(); 188 VERIFY(rcond_est >= rcond / 10 && rcond_est <= rcond * 10); 189 190 VERIFY_IS_APPROX(MatrixType(ldltlo.matrixL().transpose().conjugate()), MatrixType(ldltlo.matrixU())); 191 VERIFY_IS_APPROX(MatrixType(ldltlo.matrixU().transpose().conjugate()), MatrixType(ldltlo.matrixL())); 192 VERIFY_IS_APPROX(MatrixType(ldltup.matrixL().transpose().conjugate()), MatrixType(ldltup.matrixU())); 193 VERIFY_IS_APPROX(MatrixType(ldltup.matrixU().transpose().conjugate()), MatrixType(ldltup.matrixL())); 194 195 if(MatrixType::RowsAtCompileTime==Dynamic) 196 { 197 // note : each inplace permutation requires a small temporary vector (mask) 198 199 // check inplace solve 200 matX = matB; 201 VERIFY_EVALUATION_COUNT(matX = ldltlo.solve(matX), 0); 202 VERIFY_IS_APPROX(matX, ldltlo.solve(matB).eval()); 203 204 205 matX = matB; 206 VERIFY_EVALUATION_COUNT(matX = ldltup.solve(matX), 0); 207 VERIFY_IS_APPROX(matX, ldltup.solve(matB).eval()); 208 } 209 210 // restore 211 if(sign == -1) 212 symm = -symm; 213 214 // check matrices coming from linear constraints with Lagrange multipliers 215 if(rows>=3) 216 { 217 SquareMatrixType A = symm; 218 Index c = internal::random<Index>(0,rows-2); 219 A.bottomRightCorner(c,c).setZero(); 220 // Make sure a solution exists: 221 vecX.setRandom(); 222 vecB = A * vecX; 223 vecX.setZero(); 224 ldltlo.compute(A); 225 VERIFY_IS_APPROX(A, ldltlo.reconstructedMatrix()); 226 vecX = ldltlo.solve(vecB); 227 VERIFY_IS_APPROX(A * vecX, vecB); 228 } 229 230 // check non-full rank matrices 231 if(rows>=3) 232 { 233 Index r = internal::random<Index>(1,rows-1); 234 Matrix<Scalar,Dynamic,Dynamic> a = Matrix<Scalar,Dynamic,Dynamic>::Random(rows,r); 235 SquareMatrixType A = a * a.adjoint(); 236 // Make sure a solution exists: 237 vecX.setRandom(); 238 vecB = A * vecX; 239 vecX.setZero(); 240 ldltlo.compute(A); 241 VERIFY_IS_APPROX(A, ldltlo.reconstructedMatrix()); 242 vecX = ldltlo.solve(vecB); 243 VERIFY_IS_APPROX(A * vecX, vecB); 244 } 245 246 // check matrices with a wide spectrum 247 if(rows>=3) 248 { 249 using std::pow; 250 using std::sqrt; 251 RealScalar s = (std::min)(16,std::numeric_limits<RealScalar>::max_exponent10/8); 252 Matrix<Scalar,Dynamic,Dynamic> a = Matrix<Scalar,Dynamic,Dynamic>::Random(rows,rows); 253 Matrix<RealScalar,Dynamic,1> d = Matrix<RealScalar,Dynamic,1>::Random(rows); 254 for(Index k=0; k<rows; ++k) 255 d(k) = d(k)*pow(RealScalar(10),internal::random<RealScalar>(-s,s)); 256 SquareMatrixType A = a * d.asDiagonal() * a.adjoint(); 257 // Make sure a solution exists: 258 vecX.setRandom(); 259 vecB = A * vecX; 260 vecX.setZero(); 261 ldltlo.compute(A); 262 VERIFY_IS_APPROX(A, ldltlo.reconstructedMatrix()); 263 vecX = ldltlo.solve(vecB); 264 265 if(ldltlo.vectorD().real().cwiseAbs().minCoeff()>RealScalar(0)) 266 { 267 VERIFY_IS_APPROX(A * vecX,vecB); 268 } 269 else 270 { 271 RealScalar large_tol = sqrt(test_precision<RealScalar>()); 272 VERIFY((A * vecX).isApprox(vecB, large_tol)); 273 274 ++g_test_level; 275 VERIFY_IS_APPROX(A * vecX,vecB); 276 --g_test_level; 277 } 278 } 279 } 280 281 // update/downdate 282 CALL_SUBTEST(( test_chol_update<SquareMatrixType,LLT>(symm) )); 283 CALL_SUBTEST(( test_chol_update<SquareMatrixType,LDLT>(symm) )); 284 } 285 286 template<typename MatrixType> void cholesky_cplx(const MatrixType& m) 287 { 288 // classic test 289 cholesky(m); 290 291 // test mixing real/scalar types 292 293 Index rows = m.rows(); 294 Index cols = m.cols(); 295 296 typedef typename MatrixType::Scalar Scalar; 297 typedef typename NumTraits<Scalar>::Real RealScalar; 298 typedef Matrix<RealScalar, MatrixType::RowsAtCompileTime, MatrixType::RowsAtCompileTime> RealMatrixType; 299 typedef Matrix<Scalar, MatrixType::RowsAtCompileTime, 1> VectorType; 300 301 RealMatrixType a0 = RealMatrixType::Random(rows,cols); 302 VectorType vecB = VectorType::Random(rows), vecX(rows); 303 MatrixType matB = MatrixType::Random(rows,cols), matX(rows,cols); 304 RealMatrixType symm = a0 * a0.adjoint(); 305 // let's make sure the matrix is not singular or near singular 306 for (int k=0; k<3; ++k) 307 { 308 RealMatrixType a1 = RealMatrixType::Random(rows,cols); 309 symm += a1 * a1.adjoint(); 310 } 311 312 { 313 RealMatrixType symmLo = symm.template triangularView<Lower>(); 314 315 LLT<RealMatrixType,Lower> chollo(symmLo); 316 VERIFY_IS_APPROX(symm, chollo.reconstructedMatrix()); 317 318 check_solverbase<VectorType, VectorType>(symm, chollo, rows, rows, 1); 319 //check_solverbase<MatrixType, MatrixType>(symm, chollo, rows, cols, rows); 320 } 321 322 // LDLT 323 { 324 int sign = internal::random<int>()%2 ? 1 : -1; 325 326 if(sign == -1) 327 { 328 symm = -symm; // test a negative matrix 329 } 330 331 RealMatrixType symmLo = symm.template triangularView<Lower>(); 332 333 LDLT<RealMatrixType,Lower> ldltlo(symmLo); 334 VERIFY(ldltlo.info()==Success); 335 VERIFY_IS_APPROX(symm, ldltlo.reconstructedMatrix()); 336 337 check_solverbase<VectorType, VectorType>(symm, ldltlo, rows, rows, 1); 338 //check_solverbase<MatrixType, MatrixType>(symm, ldltlo, rows, cols, rows); 339 } 340 } 341 342 // regression test for bug 241 343 template<typename MatrixType> void cholesky_bug241(const MatrixType& m) 344 { 345 eigen_assert(m.rows() == 2 && m.cols() == 2); 346 347 typedef typename MatrixType::Scalar Scalar; 348 typedef Matrix<Scalar, MatrixType::RowsAtCompileTime, 1> VectorType; 349 350 MatrixType matA; 351 matA << 1, 1, 1, 1; 352 VectorType vecB; 353 vecB << 1, 1; 354 VectorType vecX = matA.ldlt().solve(vecB); 355 VERIFY_IS_APPROX(matA * vecX, vecB); 356 } 357 358 // LDLT is not guaranteed to work for indefinite matrices, but happens to work fine if matrix is diagonal. 359 // This test checks that LDLT reports correctly that matrix is indefinite. 360 // See http://forum.kde.org/viewtopic.php?f=74&t=106942 and bug 736 361 template<typename MatrixType> void cholesky_definiteness(const MatrixType& m) 362 { 363 eigen_assert(m.rows() == 2 && m.cols() == 2); 364 MatrixType mat; 365 LDLT<MatrixType> ldlt(2); 366 367 { 368 mat << 1, 0, 0, -1; 369 ldlt.compute(mat); 370 VERIFY(ldlt.info()==Success); 371 VERIFY(!ldlt.isNegative()); 372 VERIFY(!ldlt.isPositive()); 373 VERIFY_IS_APPROX(mat,ldlt.reconstructedMatrix()); 374 } 375 { 376 mat << 1, 2, 2, 1; 377 ldlt.compute(mat); 378 VERIFY(ldlt.info()==Success); 379 VERIFY(!ldlt.isNegative()); 380 VERIFY(!ldlt.isPositive()); 381 VERIFY_IS_APPROX(mat,ldlt.reconstructedMatrix()); 382 } 383 { 384 mat << 0, 0, 0, 0; 385 ldlt.compute(mat); 386 VERIFY(ldlt.info()==Success); 387 VERIFY(ldlt.isNegative()); 388 VERIFY(ldlt.isPositive()); 389 VERIFY_IS_APPROX(mat,ldlt.reconstructedMatrix()); 390 } 391 { 392 mat << 0, 0, 0, 1; 393 ldlt.compute(mat); 394 VERIFY(ldlt.info()==Success); 395 VERIFY(!ldlt.isNegative()); 396 VERIFY(ldlt.isPositive()); 397 VERIFY_IS_APPROX(mat,ldlt.reconstructedMatrix()); 398 } 399 { 400 mat << -1, 0, 0, 0; 401 ldlt.compute(mat); 402 VERIFY(ldlt.info()==Success); 403 VERIFY(ldlt.isNegative()); 404 VERIFY(!ldlt.isPositive()); 405 VERIFY_IS_APPROX(mat,ldlt.reconstructedMatrix()); 406 } 407 } 408 409 template<typename> 410 void cholesky_faillure_cases() 411 { 412 MatrixXd mat; 413 LDLT<MatrixXd> ldlt; 414 415 { 416 mat.resize(2,2); 417 mat << 0, 1, 1, 0; 418 ldlt.compute(mat); 419 VERIFY_IS_NOT_APPROX(mat,ldlt.reconstructedMatrix()); 420 VERIFY(ldlt.info()==NumericalIssue); 421 } 422 #if (!EIGEN_ARCH_i386) || defined(EIGEN_VECTORIZE_SSE2) 423 { 424 mat.resize(3,3); 425 mat << -1, -3, 3, 426 -3, -8.9999999999999999999, 1, 427 3, 1, 0; 428 ldlt.compute(mat); 429 VERIFY(ldlt.info()==NumericalIssue); 430 VERIFY_IS_NOT_APPROX(mat,ldlt.reconstructedMatrix()); 431 } 432 #endif 433 { 434 mat.resize(3,3); 435 mat << 1, 2, 3, 436 2, 4, 1, 437 3, 1, 0; 438 ldlt.compute(mat); 439 VERIFY(ldlt.info()==NumericalIssue); 440 VERIFY_IS_NOT_APPROX(mat,ldlt.reconstructedMatrix()); 441 } 442 443 { 444 mat.resize(8,8); 445 mat << 0.1, 0, -0.1, 0, 0, 0, 1, 0, 446 0, 4.24667, 0, 2.00333, 0, 0, 0, 0, 447 -0.1, 0, 0.2, 0, -0.1, 0, 0, 0, 448 0, 2.00333, 0, 8.49333, 0, 2.00333, 0, 0, 449 0, 0, -0.1, 0, 0.1, 0, 0, 1, 450 0, 0, 0, 2.00333, 0, 4.24667, 0, 0, 451 1, 0, 0, 0, 0, 0, 0, 0, 452 0, 0, 0, 0, 1, 0, 0, 0; 453 ldlt.compute(mat); 454 VERIFY(ldlt.info()==NumericalIssue); 455 VERIFY_IS_NOT_APPROX(mat,ldlt.reconstructedMatrix()); 456 } 457 458 // bug 1479 459 { 460 mat.resize(4,4); 461 mat << 1, 2, 0, 1, 462 2, 4, 0, 2, 463 0, 0, 0, 1, 464 1, 2, 1, 1; 465 ldlt.compute(mat); 466 VERIFY(ldlt.info()==NumericalIssue); 467 VERIFY_IS_NOT_APPROX(mat,ldlt.reconstructedMatrix()); 468 } 469 } 470 471 template<typename MatrixType> void cholesky_verify_assert() 472 { 473 MatrixType tmp; 474 475 LLT<MatrixType> llt; 476 VERIFY_RAISES_ASSERT(llt.matrixL()) 477 VERIFY_RAISES_ASSERT(llt.matrixU()) 478 VERIFY_RAISES_ASSERT(llt.solve(tmp)) 479 VERIFY_RAISES_ASSERT(llt.transpose().solve(tmp)) 480 VERIFY_RAISES_ASSERT(llt.adjoint().solve(tmp)) 481 VERIFY_RAISES_ASSERT(llt.solveInPlace(tmp)) 482 483 LDLT<MatrixType> ldlt; 484 VERIFY_RAISES_ASSERT(ldlt.matrixL()) 485 VERIFY_RAISES_ASSERT(ldlt.transpositionsP()) 486 VERIFY_RAISES_ASSERT(ldlt.vectorD()) 487 VERIFY_RAISES_ASSERT(ldlt.isPositive()) 488 VERIFY_RAISES_ASSERT(ldlt.isNegative()) 489 VERIFY_RAISES_ASSERT(ldlt.solve(tmp)) 490 VERIFY_RAISES_ASSERT(ldlt.transpose().solve(tmp)) 491 VERIFY_RAISES_ASSERT(ldlt.adjoint().solve(tmp)) 492 VERIFY_RAISES_ASSERT(ldlt.solveInPlace(tmp)) 493 } 494 495 EIGEN_DECLARE_TEST(cholesky) 496 { 497 int s = 0; 498 for(int i = 0; i < g_repeat; i++) { 499 CALL_SUBTEST_1( cholesky(Matrix<double,1,1>()) ); 500 CALL_SUBTEST_3( cholesky(Matrix2d()) ); 501 CALL_SUBTEST_3( cholesky_bug241(Matrix2d()) ); 502 CALL_SUBTEST_3( cholesky_definiteness(Matrix2d()) ); 503 CALL_SUBTEST_4( cholesky(Matrix3f()) ); 504 CALL_SUBTEST_5( cholesky(Matrix4d()) ); 505 506 s = internal::random<int>(1,EIGEN_TEST_MAX_SIZE); 507 CALL_SUBTEST_2( cholesky(MatrixXd(s,s)) ); 508 TEST_SET_BUT_UNUSED_VARIABLE(s) 509 510 s = internal::random<int>(1,EIGEN_TEST_MAX_SIZE/2); 511 CALL_SUBTEST_6( cholesky_cplx(MatrixXcd(s,s)) ); 512 TEST_SET_BUT_UNUSED_VARIABLE(s) 513 } 514 // empty matrix, regression test for Bug 785: 515 CALL_SUBTEST_2( cholesky(MatrixXd(0,0)) ); 516 517 // This does not work yet: 518 // CALL_SUBTEST_2( cholesky(Matrix<double,0,0>()) ); 519 520 CALL_SUBTEST_4( cholesky_verify_assert<Matrix3f>() ); 521 CALL_SUBTEST_7( cholesky_verify_assert<Matrix3d>() ); 522 CALL_SUBTEST_8( cholesky_verify_assert<MatrixXf>() ); 523 CALL_SUBTEST_2( cholesky_verify_assert<MatrixXd>() ); 524 525 // Test problem size constructors 526 CALL_SUBTEST_9( LLT<MatrixXf>(10) ); 527 CALL_SUBTEST_9( LDLT<MatrixXf>(10) ); 528 529 CALL_SUBTEST_2( cholesky_faillure_cases<void>() ); 530 531 TEST_SET_BUT_UNUSED_VARIABLE(nb_temporaries) 532 }