![]() |
|
Standard Cholesky decomposition (LL^T) of a matrix and associated features. More...
#include <LLT.h>
Public Types | |
enum | { MaxColsAtCompileTime = MatrixType::MaxColsAtCompileTime } |
enum | { PacketSize = internal::packet_traits<Scalar>::size , AlignmentMask = int(PacketSize) - 1 , UpLo = UpLo_ } |
typedef MatrixType_ | MatrixType |
typedef SolverBase< LLT > | Base |
typedef internal::LLT_Traits< MatrixType, UpLo > | Traits |
![]() | |
enum | |
typedef EigenBase< LLT< MatrixType_, UpLo_ > > | Base |
typedef internal::traits< LLT< MatrixType_, UpLo_ > >::Scalar | Scalar |
typedef Scalar | CoeffReturnType |
typedef Transpose< const LLT< MatrixType_, UpLo_ > > | ConstTransposeReturnType |
typedef std::conditional_t< NumTraits< Scalar >::IsComplex, CwiseUnaryOp< internal::scalar_conjugate_op< Scalar >, const ConstTransposeReturnType >, const ConstTransposeReturnType > | AdjointReturnType |
![]() | |
typedef Eigen::Index | Index |
The interface type of indices. More... | |
typedef internal::traits< Derived >::StorageKind | StorageKind |
Public Member Functions | |
LLT () | |
Default Constructor. More... | |
LLT (Index size) | |
Default Constructor with memory preallocation. More... | |
template<typename InputType > | |
LLT (const EigenBase< InputType > &matrix) | |
template<typename InputType > | |
LLT (EigenBase< InputType > &matrix) | |
Constructs a LLT factorization from a given matrix. More... | |
Traits::MatrixU | matrixU () const |
Traits::MatrixL | matrixL () const |
template<typename Derived > | |
void | solveInPlace (const MatrixBase< Derived > &bAndX) const |
template<typename InputType > | |
LLT & | compute (const EigenBase< InputType > &matrix) |
RealScalar | rcond () const |
const MatrixType & | matrixLLT () const |
MatrixType | reconstructedMatrix () const |
ComputationInfo | info () const |
Reports whether previous computation was successful. More... | |
const LLT & | adjoint () const EIGEN_NOEXCEPT |
EIGEN_CONSTEXPR Index | rows () const EIGEN_NOEXCEPT |
EIGEN_CONSTEXPR Index | cols () const EIGEN_NOEXCEPT |
template<typename VectorType > | |
LLT & | rankUpdate (const VectorType &vec, const RealScalar &sigma=1) |
template<typename RhsType , typename DstType > | |
void | _solve_impl (const RhsType &rhs, DstType &dst) const |
template<bool Conjugate, typename RhsType , typename DstType > | |
void | _solve_impl_transposed (const RhsType &rhs, DstType &dst) const |
template<typename InputType > | |
LLT< MatrixType, UpLo_ > & | compute (const EigenBase< InputType > &a) |
template<typename VectorType > | |
LLT< MatrixType_, UpLo_ > & | rankUpdate (const VectorType &v, const RealScalar &sigma) |
![]() | |
SolverBase () | |
~SolverBase () | |
const Solve< LLT< MatrixType_, UpLo_ >, Rhs > | solve (const MatrixBase< Rhs > &b) const |
const ConstTransposeReturnType | transpose () const |
const AdjointReturnType | adjoint () const |
constexpr EIGEN_DEVICE_FUNC LLT< MatrixType_, UpLo_ > & | derived () |
constexpr EIGEN_DEVICE_FUNC const LLT< MatrixType_, UpLo_ > & | derived () const |
![]() | |
constexpr EIGEN_DEVICE_FUNC Derived & | derived () |
constexpr EIGEN_DEVICE_FUNC const Derived & | derived () const |
EIGEN_DEVICE_FUNC Derived & | const_cast_derived () const |
EIGEN_DEVICE_FUNC const Derived & | const_derived () const |
EIGEN_DEVICE_FUNC EIGEN_CONSTEXPR Index | rows () const EIGEN_NOEXCEPT |
EIGEN_DEVICE_FUNC EIGEN_CONSTEXPR Index | cols () const EIGEN_NOEXCEPT |
EIGEN_DEVICE_FUNC EIGEN_CONSTEXPR Index | size () const EIGEN_NOEXCEPT |
template<typename Dest > | |
EIGEN_DEVICE_FUNC void | evalTo (Dest &dst) const |
template<typename Dest > | |
EIGEN_DEVICE_FUNC void | addTo (Dest &dst) const |
template<typename Dest > | |
EIGEN_DEVICE_FUNC void | subTo (Dest &dst) const |
template<typename Dest > | |
EIGEN_DEVICE_FUNC void | applyThisOnTheRight (Dest &dst) const |
template<typename Dest > | |
EIGEN_DEVICE_FUNC void | applyThisOnTheLeft (Dest &dst) const |
template<typename Device > | |
EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE DeviceWrapper< Derived, Device > | device (Device &device) |
template<typename Device > | |
EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE DeviceWrapper< const Derived, Device > | device (Device &device) const |
Protected Attributes | |
MatrixType | m_matrix |
RealScalar | m_l1_norm |
bool | m_isInitialized |
ComputationInfo | m_info |
Friends | |
class | SolverBase< LLT > |
Additional Inherited Members | |
![]() | |
void | _check_solve_assertion (const Rhs &b) const |
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 decomposition: 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:
Performance: for best performance, it is recommended to use a column-major storage format with the Lower triangular part (the default), or, equivalently, a row-major storage format with the Upper triangular part. Otherwise, you might get a 20% slowdown for the full factorization step, and rank-updates can be up to 3 times slower.
This class supports the inplace decomposition mechanism.
Note that during the decomposition, only the lower (or upper, as defined by UpLo_) triangular part of A is considered. Therefore, the strict lower part does not have to store correct values.
typedef SolverBase<LLT> Eigen::LLT< MatrixType_, UpLo_ >::Base |
typedef MatrixType_ Eigen::LLT< MatrixType_, UpLo_ >::MatrixType |
typedef internal::LLT_Traits<MatrixType, UpLo> Eigen::LLT< MatrixType_, UpLo_ >::Traits |
anonymous enum |
Enumerator | |
---|---|
PacketSize | |
AlignmentMask | |
UpLo |
|
inline |
Default Constructor.
The default constructor is useful in cases in which the user intends to perform decompositions via LLT::compute(const MatrixType&).
|
inlineexplicit |
Default Constructor with memory preallocation.
Like the default constructor but with preallocation of the internal data according to the specified problem size.
|
inlineexplicit |
References Eigen::LLT< MatrixType_, UpLo_ >::compute(), and matrix().
|
inlineexplicit |
Constructs a LLT factorization from a given matrix.
This overloaded constructor is provided for inplace decomposition when MatrixType
is a Eigen::Ref.
References Eigen::LLT< MatrixType_, UpLo_ >::compute(), and matrix().
void Eigen::LLT< MatrixType_, UpLo_ >::_solve_impl | ( | const RhsType & | rhs, |
DstType & | dst | ||
) | const |
void Eigen::LLT< MatrixType_, UpLo_ >::_solve_impl_transposed | ( | const RhsType & | rhs, |
DstType & | dst | ||
) | const |
|
inline |
*this
, that is, a const reference to the decomposition itself as the underlying matrix is self-adjoint.This method is provided for compatibility with other matrix decompositions, thus enabling generic code such as:
|
inline |
References Eigen::LLT< MatrixType_, UpLo_ >::m_matrix.
Referenced by gdb.printers._MatrixEntryIterator::__next__(), gdb.printers.EigenMatrixPrinter::children(), gdb.printers.EigenSparseMatrixPrinter::children(), gdb.printers.EigenMatrixPrinter::to_string(), and gdb.printers.EigenSparseMatrixPrinter::to_string().
LLT<MatrixType, UpLo_>& Eigen::LLT< MatrixType_, UpLo_ >::compute | ( | const EigenBase< InputType > & | a | ) |
Computes / recomputes the Cholesky decomposition A = LL^* = U^*U of matrix
Example:
Output:
References a, col(), eigen_assert, Eigen::internal::is_same_dense(), Eigen::Lower, Eigen::LLT< MatrixType_, UpLo_ >::m_info, Eigen::LLT< MatrixType_, UpLo_ >::m_isInitialized, Eigen::LLT< MatrixType_, UpLo_ >::m_l1_norm, Eigen::LLT< MatrixType_, UpLo_ >::m_matrix, Eigen::NumericalIssue, Eigen::EigenBase< Derived >::size(), and Eigen::Success.
LLT& Eigen::LLT< MatrixType_, UpLo_ >::compute | ( | const EigenBase< InputType > & | matrix | ) |
Referenced by cholesky(), ctms_decompositions(), and Eigen::LLT< MatrixType_, UpLo_ >::LLT().
|
inline |
Reports whether previous computation was successful.
Success
if computation was successful, NumericalIssue
if the matrix.appears not to be positive definite. References eigen_assert, Eigen::LLT< MatrixType_, UpLo_ >::m_info, and Eigen::LLT< MatrixType_, UpLo_ >::m_isInitialized.
Referenced by cholesky().
|
inline |
References eigen_assert, Eigen::LLT< MatrixType_, UpLo_ >::m_isInitialized, and Eigen::LLT< MatrixType_, UpLo_ >::m_matrix.
Referenced by __attribute__(), Eigen::LLT< MatrixType_, UpLo_ >::_solve_impl_transposed(), cholesky(), Eigen::GeneralizedSelfAdjointEigenSolver< MatrixType_ >::compute(), main(), Eigen::LLT< MatrixType_, UpLo_ >::reconstructedMatrix(), and Eigen::LLT< MatrixType_, UpLo_ >::solveInPlace().
|
inline |
TODO: document the storage layout
References eigen_assert, Eigen::LLT< MatrixType_, UpLo_ >::m_isInitialized, and Eigen::LLT< MatrixType_, UpLo_ >::m_matrix.
|
inline |
References eigen_assert, Eigen::LLT< MatrixType_, UpLo_ >::m_isInitialized, and Eigen::LLT< MatrixType_, UpLo_ >::m_matrix.
Referenced by Eigen::LLT< MatrixType_, UpLo_ >::_solve_impl_transposed(), cholesky(), Eigen::GeneralizedSelfAdjointEigenSolver< MatrixType_ >::compute(), and Eigen::LLT< MatrixType_, UpLo_ >::solveInPlace().
LLT<MatrixType_, UpLo_>& Eigen::LLT< MatrixType_, UpLo_ >::rankUpdate | ( | const VectorType & | v, |
const RealScalar & | sigma | ||
) |
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::LLT< MatrixType_, UpLo_ >::m_info, Eigen::LLT< MatrixType_, UpLo_ >::m_isInitialized, Eigen::LLT< MatrixType_, UpLo_ >::m_matrix, Eigen::NumericalIssue, calibrate::sigma, Eigen::Success, and v.
LLT& Eigen::LLT< MatrixType_, UpLo_ >::rankUpdate | ( | const VectorType & | vec, |
const RealScalar & | sigma = 1 |
||
) |
|
inline |
*this
is the Cholesky decomposition. References eigen_assert, Eigen::LLT< MatrixType_, UpLo_ >::m_info, Eigen::LLT< MatrixType_, UpLo_ >::m_isInitialized, Eigen::LLT< MatrixType_, UpLo_ >::m_l1_norm, Eigen::internal::rcond_estimate_helper(), and Eigen::Success.
Referenced by cholesky().
MatrixType Eigen::LLT< MatrixType, UpLo_ >::reconstructedMatrix |
References eigen_assert, Eigen::LLT< MatrixType_, UpLo_ >::m_isInitialized, and Eigen::LLT< MatrixType_, UpLo_ >::matrixL().
Referenced by cholesky(), and cholesky_cplx().
|
inline |
References Eigen::LLT< MatrixType_, UpLo_ >::m_matrix.
Referenced by gdb.printers._MatrixEntryIterator::__next__(), gdb.printers.EigenMatrixPrinter::children(), gdb.printers.EigenSparseMatrixPrinter::children(), gdb.printers.EigenMatrixPrinter::to_string(), and gdb.printers.EigenSparseMatrixPrinter::to_string().
void Eigen::LLT< MatrixType, UpLo_ >::solveInPlace | ( | const MatrixBase< Derived > & | bAndX | ) | const |
use x = llt_object.solve(x);
This is the in-place version of solve().
bAndX | represents both the right-hand side matrix b and result x. |
This version avoids a copy when the right hand side matrix b is not needed anymore.
References eigen_assert, Eigen::LLT< MatrixType_, UpLo_ >::m_isInitialized, Eigen::LLT< MatrixType_, UpLo_ >::m_matrix, Eigen::LLT< MatrixType_, UpLo_ >::matrixL(), and Eigen::LLT< MatrixType_, UpLo_ >::matrixU().
|
friend |
|
protected |
|
protected |
Referenced by Eigen::LLT< MatrixType_, UpLo_ >::compute(), Eigen::LLT< MatrixType_, UpLo_ >::info(), Eigen::LLT< MatrixType_, UpLo_ >::matrixL(), Eigen::LLT< MatrixType_, UpLo_ >::matrixLLT(), Eigen::LLT< MatrixType_, UpLo_ >::matrixU(), Eigen::LLT< MatrixType_, UpLo_ >::rankUpdate(), Eigen::LLT< MatrixType_, UpLo_ >::rcond(), Eigen::LLT< MatrixType_, UpLo_ >::reconstructedMatrix(), and Eigen::LLT< MatrixType_, UpLo_ >::solveInPlace().
|
protected |
Referenced by Eigen::LLT< MatrixType_, UpLo_ >::compute(), and Eigen::LLT< MatrixType_, UpLo_ >::rcond().
|
protected |
Used to compute and store L The strict upper part is not used and even not initialized.
Referenced by Eigen::LLT< MatrixType_, UpLo_ >::cols(), Eigen::LLT< MatrixType_, UpLo_ >::compute(), Eigen::SelfAdjointView< MatrixType_, UpLo >::llt(), Eigen::LLT< MatrixType_, UpLo_ >::matrixL(), Eigen::LLT< MatrixType_, UpLo_ >::matrixLLT(), Eigen::LLT< MatrixType_, UpLo_ >::matrixU(), Eigen::LLT< MatrixType_, UpLo_ >::rankUpdate(), Eigen::LLT< MatrixType_, UpLo_ >::rows(), and Eigen::LLT< MatrixType_, UpLo_ >::solveInPlace().