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

Searching...

No Matches

#include <Eigen/src/Cholesky/LLT.h>

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:

// SPDX-FileCopyrightText: The Eigen Authors
// SPDX-License-Identifier: MPL-2.0
 
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;

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

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 upper part (or the strict lower part when UpLo_ is Upper) does not have to store correct values.

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

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

Public Member Functions
RealScalar	absDeterminant () const

const LLT &	adjoint () const noexcept

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

Scalar	determinant () const

ComputationInfo	info () const
	Reports whether previous computation was successful.

Inverse< LLT >	inverse () const

	LLT ()
	Default Constructor.

template<typename InputType>
	LLT (EigenBase< InputType > &matrix)
	Constructs a LLT factorization from a given matrix.

	LLT (Index size)
	Default Constructor with memory preallocation.

RealScalar	logAbsDeterminant () const

Traits::MatrixL	matrixL () const

const MatrixType &	matrixLLT () const

Traits::MatrixU	matrixU () const

template<typename VectorType>
LLT< MatrixType_, UpLo_ > &	rankUpdate (const VectorType &v, const RealScalar &sigma)

RealScalar	rcond () const

MatrixType	reconstructedMatrix () const

Scalar	signDeterminant () const

template<typename Rhs>
Solve< LLT, Rhs >	solve (const MatrixBase< Rhs > &b) const

Public Member Functions inherited from Eigen::SolverBase< LLT< MatrixType_, UpLo_ > >
const AdjointReturnType	adjoint () const

constexpr LLT< MatrixType_, UpLo_ > &	derived ()

constexpr const LLT< MatrixType_, UpLo_ > &	derived () const

Solve< LLT< MatrixType_, UpLo_ >, Rhs >	solve (const MatrixBase< Rhs > &b) const

	SolverBase ()=default

const ConstTransposeReturnType	transpose () const

Public Member Functions inherited from Eigen::EigenBase< LLT< MatrixType_, UpLo_ > >
constexpr Index	cols () const noexcept

constexpr LLT< MatrixType_, UpLo_ > &	derived ()

constexpr const LLT< MatrixType_, UpLo_ > &	derived () const

constexpr Index	rows () const noexcept

constexpr Index	size () const noexcept

Additional Inherited Members
Public Types inherited from Eigen::EigenBase< LLT< MatrixType_, UpLo_ > >
using	Index
	The interface type of indices.

Constructor & Destructor Documentation

◆ LLT() [1/3]

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

◆ LLT() [2/3]

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

◆ LLT() [3/3]

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&)

Member Function Documentation

◆ absDeterminant()

template<typename MatrixType_, int UpLo_>

LLT< MatrixType_, UpLo_ >::RealScalar Eigen::LLT< MatrixType_, UpLo_ >::absDeterminant ( ) const

Returns: the absolute value of the determinant of the matrix of which *this is the Cholesky decomposition.

It has only linear complexity (that is, O(n) where n is the dimension of the square matrix) as the Cholesky decomposition has already been computed.

Note: The decomposed matrix is positive definite, so this is the determinant itself.

Warning: a determinant can be very big or small, so for matrices of large enough dimension, there is a risk of overflow/underflow. One way to work around that is to use logAbsDeterminant() instead.

Precondition: info() returns Success. A failed factorization does not represent the input matrix.

See also: determinant(), logAbsDeterminant(), signDeterminant(), MatrixBase::determinant()

◆ adjoint()

template<typename MatrixType_, int UpLo_>

const LLT & Eigen::LLT< MatrixType_, UpLo_ >::adjoint ( ) const

inlinenoexcept

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)

◆ compute()

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:

// SPDX-FileCopyrightText: The Eigen Authors
// SPDX-License-Identifier: MPL-2.0
 
#include <iostream>
#include <Eigen/Dense>
 
int main() {
  Eigen::Matrix2f A, b;
  Eigen::LLT<Eigen::Matrix2f> llt;
  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;
}

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

◆ determinant()

template<typename MatrixType_, int UpLo_>

LLT< MatrixType_, UpLo_ >::Scalar Eigen::LLT< MatrixType_, UpLo_ >::determinant ( ) const

Returns: the determinant of the matrix of which *this is the Cholesky decomposition.

It has only linear complexity (that is, O(n) where n is the dimension of the square matrix) as the Cholesky decomposition has already been computed.

Warning: a determinant can be very big or small, so for matrices of large enough dimension, there is a risk of overflow/underflow. One way to work around that is to use logAbsDeterminant() instead.

Precondition: info() returns Success. A failed factorization does not represent the input matrix.

See also: absDeterminant(), logAbsDeterminant(), signDeterminant(), MatrixBase::determinant()

◆ 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.

◆ inverse()

template<typename MatrixType_, int UpLo_>

Inverse< LLT > Eigen::LLT< MatrixType_, UpLo_ >::inverse ( ) const

inline

Returns: the inverse of the matrix of which *this is the Cholesky decomposition.

The result is computed as \( A^{-1} = L^{-*} L^{-1} \) by inverting the stored factor in place and squaring it, the LAPACK *POTRI sequence, for 2n^3/3 flops against the 2n^3 of solving with an explicit identity right hand side. Beyond the destination, only a block-sized scratch panel of at most 128x128 coefficients is used, whatever the matrix size.

Fewer flops is not fewer seconds at every size: for a real scalar below EIGEN_LLT_INVERSE_POTRI_THRESHOLD (32, or 256 where AVX-512 is enabled) the POTRI sequence runs in its unblocked kernels while the solve against an identity runs its 3x at the blocked TRSM rate, so that is what this method does there. Complex scalars always take the POTRI path. The result is exactly self-adjoint either way: one triangle is computed and mirrored onto the other.

An in-place decomposition (see the class documentation) may overwrite its own factor with the result, as in storage = llt.inverse(); like any other write to the referenced matrix, that leaves the decomposition unusable afterwards.

The matrix must be positive definite, that is, info() must be Success.

See also: solve(), MatrixBase::inverse()

◆ logAbsDeterminant()

template<typename MatrixType_, int UpLo_>

LLT< MatrixType_, UpLo_ >::RealScalar Eigen::LLT< MatrixType_, UpLo_ >::logAbsDeterminant ( ) const

Returns: the natural log of the absolute value of the determinant of the matrix of which *this is the Cholesky decomposition.

It has only linear complexity (that is, O(n) where n is the dimension of the square matrix) as the Cholesky decomposition has already been computed.

Note: This method is useful to work around the risk of overflow/underflow that's inherent to determinant computation.

Precondition: info() returns Success. A failed factorization does not represent the input matrix.

See also: determinant(), absDeterminant(), signDeterminant(), MatrixBase::determinant()

◆ matrixL()

template<typename MatrixType_, int UpLo_>

Traits::MatrixL Eigen::LLT< MatrixType_, UpLo_ >::matrixL ( ) const

inline

Returns: a view of the lower triangular matrix L

◆ 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

◆ matrixU()

template<typename MatrixType_, int UpLo_>

Traits::MatrixU Eigen::LLT< MatrixType_, UpLo_ >::matrixU ( ) const

inline

Returns: a view of the upper triangular matrix U

◆ rankUpdate()

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

◆ 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.

◆ reconstructedMatrix()

template<typename MatrixType, int UpLo_>

MatrixType Eigen::LLT< MatrixType, UpLo_ >::reconstructedMatrix ( ) const

Returns: the matrix represented by the decomposition, i.e., it returns the product: L L^*. This function is provided for debug purpose.

◆ signDeterminant()

template<typename MatrixType_, int UpLo_>

LLT< MatrixType_, UpLo_ >::Scalar Eigen::LLT< MatrixType_, UpLo_ >::signDeterminant ( ) const

Returns: the sign of the determinant of the matrix of which *this is the Cholesky decomposition, which is 1 since that matrix is positive definite.

This method is provided for compatibility with the other decompositions, thus enabling generic code.

Precondition: info() returns Success. A failed factorization does not represent the input matrix.

See also: determinant(), absDeterminant(), logAbsDeterminant(), MatrixBase::determinant()

◆ solve()

template<typename MatrixType_, int UpLo_>

template<typename Rhs>

Solve< LLT, Rhs > Eigen::LLT< MatrixType_, UpLo_ >::solve ( const MatrixBase< Rhs > & b ) const

inline

Returns: the solution x of \( A x = b \) using the current decomposition of A.

Since this LLT class assumes anyway that the matrix A is invertible, the solution theoretically exists and is unique regardless of b.

Example:

// SPDX-FileCopyrightText: The Eigen Authors
// SPDX-License-Identifier: MPL-2.0
 
typedef Matrix<float, Dynamic, 2> DataMatrix;
// let's generate some samples on the 3D plane of equation z = 2x+3y (with some noise)
DataMatrix samples = DataMatrix::Random(12, 2);
VectorXf elevations = 2 * samples.col(0) + 3 * samples.col(1) + VectorXf::Random(12) * 0.1;
// and let's solve samples * [x y]^T = elevations in least square sense:
Matrix<float, 2, 1> xy = (samples.adjoint() * samples).llt().solve((samples.adjoint() * elevations));
cout << xy << endl;

Output:

2.02
2.97

See also: solveInPlace(), MatrixBase::llt(), SelfAdjointView::llt()

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

/builds/libeigen/eigen/Eigen/src/Cholesky/LLT.h

Detailed Description

Public Member Functions

Additional Inherited Members

Constructor & Destructor Documentation

◆ LLT() [1/3]

◆ LLT() [2/3]

◆ LLT() [3/3]

Member Function Documentation

◆ absDeterminant()

◆ adjoint()

◆ compute()

◆ determinant()

◆ info()

◆ inverse()

◆ logAbsDeterminant()

◆ matrixL()

◆ matrixLLT()

◆ matrixU()

◆ rankUpdate()

◆ rcond()

◆ reconstructedMatrix()

◆ signDeterminant()

◆ solve()