Eigen::LLT< MatrixType_, UpLo_ > Class Template Reference

Standard Cholesky decomposition (LL^T) of a matrix and associated features. More...

#include <LLT.h>

+ Inheritance diagram for Eigen::LLT< MatrixType_, UpLo_ >:

Public Types

enum  { MaxColsAtCompileTime = MatrixType::MaxColsAtCompileTime }
 
enum  { PacketSize = internal::packet_traits<Scalar>::size , AlignmentMask = int(PacketSize) - 1 , UpLo = UpLo_ }
 
typedef MatrixType_ MatrixType
 
typedef SolverBase< LLTBase
 
typedef internal::LLT_Traits< MatrixType, UpLoTraits
 
- Public Types inherited from Eigen::SolverBase< LLT< MatrixType_, UpLo_ > >
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 ConstTransposeReturnTypeAdjointReturnType
 
- Public Types inherited from Eigen::EigenBase< Derived >
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 >
LLTcompute (const EigenBase< InputType > &matrix)
 
RealScalar rcond () const
 
const MatrixTypematrixLLT () const
 
MatrixType reconstructedMatrix () const
 
ComputationInfo info () const
 Reports whether previous computation was successful. More...
 
const LLTadjoint () const EIGEN_NOEXCEPT
 
EIGEN_CONSTEXPR Index rows () const EIGEN_NOEXCEPT
 
EIGEN_CONSTEXPR Index cols () const EIGEN_NOEXCEPT
 
template<typename VectorType >
LLTrankUpdate (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)
 
- Public Member Functions inherited from Eigen::SolverBase< LLT< MatrixType_, UpLo_ > >
 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
 
- Public Member Functions inherited from Eigen::EigenBase< Derived >
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

- Protected Member Functions inherited from Eigen::SolverBase< LLT< MatrixType_, UpLo_ > >
void _check_solve_assertion (const Rhs &b) const
 

Detailed Description

template<typename MatrixType_, int UpLo_>
class Eigen::LLT< MatrixType_, UpLo_ >

Standard Cholesky decomposition (LL^T) of a matrix and associated features.

Template Parameters
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:

MatrixXd A(3, 3);
A << 4, -1, 2, -1, 6, 0, 2, 0, 5;
cout << "The matrix A is" << endl << A << endl;
LLT<MatrixXd> lltOfA(A); // compute the Cholesky decomposition of A
MatrixXd L = lltOfA.matrixL(); // retrieve factor L in the decomposition
// The previous two lines can also be written as "L = A.llt().matrixL()"
cout << "The Cholesky factor L is" << endl << L << endl;
cout << "To check this, let us compute L * L.transpose()" << endl;
cout << L * L.transpose() << endl;
cout << "This should equal the matrix A" << endl;
MatrixXd L
Definition: LLT_example.cpp:6
A<< 4, -1, 2, -1, 6, 0, 2, 0, 5;cout<< "The matrix A is"<< endl<< A<< endl;LLT< MatrixXd > lltOfA(A)
Matrix< SCALARA, Dynamic, Dynamic, opt_A > A
Definition: bench_gemm.cpp:47

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.

See also
MatrixBase::llt(), SelfAdjointView::llt(), class LDLT

Member Typedef Documentation

◆ Base

template<typename MatrixType_ , int UpLo_>
typedef SolverBase<LLT> Eigen::LLT< MatrixType_, UpLo_ >::Base

◆ MatrixType

template<typename MatrixType_ , int UpLo_>
typedef MatrixType_ Eigen::LLT< MatrixType_, UpLo_ >::MatrixType

◆ Traits

template<typename MatrixType_ , int UpLo_>
typedef internal::LLT_Traits<MatrixType, UpLo> Eigen::LLT< MatrixType_, UpLo_ >::Traits

Member Enumeration Documentation

◆ anonymous enum

template<typename MatrixType_ , int UpLo_>
anonymous enum
Enumerator
MaxColsAtCompileTime 
77 { MaxColsAtCompileTime = MatrixType::MaxColsAtCompileTime };
@ MaxColsAtCompileTime
Definition: LLT.h:77

◆ anonymous enum

template<typename MatrixType_ , int UpLo_>
anonymous enum
Enumerator
PacketSize 
AlignmentMask 
UpLo 
@ UpLo
Definition: LLT.h:79
@ PacketSize
Definition: LLT.h:79
@ AlignmentMask
Definition: LLT.h:79
return int(ret)+1
@ size
Definition: GenericPacketMath.h:113

Constructor & Destructor Documentation

◆ LLT() [1/4]

template<typename MatrixType_ , int UpLo_>
Eigen::LLT< MatrixType_, UpLo_ >::LLT ( )
inline

Default Constructor.

The default constructor is useful in cases in which the user intends to perform decompositions via LLT::compute(const MatrixType&).

89 : m_matrix(), m_isInitialized(false) {}
MatrixType m_matrix
Definition: LLT.h:208
bool m_isInitialized
Definition: LLT.h:210

◆ LLT() [2/4]

template<typename MatrixType_ , int UpLo_>
Eigen::LLT< MatrixType_, UpLo_ >::LLT ( Index  size)
inlineexplicit

Default Constructor with memory preallocation.

Like the default constructor but with preallocation of the internal data according to the specified problem size.

See also
LLT()
97 : m_matrix(size, size), m_isInitialized(false) {}
EIGEN_DEVICE_FUNC EIGEN_CONSTEXPR Index size() const EIGEN_NOEXCEPT
Definition: EigenBase.h:64

◆ LLT() [3/4]

template<typename MatrixType_ , int UpLo_>
template<typename InputType >
Eigen::LLT< MatrixType_, UpLo_ >::LLT ( const EigenBase< InputType > &  matrix)
inlineexplicit
100  : m_matrix(matrix.rows(), matrix.cols()), m_isInitialized(false) {
101  compute(matrix.derived());
102  }
LLT & compute(const EigenBase< InputType > &matrix)
Eigen::Map< Eigen::Matrix< T, Eigen::Dynamic, Eigen::Dynamic, Eigen::ColMajor >, 0, Eigen::OuterStride<> > matrix(T *data, int rows, int cols, int stride)
Definition: common.h:85

References Eigen::LLT< MatrixType_, UpLo_ >::compute(), and matrix().

◆ LLT() [4/4]

template<typename MatrixType_ , int UpLo_>
template<typename InputType >
Eigen::LLT< MatrixType_, UpLo_ >::LLT ( EigenBase< InputType > &  matrix)
inlineexplicit

Constructs a LLT factorization from a given matrix.

This overloaded constructor is provided for inplace decomposition when MatrixType is a Eigen::Ref.

See also
LLT(const EigenBase&)
112  : m_matrix(matrix.derived()), m_isInitialized(false) {
113  compute(matrix.derived());
114  }

References Eigen::LLT< MatrixType_, UpLo_ >::compute(), and matrix().

Member Function Documentation

◆ _solve_impl()

template<typename MatrixType_ , int UpLo_>
template<typename RhsType , typename DstType >
void Eigen::LLT< MatrixType_, UpLo_ >::_solve_impl ( const RhsType &  rhs,
DstType &  dst 
) const
448  {
449  _solve_impl_transposed<true>(rhs, dst);
450 }

◆ _solve_impl_transposed()

template<typename MatrixType_ , int UpLo_>
template<bool Conjugate, typename RhsType , typename DstType >
void Eigen::LLT< MatrixType_, UpLo_ >::_solve_impl_transposed ( const RhsType &  rhs,
DstType &  dst 
) const
454  {
455  dst = rhs;
456 
457  matrixL().template conjugateIf<!Conjugate>().solveInPlace(dst);
458  matrixU().template conjugateIf<!Conjugate>().solveInPlace(dst);
459 }
Traits::MatrixU matrixU() const
Definition: LLT.h:117
Traits::MatrixL matrixL() const
Definition: LLT.h:123

References Eigen::LLT< MatrixType_, UpLo_ >::matrixL(), and Eigen::LLT< MatrixType_, UpLo_ >::matrixU().

◆ adjoint()

template<typename MatrixType_ , int UpLo_>
const LLT& Eigen::LLT< MatrixType_, UpLo_ >::adjoint ( ) const
inline
Returns
the adjoint of *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:

x = decomposition.adjoint().solve(b)
Scalar * b
Definition: benchVecAdd.cpp:17
list x
Definition: plotDoE.py:28
185 { return *this; }

◆ cols()

◆ compute() [1/2]

template<typename MatrixType_ , int UpLo_>
template<typename InputType >
LLT<MatrixType, UpLo_>& Eigen::LLT< MatrixType_, UpLo_ >::compute ( const EigenBase< InputType > &  a)

Computes / recomputes the Cholesky decomposition A = LL^* = U^*U of matrix

Returns
a reference to *this

Example:

#include <iostream>
#include <Eigen/Dense>
int main() {
Eigen::Matrix2f A, b;
A << 2, -1, -1, 3;
b << 1, 2, 3, 1;
std::cout << "Here is the matrix A:\n" << A << std::endl;
std::cout << "Here is the right hand side b:\n" << b << std::endl;
std::cout << "Computing LLT decomposition..." << std::endl;
llt.compute(A);
std::cout << "The solution is:\n" << llt.solve(b) << std::endl;
A(1, 1)++;
std::cout << "The matrix A is now:\n" << A << std::endl;
std::cout << "Computing LLT decomposition..." << std::endl;
llt.compute(A);
std::cout << "The solution is now:\n" << llt.solve(b) << std::endl;
}
int main(int argc, char *argv[])
Definition: T_protectiveWall.cpp:194
Standard Cholesky decomposition (LL^T) of a matrix and associated features.
Definition: LLT.h:70
EIGEN_DONT_INLINE void llt(const Mat &A, const Mat &B, Mat &C)
Definition: llt.cpp:4

Output:

 
399  {
400  eigen_assert(a.rows() == a.cols());
401  const Index size = a.rows();
402  m_matrix.resize(size, size);
403  if (!internal::is_same_dense(m_matrix, a.derived())) m_matrix = a.derived();
404 
405  // Compute matrix L1 norm = max abs column sum.
406  m_l1_norm = RealScalar(0);
407  // TODO move this code to SelfAdjointView
408  for (Index col = 0; col < size; ++col) {
409  RealScalar abs_col_sum;
410  if (UpLo_ == Lower)
411  abs_col_sum =
412  m_matrix.col(col).tail(size - col).template lpNorm<1>() + m_matrix.row(col).head(col).template lpNorm<1>();
413  else
414  abs_col_sum =
415  m_matrix.col(col).head(col).template lpNorm<1>() + m_matrix.row(col).tail(size - col).template lpNorm<1>();
416  if (abs_col_sum > m_l1_norm) m_l1_norm = abs_col_sum;
417  }
418 
419  m_isInitialized = true;
420  bool ok = Traits::inplace_decomposition(m_matrix);
421  m_info = ok ? Success : NumericalIssue;
422 
423  return *this;
424 }
#define eigen_assert(x)
Definition: Macros.h:910
m col(1)
NumTraits< Scalar >::Real RealScalar
Definition: bench_gemm.cpp:46
RealScalar m_l1_norm
Definition: LLT.h:209
ComputationInfo m_info
Definition: LLT.h:211
@ Lower
Definition: Constants.h:211
@ NumericalIssue
Definition: Constants.h:442
@ Success
Definition: Constants.h:440
const Scalar * a
Definition: level2_cplx_impl.h:32
EIGEN_DEVICE_FUNC bool is_same_dense(const T1 &mat1, const T2 &mat2, std::enable_if_t< possibly_same_dense< T1, T2 >::value > *=0)
Definition: XprHelper.h:869
Eigen::Index Index
The interface type of indices.
Definition: EigenBase.h:43

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.

◆ compute() [2/2]

template<typename MatrixType_ , int UpLo_>
template<typename InputType >
LLT& Eigen::LLT< MatrixType_, UpLo_ >::compute ( const EigenBase< InputType > &  matrix)

◆ info()

template<typename MatrixType_ , int UpLo_>
ComputationInfo Eigen::LLT< MatrixType_, UpLo_ >::info ( ) const
inline

Reports whether previous computation was successful.

Returns
Success if computation was successful, NumericalIssue if the matrix.appears not to be positive definite.
174  {
175  eigen_assert(m_isInitialized && "LLT is not initialized.");
176  return m_info;
177  }

References eigen_assert, Eigen::LLT< MatrixType_, UpLo_ >::m_info, and Eigen::LLT< MatrixType_, UpLo_ >::m_isInitialized.

Referenced by cholesky().

◆ matrixL()

template<typename MatrixType_ , int UpLo_>
Traits::MatrixL Eigen::LLT< MatrixType_, UpLo_ >::matrixL ( ) const
inline

◆ matrixLLT()

template<typename MatrixType_ , int UpLo_>
const MatrixType& Eigen::LLT< MatrixType_, UpLo_ >::matrixLLT ( ) const
inline
Returns
the LLT decomposition matrix

TODO: document the storage layout

162  {
163  eigen_assert(m_isInitialized && "LLT is not initialized.");
164  return m_matrix;
165  }

References eigen_assert, Eigen::LLT< MatrixType_, UpLo_ >::m_isInitialized, and Eigen::LLT< MatrixType_, UpLo_ >::m_matrix.

◆ matrixU()

template<typename MatrixType_ , int UpLo_>
Traits::MatrixU Eigen::LLT< MatrixType_, UpLo_ >::matrixU ( ) const
inline

◆ rankUpdate() [1/2]

template<typename MatrixType_ , int UpLo_>
template<typename VectorType >
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.

433  {
435  eigen_assert(v.size() == m_matrix.cols());
437  if (internal::llt_inplace<typename MatrixType::Scalar, UpLo>::rankUpdate(m_matrix, v, sigma) >= 0)
439  else
440  m_info = Success;
441 
442  return *this;
443 }
Array< int, Dynamic, 1 > v
Definition: Array_initializer_list_vector_cxx11.cpp:1
#define EIGEN_STATIC_ASSERT_VECTOR_ONLY(TYPE)
Definition: StaticAssert.h:36
int sigma
Definition: calibrate.py:179
Definition: fft_test_shared.h:66

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.

◆ rankUpdate() [2/2]

template<typename MatrixType_ , int UpLo_>
template<typename VectorType >
LLT& Eigen::LLT< MatrixType_, UpLo_ >::rankUpdate ( const VectorType vec,
const RealScalar sigma = 1 
)

◆ rcond()

template<typename MatrixType_ , int UpLo_>
RealScalar Eigen::LLT< MatrixType_, UpLo_ >::rcond ( ) const
inline
Returns
an estimate of the reciprocal condition number of the matrix of which *this is the Cholesky decomposition.
152  {
153  eigen_assert(m_isInitialized && "LLT is not initialized.");
154  eigen_assert(m_info == Success && "LLT failed because matrix appears to be negative");
156  }
Decomposition::RealScalar rcond_estimate_helper(typename Decomposition::RealScalar matrix_norm, const Decomposition &dec)
Reciprocal condition number estimator.
Definition: ConditionEstimator.h:157

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().

◆ reconstructedMatrix()

template<typename MatrixType , int UpLo_>
MatrixType Eigen::LLT< MatrixType, UpLo_ >::reconstructedMatrix
Returns
the matrix represented by the decomposition, i.e., it returns the product: L L^*. This function is provided for debug purpose.
488  {
489  eigen_assert(m_isInitialized && "LLT is not initialized.");
490  return matrixL() * matrixL().adjoint().toDenseMatrix();
491 }

References eigen_assert, Eigen::LLT< MatrixType_, UpLo_ >::m_isInitialized, and Eigen::LLT< MatrixType_, UpLo_ >::matrixL().

Referenced by cholesky(), and cholesky_cplx().

◆ rows()

◆ solveInPlace()

template<typename MatrixType , int UpLo_>
template<typename Derived >
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().

Parameters
bAndXrepresents 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.

Warning
The parameter is only marked 'const' to make the C++ compiler accept a temporary expression here. This function will const_cast it, so constness isn't honored here.
See also
LLT::solve(), MatrixBase::llt()
477  {
478  eigen_assert(m_isInitialized && "LLT is not initialized.");
479  eigen_assert(m_matrix.rows() == bAndX.rows());
480  matrixL().solveInPlace(bAndX);
481  matrixU().solveInPlace(bAndX);
482 }

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().

Friends And Related Function Documentation

◆ SolverBase< LLT >

template<typename MatrixType_ , int UpLo_>
friend class SolverBase< LLT >
friend

Member Data Documentation

◆ m_info

◆ m_isInitialized

◆ m_l1_norm

template<typename MatrixType_ , int UpLo_>
RealScalar Eigen::LLT< MatrixType_, UpLo_ >::m_l1_norm
protected

◆ m_matrix


The documentation for this class was generated from the following file: