Standard Cholesky decomposition (LL^T) of a matrix and associated features. More...
#include <LLT.h>
Public Types | |
enum | { RowsAtCompileTime, ColsAtCompileTime, Options, MaxColsAtCompileTime } |
enum | { PacketSize, AlignmentMask, UpLo } |
typedef MatrixType::Index | Index |
typedef _MatrixType | MatrixType |
typedef NumTraits< typename MatrixType::Scalar >::Real | RealScalar |
typedef MatrixType::Scalar | Scalar |
typedef internal::LLT_Traits < MatrixType, UpLo > | Traits |
Public Member Functions | |
Index | cols () const |
LLT & | compute (const MatrixType &matrix) |
ComputationInfo | info () const |
Reports whether previous computation was successful. | |
LLT () | |
Default Constructor. | |
LLT (Index size) | |
Default Constructor with memory preallocation. | |
LLT (const MatrixType &matrix) | |
Traits::MatrixL | matrixL () const |
const MatrixType & | matrixLLT () const |
Traits::MatrixU | matrixU () const |
template<typename VectorType > | |
LLT | rankUpdate (const VectorType &vec, const RealScalar &sigma=1) |
MatrixType | reconstructedMatrix () const |
Index | rows () const |
template<typename Rhs > | |
const internal::solve_retval < LLT, Rhs > | solve (const MatrixBase< Rhs > &b) const |
template<typename Derived > | |
void | solveInPlace (MatrixBase< Derived > &bAndX) const |
Protected Attributes | |
ComputationInfo | m_info |
bool | m_isInitialized |
MatrixType | m_matrix |
Standard Cholesky decomposition (LL^T) of a matrix and associated features.
MatrixType | the type of the matrix of which we are computing the LL^T Cholesky decomposition |
UpLo | the triangular part that will be used for the decompositon: Lower (default) or Upper. The other triangular part won't be read. |
This class performs a LL^T Cholesky decomposition of a symmetric, positive definite matrix A such that A = LL^* = U^*U, where L is lower triangular.
While the Cholesky decomposition is particularly useful to solve selfadjoint problems like D^*D x = b, for that purpose, we recommend the Cholesky decomposition without square root which is more stable and even faster. Nevertheless, this standard Cholesky decomposition remains useful in many other situations like generalised eigen problems with hermitian matrices.
Remember that Cholesky decompositions are not rank-revealing. This LLT decomposition is only stable on positive definite matrices, use LDLT instead for the semidefinite case. Also, do not use a Cholesky decomposition to determine whether a system of equations has a solution.
Example:
Output:
The matrix A is 4 -1 2 -1 6 0 2 0 5 The Cholesky factor L is 2 0 0 -0.5 2.4 0 1 0.209 1.99 To check this, let us compute L * L.transpose() 4 -1 2 -1 6 0 2 0 5 This should equal the matrix A
typedef MatrixType::Index Index |
typedef _MatrixType MatrixType |
typedef NumTraits<typename MatrixType::Scalar>::Real RealScalar |
typedef MatrixType::Scalar Scalar |
typedef internal::LLT_Traits<MatrixType,UpLo> Traits |
|
inline |
Default Constructor.
The default constructor is useful in cases in which the user intends to perform decompositions via LLT::compute(const MatrixType&).
Default Constructor with memory preallocation.
Like the default constructor but with preallocation of the internal data according to the specified problem size.
|
inline |
References LLT< _MatrixType, _UpLo >::compute().
|
inline |
References LLT< _MatrixType, _UpLo >::m_matrix.
LLT< MatrixType, _UpLo > & compute | ( | const MatrixType & | a | ) |
Computes / recomputes the Cholesky decomposition A = LL^* = U^*U of matrix
Example:
Output:
Here is the matrix A: 2 -1 -1 3 Here is the right hand side b: 1 2 3 1 Computing LLT decomposition... The solution is: 1.2 1.4 1.4 0.8 The matrix A is now: 2 -1 -1 4 Computing LLT decomposition... The solution is now: 1 1.29 1 0.571
References eigen_assert, Eigen::NumericalIssue, and Eigen::Success.
Referenced by LLT< _MatrixType, _UpLo >::LLT().
|
inline |
Reports whether previous computation was successful.
Success
if computation was succesful, NumericalIssue
if the matrix.appears to be negative. References eigen_assert, LLT< _MatrixType, _UpLo >::m_info, and LLT< _MatrixType, _UpLo >::m_isInitialized.
|
inline |
References eigen_assert, LLT< _MatrixType, _UpLo >::m_isInitialized, and LLT< _MatrixType, _UpLo >::m_matrix.
Referenced by GeneralizedSelfAdjointEigenSolver< _MatrixType >::compute().
|
inline |
TODO: document the storage layout
References eigen_assert, LLT< _MatrixType, _UpLo >::m_isInitialized, and LLT< _MatrixType, _UpLo >::m_matrix.
|
inline |
References eigen_assert, LLT< _MatrixType, _UpLo >::m_isInitialized, and LLT< _MatrixType, _UpLo >::m_matrix.
Referenced by GeneralizedSelfAdjointEigenSolver< _MatrixType >::compute().
LLT< _MatrixType, _UpLo > rankUpdate | ( | const VectorType & | v, |
const RealScalar & | sigma = 1 |
||
) |
Performs a rank one update (or dowdate) of the current decomposition. If A = LL^* before the rank one update, then after it we have LL^* = A + sigma * v v^* where v must be a vector of same dimension.
References eigen_assert, EIGEN_STATIC_ASSERT_VECTOR_ONLY, Eigen::NumericalIssue, and Eigen::Success.
MatrixType reconstructedMatrix | ( | ) | const |
References eigen_assert.
|
inline |
References LLT< _MatrixType, _UpLo >::m_matrix.
|
inline |
Since this LLT class assumes anyway that the matrix A is invertible, the solution theoretically exists and is unique regardless of b.
Example:
Output:
2.02 2.97
References eigen_assert, LLT< _MatrixType, _UpLo >::m_isInitialized, and LLT< _MatrixType, _UpLo >::m_matrix.
void solveInPlace | ( | MatrixBase< Derived > & | bAndX | ) | const |
References eigen_assert.
|
protected |
Referenced by LLT< _MatrixType, _UpLo >::info().
|
protected |
|
protected |