MOAB  4.9.3pre
Eigen::LDLT< _MatrixType, _UpLo > Class Template Reference

Robust Cholesky decomposition of a matrix with pivoting. More...

#include <LDLT.h>

Collaboration diagram for Eigen::LDLT< _MatrixType, _UpLo >:

List of all members.

Public Types

enum  {
  RowsAtCompileTime = MatrixType::RowsAtCompileTime, ColsAtCompileTime = MatrixType::ColsAtCompileTime, Options = MatrixType::Options & ~RowMajorBit, MaxRowsAtCompileTime = MatrixType::MaxRowsAtCompileTime,
  MaxColsAtCompileTime = MatrixType::MaxColsAtCompileTime, UpLo = _UpLo
}
typedef _MatrixType MatrixType
typedef MatrixType::Scalar Scalar
typedef NumTraits< typename
MatrixType::Scalar >::Real 
RealScalar
typedef Eigen::Index Index
typedef MatrixType::StorageIndex StorageIndex
typedef Matrix< Scalar,
RowsAtCompileTime, 1, Options,
MaxRowsAtCompileTime, 1 > 
TmpMatrixType
typedef Transpositions
< RowsAtCompileTime,
MaxRowsAtCompileTime
TranspositionType
typedef PermutationMatrix
< RowsAtCompileTime,
MaxRowsAtCompileTime
PermutationType
typedef internal::LDLT_Traits
< MatrixType, UpLo
Traits

Public Member Functions

 LDLT ()
 Default Constructor.
 LDLT (Index size)
 Default Constructor with memory preallocation.
template<typename InputType >
 LDLT (const EigenBase< InputType > &matrix)
 Constructor with decomposition.
void setZero ()
Traits::MatrixU matrixU () const
Traits::MatrixL matrixL () const
const TranspositionTypetranspositionsP () const
Diagonal< const MatrixTypevectorD () const
bool isPositive () const
bool isNegative (void) const
template<typename Rhs >
const Solve< LDLT, Rhs > solve (const MatrixBase< Rhs > &b) const
template<typename Derived >
bool solveInPlace (MatrixBase< Derived > &bAndX) const
template<typename InputType >
LDLTcompute (const EigenBase< InputType > &matrix)
template<typename Derived >
LDLTrankUpdate (const MatrixBase< Derived > &w, const RealScalar &alpha=1)
const MatrixTypematrixLDLT () const
MatrixType reconstructedMatrix () const
Index rows () const
Index cols () const
ComputationInfo info () const
 Reports whether previous computation was successful.
template<typename RhsType , typename DstType >
EIGEN_DEVICE_FUNC void _solve_impl (const RhsType &rhs, DstType &dst) const
template<typename Derived >
LDLT< MatrixType, _UpLo > & rankUpdate (const MatrixBase< Derived > &w, const typename LDLT< MatrixType, _UpLo >::RealScalar &sigma)

Static Protected Member Functions

static void check_template_parameters ()

Protected Attributes

MatrixType m_matrix
TranspositionType m_transpositions
TmpMatrixType m_temporary
internal::SignMatrix m_sign
bool m_isInitialized

Detailed Description

template<typename _MatrixType, int _UpLo>
class Eigen::LDLT< _MatrixType, _UpLo >

Robust Cholesky decomposition of a matrix with pivoting.

Template Parameters:
_MatrixTypethe type of the matrix of which to compute the LDL^T Cholesky decomposition
_UpLothe triangular part that will be used for the decompositon: Lower (default) or Upper. The other triangular part won't be read.

Perform a robust Cholesky decomposition of a positive semidefinite or negative semidefinite matrix $ A $ such that $ A = P^TLDL^*P $, where P is a permutation matrix, L is lower triangular with a unit diagonal and D is a diagonal matrix.

The decomposition uses pivoting to ensure stability, so that L will have zeros in the bottom right rank(A) - n submatrix. Avoiding the square root on D also stabilizes the computation.

Remember that Cholesky decompositions are not rank-revealing. Also, do not use a Cholesky decomposition to determine whether a system of equations has a solution.

See also:
MatrixBase::ldlt(), SelfAdjointView::ldlt(), class LLT

Definition at line 48 of file LDLT.h.


Member Typedef Documentation

template<typename _MatrixType, int _UpLo>
typedef Eigen::Index Eigen::LDLT< _MatrixType, _UpLo >::Index
Deprecated:
since Eigen 3.3

Definition at line 62 of file LDLT.h.

template<typename _MatrixType, int _UpLo>
typedef _MatrixType Eigen::LDLT< _MatrixType, _UpLo >::MatrixType

Definition at line 51 of file LDLT.h.

template<typename _MatrixType, int _UpLo>
typedef PermutationMatrix<RowsAtCompileTime, MaxRowsAtCompileTime> Eigen::LDLT< _MatrixType, _UpLo >::PermutationType

Definition at line 67 of file LDLT.h.

template<typename _MatrixType, int _UpLo>
typedef NumTraits<typename MatrixType::Scalar>::Real Eigen::LDLT< _MatrixType, _UpLo >::RealScalar

Definition at line 61 of file LDLT.h.

template<typename _MatrixType, int _UpLo>
typedef MatrixType::Scalar Eigen::LDLT< _MatrixType, _UpLo >::Scalar

Definition at line 60 of file LDLT.h.

template<typename _MatrixType, int _UpLo>
typedef MatrixType::StorageIndex Eigen::LDLT< _MatrixType, _UpLo >::StorageIndex

Definition at line 63 of file LDLT.h.

template<typename _MatrixType, int _UpLo>
typedef Matrix<Scalar, RowsAtCompileTime, 1, Options, MaxRowsAtCompileTime, 1> Eigen::LDLT< _MatrixType, _UpLo >::TmpMatrixType

Definition at line 64 of file LDLT.h.

template<typename _MatrixType, int _UpLo>
typedef internal::LDLT_Traits<MatrixType,UpLo> Eigen::LDLT< _MatrixType, _UpLo >::Traits

Definition at line 69 of file LDLT.h.

template<typename _MatrixType, int _UpLo>
typedef Transpositions<RowsAtCompileTime, MaxRowsAtCompileTime> Eigen::LDLT< _MatrixType, _UpLo >::TranspositionType

Definition at line 66 of file LDLT.h.


Member Enumeration Documentation

template<typename _MatrixType, int _UpLo>
anonymous enum
Enumerator:
RowsAtCompileTime 
ColsAtCompileTime 
Options 
MaxRowsAtCompileTime 
MaxColsAtCompileTime 
UpLo 

Definition at line 52 of file LDLT.h.


Constructor & Destructor Documentation

template<typename _MatrixType, int _UpLo>
Eigen::LDLT< _MatrixType, _UpLo >::LDLT ( ) [inline]

Default Constructor.

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

Definition at line 76 of file LDLT.h.

template<typename _MatrixType, int _UpLo>
Eigen::LDLT< _MatrixType, _UpLo >::LDLT ( Index  size) [inline, explicit]

Default Constructor with memory preallocation.

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

See also:
LDLT()

Definition at line 89 of file LDLT.h.

template<typename _MatrixType, int _UpLo>
template<typename InputType >
Eigen::LDLT< _MatrixType, _UpLo >::LDLT ( const EigenBase< InputType > &  matrix) [inline, explicit]

Constructor with decomposition.

This calculates the decomposition for the input matrix.

See also:
LDLT(Index size)

Definition at line 103 of file LDLT.h.

      : m_matrix(matrix.rows(), matrix.cols()),
        m_transpositions(matrix.rows()),
        m_temporary(matrix.rows()),
        m_sign(internal::ZeroSign),
        m_isInitialized(false)
    {
      compute(matrix.derived());
    }

Member Function Documentation

template<typename _MatrixType , int _UpLo>
template<typename RhsType , typename DstType >
void Eigen::LDLT< _MatrixType, _UpLo >::_solve_impl ( const RhsType &  rhs,
DstType &  dst 
) const

Definition at line 488 of file LDLT.h.

{
  eigen_assert(rhs.rows() == rows());
  // dst = P b
  dst = m_transpositions * rhs;

  // dst = L^-1 (P b)
  matrixL().solveInPlace(dst);

  // dst = D^-1 (L^-1 P b)
  // more precisely, use pseudo-inverse of D (see bug 241)
  using std::abs;
  const typename Diagonal<const MatrixType>::RealReturnType vecD(vectorD());
  // In some previous versions, tolerance was set to the max of 1/highest and the maximal diagonal entry * epsilon
  // as motivated by LAPACK's xGELSS:
  // RealScalar tolerance = numext::maxi(vecD.array().abs().maxCoeff() * NumTraits<RealScalar>::epsilon(),RealScalar(1) / NumTraits<RealScalar>::highest());
  // However, LDLT is not rank revealing, and so adjusting the tolerance wrt to the highest
  // diagonal element is not well justified and leads to numerical issues in some cases.
  // Moreover, Lapack's xSYTRS routines use 0 for the tolerance.
  RealScalar tolerance = RealScalar(1) / NumTraits<RealScalar>::highest();
  
  for (Index i = 0; i < vecD.size(); ++i)
  {
    if(abs(vecD(i)) > tolerance)
      dst.row(i) /= vecD(i);
    else
      dst.row(i).setZero();
  }

  // dst = L^-T (D^-1 L^-1 P b)
  matrixU().solveInPlace(dst);

  // dst = P^-1 (L^-T D^-1 L^-1 P b) = A^-1 b
  dst = m_transpositions.transpose() * dst;
}
template<typename _MatrixType, int _UpLo>
static void Eigen::LDLT< _MatrixType, _UpLo >::check_template_parameters ( ) [inline, static, protected]

Definition at line 232 of file LDLT.h.

template<typename _MatrixType, int _UpLo>
Index Eigen::LDLT< _MatrixType, _UpLo >::cols ( ) const [inline]

Definition at line 211 of file LDLT.h.

{ return m_matrix.cols(); }
template<typename MatrixType , int _UpLo>
template<typename InputType >
LDLT< MatrixType, _UpLo > & Eigen::LDLT< MatrixType, _UpLo >::compute ( const EigenBase< InputType > &  a)

Compute / recompute the LDLT decomposition A = L D L^* = U^* D U of matrix

Definition at line 433 of file LDLT.h.

{
  check_template_parameters();
  
  eigen_assert(a.rows()==a.cols());
  const Index size = a.rows();

  m_matrix = a.derived();

  m_transpositions.resize(size);
  m_isInitialized = false;
  m_temporary.resize(size);
  m_sign = internal::ZeroSign;

  internal::ldlt_inplace<UpLo>::unblocked(m_matrix, m_transpositions, m_temporary, m_sign);

  m_isInitialized = true;
  return *this;
}
template<typename _MatrixType, int _UpLo>
ComputationInfo Eigen::LDLT< _MatrixType, _UpLo >::info ( ) const [inline]

Reports whether previous computation was successful.

Returns:
Success if computation was succesful, NumericalIssue if the matrix.appears to be negative.

Definition at line 218 of file LDLT.h.

    {
      eigen_assert(m_isInitialized && "LDLT is not initialized.");
      return Success;
    }
template<typename _MatrixType, int _UpLo>
bool Eigen::LDLT< _MatrixType, _UpLo >::isNegative ( void  ) const [inline]
Returns:
true if the matrix is negative (semidefinite)

Definition at line 158 of file LDLT.h.

    {
      eigen_assert(m_isInitialized && "LDLT is not initialized.");
      return m_sign == internal::NegativeSemiDef || m_sign == internal::ZeroSign;
    }
template<typename _MatrixType, int _UpLo>
bool Eigen::LDLT< _MatrixType, _UpLo >::isPositive ( ) const [inline]
Returns:
true if the matrix is positive (semidefinite)

Definition at line 151 of file LDLT.h.

    {
      eigen_assert(m_isInitialized && "LDLT is not initialized.");
      return m_sign == internal::PositiveSemiDef || m_sign == internal::ZeroSign;
    }
template<typename _MatrixType, int _UpLo>
Traits::MatrixL Eigen::LDLT< _MatrixType, _UpLo >::matrixL ( ) const [inline]
Returns:
a view of the lower triangular matrix L

Definition at line 129 of file LDLT.h.

    {
      eigen_assert(m_isInitialized && "LDLT is not initialized.");
      return Traits::getL(m_matrix);
    }
template<typename _MatrixType, int _UpLo>
const MatrixType& Eigen::LDLT< _MatrixType, _UpLo >::matrixLDLT ( ) const [inline]
Returns:
the internal LDLT decomposition matrix

TODO: document the storage layout

Definition at line 202 of file LDLT.h.

    {
      eigen_assert(m_isInitialized && "LDLT is not initialized.");
      return m_matrix;
    }
template<typename _MatrixType, int _UpLo>
Traits::MatrixU Eigen::LDLT< _MatrixType, _UpLo >::matrixU ( ) const [inline]
Returns:
a view of the upper triangular matrix U

Definition at line 122 of file LDLT.h.

    {
      eigen_assert(m_isInitialized && "LDLT is not initialized.");
      return Traits::getU(m_matrix);
    }
template<typename _MatrixType, int _UpLo>
template<typename Derived >
LDLT& Eigen::LDLT< _MatrixType, _UpLo >::rankUpdate ( const MatrixBase< Derived > &  w,
const RealScalar alpha = 1 
)
template<typename _MatrixType, int _UpLo>
template<typename Derived >
LDLT<MatrixType,_UpLo>& Eigen::LDLT< _MatrixType, _UpLo >::rankUpdate ( const MatrixBase< Derived > &  w,
const typename LDLT< MatrixType, _UpLo >::RealScalar sigma 
)

Update the LDLT decomposition: given A = L D L^T, efficiently compute the decomposition of A + sigma w w^T.

Parameters:
wa vector to be incorporated into the decomposition.
sigmaa scalar, +1 for updates and -1 for "downdates," which correspond to removing previously-added column vectors. Optional; default value is +1.
See also:
setZero()

Definition at line 460 of file LDLT.h.

{
  typedef typename TranspositionType::StorageIndex IndexType;
  const Index size = w.rows();
  if (m_isInitialized)
  {
    eigen_assert(m_matrix.rows()==size);
  }
  else
  {    
    m_matrix.resize(size,size);
    m_matrix.setZero();
    m_transpositions.resize(size);
    for (Index i = 0; i < size; i++)
      m_transpositions.coeffRef(i) = IndexType(i);
    m_temporary.resize(size);
    m_sign = sigma>=0 ? internal::PositiveSemiDef : internal::NegativeSemiDef;
    m_isInitialized = true;
  }

  internal::ldlt_inplace<UpLo>::update(m_matrix, m_transpositions, m_temporary, w, sigma);

  return *this;
}
template<typename MatrixType , int _UpLo>
MatrixType Eigen::LDLT< MatrixType, _UpLo >::reconstructedMatrix ( ) const
Returns:
the matrix represented by the decomposition, i.e., it returns the product: P^T L D L^* P. This function is provided for debug purpose.

Definition at line 554 of file LDLT.h.

{
  eigen_assert(m_isInitialized && "LDLT is not initialized.");
  const Index size = m_matrix.rows();
  MatrixType res(size,size);

  // P
  res.setIdentity();
  res = transpositionsP() * res;
  // L^* P
  res = matrixU() * res;
  // D(L^*P)
  res = vectorD().real().asDiagonal() * res;
  // L(DL^*P)
  res = matrixL() * res;
  // P^T (LDL^*P)
  res = transpositionsP().transpose() * res;

  return res;
}
template<typename _MatrixType, int _UpLo>
Index Eigen::LDLT< _MatrixType, _UpLo >::rows ( ) const [inline]

Definition at line 210 of file LDLT.h.

{ return m_matrix.rows(); }
template<typename _MatrixType, int _UpLo>
void Eigen::LDLT< _MatrixType, _UpLo >::setZero ( ) [inline]

Clear any existing decomposition

See also:
rankUpdate(w,sigma)

Definition at line 116 of file LDLT.h.

    {
      m_isInitialized = false;
    }
template<typename _MatrixType, int _UpLo>
template<typename Rhs >
const Solve<LDLT, Rhs> Eigen::LDLT< _MatrixType, _UpLo >::solve ( const MatrixBase< Rhs > &  b) const [inline]
Returns:
a solution x of $ A x = b $ using the current decomposition of A.

This function also supports in-place solves using the syntax x = decompositionObject.solve(x) .

More precisely, this method solves $ A x = b $ using the decomposition $ A = P^T L D L^* P $ by solving the systems $ P^T y_1 = b $, $ L y_2 = y_1 $, $ D y_3 = y_2 $, $ L^* y_4 = y_3 $ and $ P x = y_4 $ in succession. If the matrix $ A $ is singular, then $ D $ will also be singular (all the other matrices are invertible). In that case, the least-square solution of $ D y_3 = y_2 $ is computed. This does not mean that this function computes the least-square solution of $ A x = b $ is $ A $ is singular.

See also:
MatrixBase::ldlt(), SelfAdjointView::ldlt()

Definition at line 181 of file LDLT.h.

    {
      eigen_assert(m_isInitialized && "LDLT is not initialized.");
      eigen_assert(m_matrix.rows()==b.rows()
                && "LDLT::solve(): invalid number of rows of the right hand side matrix b");
      return Solve<LDLT, Rhs>(*this, b.derived());
    }
template<typename MatrixType , int _UpLo>
template<typename Derived >
bool Eigen::LDLT< MatrixType, _UpLo >::solveInPlace ( MatrixBase< Derived > &  bAndX) const

use x = ldlt_object.solve(x);

This is the in-place version of solve().

Parameters:
bAndXrepresents both the right-hand side matrix b and result x.
Returns:
true always! If you need to check for existence of solutions, use another decomposition like LU, QR, or SVD.

This version avoids a copy when the right hand side matrix b is not needed anymore.

See also:
LDLT::solve(), MatrixBase::ldlt()

Definition at line 540 of file LDLT.h.

{
  eigen_assert(m_isInitialized && "LDLT is not initialized.");
  eigen_assert(m_matrix.rows() == bAndX.rows());

  bAndX = this->solve(bAndX);

  return true;
}
template<typename _MatrixType, int _UpLo>
const TranspositionType& Eigen::LDLT< _MatrixType, _UpLo >::transpositionsP ( ) const [inline]
Returns:
the permutation matrix P as a transposition sequence.

Definition at line 137 of file LDLT.h.

    {
      eigen_assert(m_isInitialized && "LDLT is not initialized.");
      return m_transpositions;
    }
template<typename _MatrixType, int _UpLo>
Diagonal<const MatrixType> Eigen::LDLT< _MatrixType, _UpLo >::vectorD ( ) const [inline]
Returns:
the coefficients of the diagonal matrix D

Definition at line 144 of file LDLT.h.

    {
      eigen_assert(m_isInitialized && "LDLT is not initialized.");
      return m_matrix.diagonal();
    }

Member Data Documentation

template<typename _MatrixType, int _UpLo>
bool Eigen::LDLT< _MatrixType, _UpLo >::m_isInitialized [protected]

Definition at line 247 of file LDLT.h.

template<typename _MatrixType, int _UpLo>
MatrixType Eigen::LDLT< _MatrixType, _UpLo >::m_matrix [protected]

Used to compute and store the Cholesky decomposition A = L D L^* = U^* D U. The strict upper part is used during the decomposition, the strict lower part correspond to the coefficients of L (its diagonal is equal to 1 and is not stored), and the diagonal entries correspond to D.

Definition at line 243 of file LDLT.h.

template<typename _MatrixType, int _UpLo>
internal::SignMatrix Eigen::LDLT< _MatrixType, _UpLo >::m_sign [protected]

Definition at line 246 of file LDLT.h.

template<typename _MatrixType, int _UpLo>
TmpMatrixType Eigen::LDLT< _MatrixType, _UpLo >::m_temporary [protected]

Definition at line 245 of file LDLT.h.

template<typename _MatrixType, int _UpLo>
TranspositionType Eigen::LDLT< _MatrixType, _UpLo >::m_transpositions [protected]

Definition at line 244 of file LDLT.h.


The documentation for this class was generated from the following file:
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines