Eigen::SimplicialCholeskyBase< Derived > Class Template Reference

A base class for direct sparse Cholesky factorizations. More...

Inheritance diagram for Eigen::SimplicialCholeskyBase< Derived >:

Classes
struct	keep_diag

Public Types
enum	{ UpLo }
enum	{   ColsAtCompileTime ,   MaxColsAtCompileTime }
typedef SparseMatrix< Scalar, ColMajor, StorageIndex >	CholMatrixType
typedef CholMatrixType const *	ConstCholMatrixPtr
typedef internal::traits< Derived >::MatrixType	MatrixType
typedef internal::traits< Derived >::OrderingType	OrderingType
typedef MatrixType::RealScalar	RealScalar
typedef MatrixType::Scalar	Scalar
typedef MatrixType::StorageIndex	StorageIndex
typedef Matrix< StorageIndex, Dynamic, 1 >	VectorI
typedef Matrix< Scalar, Dynamic, 1 >	VectorType

Public Member Functions
Index	cols () const
Derived &	derived ()
Derived &	derived ()
const Derived &	derived () const
const Derived &	derived () const
ComputationInfo	info () const
	Reports whether previous computation was successful. More...
const PermutationMatrix< Dynamic, Dynamic, StorageIndex > &	permutationP () const
const PermutationMatrix< Dynamic, Dynamic, StorageIndex > &	permutationPinv () const
Index	rows () const
Derived &	setShift (const RealScalar &offset, const RealScalar &scale=1)
	SimplicialCholeskyBase ()
	SimplicialCholeskyBase (const MatrixType &matrix)
	~SimplicialCholeskyBase ()
Public Member Functions inherited from Eigen::SparseSolverBase< Derived >
Derived &	derived ()
const Derived &	derived () const
template<typename Rhs >
const Solve< Derived, Rhs >	solve (const MatrixBase< Rhs > &b) const
template<typename Rhs >
const Solve< Derived, Rhs >	solve (const SparseMatrixBase< Rhs > &b) const
	SparseSolverBase ()
	SparseSolverBase (SparseSolverBase &&other)
	~SparseSolverBase ()

Protected Member Functions
void	analyzePattern (const MatrixType &a, bool doLDLT)
void	analyzePattern_preordered (const CholMatrixType &a, bool doLDLT)
template<bool DoLDLT>
void	compute (const MatrixType &matrix)
template<bool DoLDLT>
void	factorize (const MatrixType &a)
template<bool DoLDLT>
void	factorize_preordered (const CholMatrixType &a)
void	ordering (const MatrixType &a, ConstCholMatrixPtr &pmat, CholMatrixType &ap)

Protected Attributes
bool	m_analysisIsOk
VectorType	m_diag
bool	m_factorizationIsOk
ComputationInfo	m_info
CholMatrixType	m_matrix
VectorI	m_nonZerosPerCol
PermutationMatrix< Dynamic, Dynamic, StorageIndex >	m_P
VectorI	m_parent
PermutationMatrix< Dynamic, Dynamic, StorageIndex >	m_Pinv
RealScalar	m_shiftOffset
RealScalar	m_shiftScale
Protected Attributes inherited from Eigen::SparseSolverBase< Derived >
bool	m_isInitialized

Private Types
typedef SparseSolverBase< Derived >	Base

Private Attributes
bool	m_isInitialized

Detailed Description

template<typename Derived>
class Eigen::SimplicialCholeskyBase< Derived >

A base class for direct sparse Cholesky factorizations.

This is a base class for LL^T and LDL^T Cholesky factorizations of sparse matrices that are selfadjoint and positive definite. These factorizations allow for solving A.X = B where X and B can be either dense or sparse.

In order to reduce the fill-in, a symmetric permutation P is applied prior to the factorization such that the factorized matrix is P A P^-1.

Template Parameters

Derived

the type of the derived class, that is the actual factorization type.

Definition at line 57 of file SimplicialCholesky.h.

Member Typedef Documentation

Base

template<typename Derived >

typedef SparseSolverBase<Derived> Eigen::SimplicialCholeskyBase< Derived >::Base

private

Definition at line 59 of file SimplicialCholesky.h.

CholMatrixType

template<typename Derived >

typedef SparseMatrix<Scalar,ColMajor,StorageIndex> Eigen::SimplicialCholeskyBase< Derived >::CholMatrixType

Definition at line 69 of file SimplicialCholesky.h.

ConstCholMatrixPtr

template<typename Derived >

typedef CholMatrixType const* Eigen::SimplicialCholeskyBase< Derived >::ConstCholMatrixPtr

Definition at line 70 of file SimplicialCholesky.h.

MatrixType

template<typename Derived >

typedef internal::traits<Derived>::MatrixType Eigen::SimplicialCholeskyBase< Derived >::MatrixType

Definition at line 63 of file SimplicialCholesky.h.

OrderingType

template<typename Derived >

typedef internal::traits<Derived>::OrderingType Eigen::SimplicialCholeskyBase< Derived >::OrderingType

Definition at line 64 of file SimplicialCholesky.h.

RealScalar

template<typename Derived >

typedef MatrixType::RealScalar Eigen::SimplicialCholeskyBase< Derived >::RealScalar

Definition at line 67 of file SimplicialCholesky.h.

Scalar

template<typename Derived >

typedef MatrixType::Scalar Eigen::SimplicialCholeskyBase< Derived >::Scalar

Definition at line 66 of file SimplicialCholesky.h.

StorageIndex

template<typename Derived >

typedef MatrixType::StorageIndex Eigen::SimplicialCholeskyBase< Derived >::StorageIndex

Definition at line 68 of file SimplicialCholesky.h.

VectorI

template<typename Derived >

typedef Matrix<StorageIndex,Dynamic,1> Eigen::SimplicialCholeskyBase< Derived >::VectorI

Definition at line 72 of file SimplicialCholesky.h.

VectorType

template<typename Derived >

typedef Matrix<Scalar,Dynamic,1> Eigen::SimplicialCholeskyBase< Derived >::VectorType

Definition at line 71 of file SimplicialCholesky.h.

Member Enumeration Documentation

anonymous enum

template<typename Derived >

anonymous enum

Enumerator
UpLo

Definition at line 65 of file SimplicialCholesky.h.

{ UpLo = internal::traits<Derived>::UpLo };

anonymous enum

template<typename Derived >

anonymous enum

Enumerator
ColsAtCompileTime
MaxColsAtCompileTime

Definition at line 74 of file SimplicialCholesky.h.

         {
      ColsAtCompileTime = MatrixType::ColsAtCompileTime,
      MaxColsAtCompileTime = MatrixType::MaxColsAtCompileTime
    };

Constructor & Destructor Documentation

SimplicialCholeskyBase() [1/2]

template<typename Derived >

Eigen::SimplicialCholeskyBase< Derived >::SimplicialCholeskyBase ( )

inline

Default constructor

Definition at line 84 of file SimplicialCholesky.h.

      : m_info(Success),
        m_factorizationIsOk(false),
        m_analysisIsOk(false),
        m_shiftOffset(0),
        m_shiftScale(1)
    {}

SimplicialCholeskyBase() [2/2]

template<typename Derived >

Eigen::SimplicialCholeskyBase< Derived >::SimplicialCholeskyBase ( const MatrixType &  matrix )

inlineexplicit

Definition at line 92 of file SimplicialCholesky.h.

      : m_info(Success),
        m_factorizationIsOk(false),
        m_analysisIsOk(false),
        m_shiftOffset(0),
        m_shiftScale(1)
    {
      derived().compute(matrix);
    }

~SimplicialCholeskyBase()

template<typename Derived >

Eigen::SimplicialCholeskyBase< Derived >::~SimplicialCholeskyBase ( )

inline

Definition at line 102 of file SimplicialCholesky.h.

    {
    }

Member Function Documentation

analyzePattern()

template<typename Derived >

void Eigen::SimplicialCholeskyBase< Derived >::analyzePattern ( const MatrixType &  a, bool  doLDLT  )

inlineprotected

Definition at line 240 of file SimplicialCholesky.h.

    {
      eigen_assert(a.rows()==a.cols());
      Index size = a.cols();
      CholMatrixType tmp(size,size);
      ConstCholMatrixPtr pmat;
      ordering(a, pmat, tmp);
      analyzePattern_preordered(*pmat,doLDLT);
    }

analyzePattern_preordered()

template<typename Derived >

void Eigen::SimplicialCholeskyBase< Derived >::analyzePattern_preordered ( const CholMatrixType &  a, bool  doLDLT  )

protected

Definition at line 28 of file SimplicialCholesky_impl.h.

{
  const StorageIndex size = StorageIndex(ap.rows());
  m_matrix.resize(size, size);
  m_parent.resize(size);
  m_nonZerosPerCol.resize(size);
 
  ei_declare_aligned_stack_constructed_variable(StorageIndex, tags, size, 0);
 
  for(StorageIndex k = 0; k < size; ++k)
  {
    /* L(k,:) pattern: all nodes reachable in etree from nz in A(0:k-1,k) */
    m_parent[k] = -1;             /* parent of k is not yet known */
    tags[k] = k;                  /* mark node k as visited */
    m_nonZerosPerCol[k] = 0;      /* count of nonzeros in column k of L */
    for(typename CholMatrixType::InnerIterator it(ap,k); it; ++it)
    {
      StorageIndex i = it.index();
      if(i < k)
      {
        /* follow path from i to root of etree, stop at flagged node */
        for(; tags[i] != k; i = m_parent[i])
        {
          /* find parent of i if not yet determined */
          if (m_parent[i] == -1)
            m_parent[i] = k;
          m_nonZerosPerCol[i]++;        /* L (k,i) is nonzero */
          tags[i] = k;                  /* mark i as visited */
        }
      }
    }
  }
 
  /* construct Lp index array from m_nonZerosPerCol column counts */
  StorageIndex* Lp = m_matrix.outerIndexPtr();
  Lp[0] = 0;
  for(StorageIndex k = 0; k < size; ++k)
    Lp[k+1] = Lp[k] + m_nonZerosPerCol[k] + (doLDLT ? 0 : 1);
 
  m_matrix.resizeNonZeros(Lp[size]);
 
  m_isInitialized     = true;
  m_info              = Success;
  m_analysisIsOk      = true;
  m_factorizationIsOk = false;
}

cols()

template<typename Derived >

Index Eigen::SimplicialCholeskyBase< Derived >::cols ( void  ) const

inline

Definition at line 109 of file SimplicialCholesky.h.

{ return m_matrix.cols(); }

compute()

template<typename Derived >

template<bool DoLDLT>

void Eigen::SimplicialCholeskyBase< Derived >::compute ( const MatrixType &  matrix )

inlineprotected

Computes the sparse Cholesky decomposition of matrix

Definition at line 204 of file SimplicialCholesky.h.

    {
      eigen_assert(matrix.rows()==matrix.cols());
      Index size = matrix.cols();
      CholMatrixType tmp(size,size);
      ConstCholMatrixPtr pmat;
      ordering(matrix, pmat, tmp);
      analyzePattern_preordered(*pmat, DoLDLT);
      factorize_preordered<DoLDLT>(*pmat);
    }

derived() [1/4]

template<typename Derived >

Derived& Eigen::SparseSolverBase< Derived >::derived

inline

Definition at line 83 of file SparseSolverBase.h.

{ return *static_cast<Derived*>(this); }

derived() [2/4]

template<typename Derived >

Derived& Eigen::SimplicialCholeskyBase< Derived >::derived ( )

inline

Definition at line 106 of file SimplicialCholesky.h.

{ return *static_cast<Derived*>(this); }

derived() [3/4]

template<typename Derived >

const Derived& Eigen::SparseSolverBase< Derived >::derived

inline

Definition at line 84 of file SparseSolverBase.h.

{ return *static_cast<const Derived*>(this); }

derived() [4/4]

template<typename Derived >

const Derived& Eigen::SimplicialCholeskyBase< Derived >::derived ( ) const

inline

Definition at line 107 of file SimplicialCholesky.h.

{ return *static_cast<const Derived*>(this); }

factorize()

template<typename Derived >

template<bool DoLDLT>

void Eigen::SimplicialCholeskyBase< Derived >::factorize ( const MatrixType &  a )

inlineprotected

Definition at line 216 of file SimplicialCholesky.h.

    {
      eigen_assert(a.rows()==a.cols());
      Index size = a.cols();
      CholMatrixType tmp(size,size);
      ConstCholMatrixPtr pmat;
      
      if(m_P.size() == 0 && (int(UpLo) & int(Upper)) == Upper)
      {
        // If there is no ordering, try to directly use the input matrix without any copy
        internal::simplicial_cholesky_grab_input<CholMatrixType,MatrixType>::run(a, pmat, tmp);
      }
      else
      {
        tmp.template selfadjointView<Upper>() = a.template selfadjointView<UpLo>().twistedBy(m_P);
        pmat = &tmp;
      }
      
      factorize_preordered<DoLDLT>(*pmat);
    }

factorize_preordered()

template<typename Derived >

template<bool DoLDLT>

void Eigen::SimplicialCholeskyBase< Derived >::factorize_preordered ( const CholMatrixType &  a )

protected

Definition at line 78 of file SimplicialCholesky_impl.h.

{
  using std::sqrt;
 
  eigen_assert(m_analysisIsOk && "You must first call analyzePattern()");
  eigen_assert(ap.rows()==ap.cols());
  eigen_assert(m_parent.size()==ap.rows());
  eigen_assert(m_nonZerosPerCol.size()==ap.rows());
 
  const StorageIndex size = StorageIndex(ap.rows());
  const StorageIndex* Lp = m_matrix.outerIndexPtr();
  StorageIndex* Li = m_matrix.innerIndexPtr();
  Scalar* Lx = m_matrix.valuePtr();
 
  ei_declare_aligned_stack_constructed_variable(Scalar, y, size, 0);
  ei_declare_aligned_stack_constructed_variable(StorageIndex,  pattern, size, 0);
  ei_declare_aligned_stack_constructed_variable(StorageIndex,  tags, size, 0);
 
  bool ok = true;
  m_diag.resize(DoLDLT ? size : 0);
 
  for(StorageIndex k = 0; k < size; ++k)
  {
    // compute nonzero pattern of kth row of L, in topological order
    y[k] = Scalar(0);                     // Y(0:k) is now all zero
    StorageIndex top = size;               // stack for pattern is empty
    tags[k] = k;                    // mark node k as visited
    m_nonZerosPerCol[k] = 0;        // count of nonzeros in column k of L
    for(typename CholMatrixType::InnerIterator it(ap,k); it; ++it)
    {
      StorageIndex i = it.index();
      if(i <= k)
      {
        y[i] += numext::conj(it.value());            /* scatter A(i,k) into Y (sum duplicates) */
        Index len;
        for(len = 0; tags[i] != k; i = m_parent[i])
        {
          pattern[len++] = i;     /* L(k,i) is nonzero */
          tags[i] = k;            /* mark i as visited */
        }
        while(len > 0)
          pattern[--top] = pattern[--len];
      }
    }
 
    /* compute numerical values kth row of L (a sparse triangular solve) */
 
    RealScalar d = numext::real(y[k]) * m_shiftScale + m_shiftOffset;    // get D(k,k), apply the shift function, and clear Y(k)
    y[k] = Scalar(0);
    for(; top < size; ++top)
    {
      Index i = pattern[top];       /* pattern[top:n-1] is pattern of L(:,k) */
      Scalar yi = y[i];             /* get and clear Y(i) */
      y[i] = Scalar(0);
 
      /* the nonzero entry L(k,i) */
      Scalar l_ki;
      if(DoLDLT)
        l_ki = yi / numext::real(m_diag[i]);
      else
        yi = l_ki = yi / Lx[Lp[i]];
 
      Index p2 = Lp[i] + m_nonZerosPerCol[i];
      Index p;
      for(p = Lp[i] + (DoLDLT ? 0 : 1); p < p2; ++p)
        y[Li[p]] -= numext::conj(Lx[p]) * yi;
      d -= numext::real(l_ki * numext::conj(yi));
      Li[p] = k;                          /* store L(k,i) in column form of L */
      Lx[p] = l_ki;
      ++m_nonZerosPerCol[i];              /* increment count of nonzeros in col i */
    }
    if(DoLDLT)
    {
      m_diag[k] = d;
      if(d == RealScalar(0))
      {
        ok = false;                         /* failure, D(k,k) is zero */
        break;
      }
    }
    else
    {
      Index p = Lp[k] + m_nonZerosPerCol[k]++;
      Li[p] = k ;                /* store L(k,k) = sqrt (d) in column k */
      if(d <= RealScalar(0)) {
        ok = false;              /* failure, matrix is not positive definite */
        break;
      }
      Lx[p] = sqrt(d) ;
    }
  }
 
  m_info = ok ? Success : NumericalIssue;
  m_factorizationIsOk = true;
}

info()

template<typename Derived >

ComputationInfo Eigen::SimplicialCholeskyBase< Derived >::info ( ) const

inline

Reports whether previous computation was successful.

Returns

Success if computation was successful, NumericalIssue if the matrix.appears to be negative.

Definition at line 117 of file SimplicialCholesky.h.

    {
      eigen_assert(m_isInitialized && "Decomposition is not initialized.");
      return m_info;
    }

ordering()

template<typename Derived >

void Eigen::SimplicialCholeskyBase< Derived >::ordering ( const MatrixType &  a, ConstCholMatrixPtr &  pmat, CholMatrixType &  ap  )

protected

Definition at line 660 of file SimplicialCholesky.h.

{
  eigen_assert(a.rows()==a.cols());
  const Index size = a.rows();
  pmat = ≈
  // Note that ordering methods compute the inverse permutation
  if(!internal::is_same<OrderingType,NaturalOrdering<Index> >::value)
  {
    {
      CholMatrixType C;
      C = a.template selfadjointView<UpLo>();
      
      OrderingType ordering;
      ordering(C,m_Pinv);
    }
 
    if(m_Pinv.size()>0) m_P = m_Pinv.inverse();
    else                m_P.resize(0);
    
    ap.resize(size,size);
    ap.template selfadjointView<Upper>() = a.template selfadjointView<UpLo>().twistedBy(m_P);
  }
  else
  {
    m_Pinv.resize(0);
    m_P.resize(0);
    if(int(UpLo)==int(Lower) || MatrixType::IsRowMajor)
    {
      // we have to transpose the lower part to to the upper one
      ap.resize(size,size);
      ap.template selfadjointView<Upper>() = a.template selfadjointView<UpLo>();
    }
    else
      internal::simplicial_cholesky_grab_input<CholMatrixType,MatrixType>::run(a, pmat, ap);
  }  
}

permutationP()

template<typename Derived >

const PermutationMatrix<Dynamic,Dynamic,StorageIndex>& Eigen::SimplicialCholeskyBase< Derived >::permutationP ( ) const

inline

Returns

the permutation P

See also

permutationPinv()

Definition at line 125 of file SimplicialCholesky.h.

    { return m_P; }

permutationPinv()

template<typename Derived >

const PermutationMatrix<Dynamic,Dynamic,StorageIndex>& Eigen::SimplicialCholeskyBase< Derived >::permutationPinv ( ) const

inline

Returns

the inverse P^-1 of the permutation P

See also

permutationP()

Definition at line 130 of file SimplicialCholesky.h.

    { return m_Pinv; }

rows()

template<typename Derived >

Index Eigen::SimplicialCholeskyBase< Derived >::rows ( void  ) const

inline

Definition at line 110 of file SimplicialCholesky.h.

{ return m_matrix.rows(); }

setShift()

template<typename Derived >

Derived& Eigen::SimplicialCholeskyBase< Derived >::setShift ( const RealScalar &  offset, const RealScalar &  scale = 1  )

inline

Sets the shift parameters that will be used to adjust the diagonal coefficients during the numerical factorization.

During the numerical factorization, the diagonal coefficients are transformed by the following linear model:
d_ii = offset + scale * d_ii

The default is the identity transformation with offset=0, and scale=1.

Returns

a reference to *this.

Definition at line 142 of file SimplicialCholesky.h.

    {
      m_shiftOffset = offset;
      m_shiftScale = scale;
      return derived();
    }

Member Data Documentation

m_analysisIsOk

template<typename Derived >

bool Eigen::SimplicialCholeskyBase< Derived >::m_analysisIsOk

protected

Definition at line 263 of file SimplicialCholesky.h.

m_diag

template<typename Derived >

VectorType Eigen::SimplicialCholeskyBase< Derived >::m_diag

protected

Definition at line 266 of file SimplicialCholesky.h.

m_factorizationIsOk

template<typename Derived >

bool Eigen::SimplicialCholeskyBase< Derived >::m_factorizationIsOk

protected

Definition at line 262 of file SimplicialCholesky.h.

m_info

template<typename Derived >

ComputationInfo Eigen::SimplicialCholeskyBase< Derived >::m_info

mutableprotected

Definition at line 261 of file SimplicialCholesky.h.

m_isInitialized

template<typename Derived >

bool Eigen::SparseSolverBase< Derived >::m_isInitialized

mutableprivate

Definition at line 123 of file SparseSolverBase.h.

m_matrix

template<typename Derived >

CholMatrixType Eigen::SimplicialCholeskyBase< Derived >::m_matrix

protected

Definition at line 265 of file SimplicialCholesky.h.

m_nonZerosPerCol

template<typename Derived >

VectorI Eigen::SimplicialCholeskyBase< Derived >::m_nonZerosPerCol

protected

Definition at line 268 of file SimplicialCholesky.h.

m_P

template<typename Derived >

PermutationMatrix<Dynamic,Dynamic,StorageIndex> Eigen::SimplicialCholeskyBase< Derived >::m_P

protected

Definition at line 269 of file SimplicialCholesky.h.

m_parent

template<typename Derived >

VectorI Eigen::SimplicialCholeskyBase< Derived >::m_parent

protected

Definition at line 267 of file SimplicialCholesky.h.

m_Pinv

template<typename Derived >

PermutationMatrix<Dynamic,Dynamic,StorageIndex> Eigen::SimplicialCholeskyBase< Derived >::m_Pinv

protected

Definition at line 270 of file SimplicialCholesky.h.

m_shiftOffset

template<typename Derived >

RealScalar Eigen::SimplicialCholeskyBase< Derived >::m_shiftOffset

protected

Definition at line 272 of file SimplicialCholesky.h.

m_shiftScale

template<typename Derived >

RealScalar Eigen::SimplicialCholeskyBase< Derived >::m_shiftScale

protected

Definition at line 273 of file SimplicialCholesky.h.

---

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

• SimplicialCholesky.h
• SimplicialCholesky_impl.h