Prusa Slicer 2.6.0
Loading...
Searching...
No Matches
Eigen::CholmodDecomposition< _MatrixType, _UpLo > Class Template Reference

A general Cholesky factorization and solver based on Cholmod. More...

#include <src/eigen/Eigen/src/CholmodSupport/CholmodSupport.h>

+ Inheritance diagram for Eigen::CholmodDecomposition< _MatrixType, _UpLo >:
+ Collaboration diagram for Eigen::CholmodDecomposition< _MatrixType, _UpLo >:

Public Types

typedef _MatrixType MatrixType
 
enum  { UpLo = _UpLo }
 
enum  { ColsAtCompileTime = MatrixType::ColsAtCompileTime , MaxColsAtCompileTime = MatrixType::MaxColsAtCompileTime }
 
typedef MatrixType::Scalar Scalar
 
typedef MatrixType::RealScalar RealScalar
 
typedef MatrixType CholMatrixType
 
typedef MatrixType::StorageIndex StorageIndex
 

Public Member Functions

 CholmodDecomposition ()
 
 CholmodDecomposition (const MatrixType &matrix)
 
 ~CholmodDecomposition ()
 
void setMode (CholmodMode mode)
 
StorageIndex cols () const
 
StorageIndex rows () const
 
ComputationInfo info () const
 Reports whether previous computation was successful.
 
Derived & compute (const MatrixType &matrix)
 
void analyzePattern (const MatrixType &matrix)
 
void factorize (const MatrixType &matrix)
 
cholmod_common & cholmod ()
 
template<typename Rhs , typename Dest >
void _solve_impl (const MatrixBase< Rhs > &b, MatrixBase< Dest > &dest) const
 
template<typename RhsDerived , typename DestDerived >
void _solve_impl (const SparseMatrixBase< RhsDerived > &b, SparseMatrixBase< DestDerived > &dest) const
 
template<typename Rhs , typename Dest >
void _solve_impl (const SparseMatrixBase< Rhs > &b, SparseMatrixBase< Dest > &dest) const
 
Derived & setShift (const RealScalar &offset)
 
Scalar determinant () const
 
Scalar logDeterminant () const
 
template<typename Stream >
void dumpMemory (Stream &)
 
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
 

Protected Member Functions

void init ()
 
Derived & derived ()
 
const Derived & derived () const
 

Protected Attributes

cholmod_factor * m_cholmodFactor
 
double m_shiftOffset [2]
 
ComputationInfo m_info
 
int m_factorizationIsOk
 
int m_analysisIsOk
 
bool m_isInitialized
 

Private Types

typedef CholmodBase< _MatrixType, _UpLo, CholmodDecompositionBase
 

Private Attributes

cholmod_common m_cholmod
 

Detailed Description

template<typename _MatrixType, int _UpLo = Lower>
class Eigen::CholmodDecomposition< _MatrixType, _UpLo >

A general Cholesky factorization and solver based on Cholmod.

This class allows to solve for A.X = B sparse linear problems via a LL^T or LDL^T Cholesky factorization using the Cholmod library. The sparse matrix A must be selfadjoint and positive definite. The vectors or matrices X and B can be either dense or sparse.

This variant permits to change the underlying Cholesky method at runtime. On the other hand, it does not provide access to the result of the factorization. The default is to let Cholmod automatically choose between a simplicial and supernodal factorization.

Template Parameters
_MatrixTypethe type of the sparse matrix A, it must be a SparseMatrix<>
_UpLothe triangular part that will be used for the computations. It can be Lower or Upper. Default is Lower.

\implsparsesolverconcept

This class supports all kind of SparseMatrix<>: row or column major; upper, lower, or both; compressed or non compressed.

Warning
Only double precision real and complex scalar types are supported by Cholmod.
See also
TutorialSparseSolverConcept

Member Typedef Documentation

◆ Base

template<typename _MatrixType , int _UpLo = Lower>
typedef CholmodBase<_MatrixType, _UpLo, CholmodDecomposition> Eigen::CholmodDecomposition< _MatrixType, _UpLo >::Base
private

◆ CholMatrixType

template<typename _MatrixType , int _UpLo, typename Derived >
typedef MatrixType Eigen::CholmodBase< _MatrixType, _UpLo, Derived >::CholMatrixType
inherited

◆ MatrixType

template<typename _MatrixType , int _UpLo = Lower>
typedef _MatrixType Eigen::CholmodDecomposition< _MatrixType, _UpLo >::MatrixType

◆ RealScalar

template<typename _MatrixType , int _UpLo, typename Derived >
typedef MatrixType::RealScalar Eigen::CholmodBase< _MatrixType, _UpLo, Derived >::RealScalar
inherited

◆ Scalar

template<typename _MatrixType , int _UpLo, typename Derived >
typedef MatrixType::Scalar Eigen::CholmodBase< _MatrixType, _UpLo, Derived >::Scalar
inherited

◆ StorageIndex

template<typename _MatrixType , int _UpLo, typename Derived >
typedef MatrixType::StorageIndex Eigen::CholmodBase< _MatrixType, _UpLo, Derived >::StorageIndex
inherited

Member Enumeration Documentation

◆ anonymous enum

template<typename _MatrixType , int _UpLo, typename Derived >
anonymous enum
inherited
Enumerator
UpLo 
181{ UpLo = _UpLo };
@ UpLo
Definition CholmodSupport.h:181

◆ anonymous enum

template<typename _MatrixType , int _UpLo, typename Derived >
anonymous enum
inherited
Enumerator
ColsAtCompileTime 
MaxColsAtCompileTime 
186 {
187 ColsAtCompileTime = MatrixType::ColsAtCompileTime,
188 MaxColsAtCompileTime = MatrixType::MaxColsAtCompileTime
189 };
@ MaxColsAtCompileTime
Definition CholmodSupport.h:188
@ ColsAtCompileTime
Definition CholmodSupport.h:187

Constructor & Destructor Documentation

◆ CholmodDecomposition() [1/2]

template<typename _MatrixType , int _UpLo = Lower>
Eigen::CholmodDecomposition< _MatrixType, _UpLo >::CholmodDecomposition ( )
inline
594: Base() { init(); }
void init()
Definition CholmodSupport.h:630
CholmodBase< _MatrixType, _UpLo, CholmodDecomposition > Base
Definition CholmodSupport.h:587

References Eigen::CholmodDecomposition< _MatrixType, _UpLo >::init().

+ Here is the call graph for this function:

◆ CholmodDecomposition() [2/2]

template<typename _MatrixType , int _UpLo = Lower>
Eigen::CholmodDecomposition< _MatrixType, _UpLo >::CholmodDecomposition ( const MatrixType matrix)
inline
596 : Base()
597 {
598 init();
599 this->compute(matrix);
600 }
Derived & compute(const MatrixType &matrix)
Definition CholmodSupport.h:232

References Eigen::CholmodBase< _MatrixType, _UpLo, Derived >::compute(), and Eigen::CholmodDecomposition< _MatrixType, _UpLo >::init().

+ Here is the call graph for this function:

◆ ~CholmodDecomposition()

template<typename _MatrixType , int _UpLo = Lower>
Eigen::CholmodDecomposition< _MatrixType, _UpLo >::~CholmodDecomposition ( )
inline
602{}

Member Function Documentation

◆ _solve_impl() [1/3]

template<typename _MatrixType , int _UpLo, typename Derived >
template<typename Rhs , typename Dest >
void Eigen::CholmodBase< _MatrixType, _UpLo, Derived >::_solve_impl ( const MatrixBase< Rhs > &  b,
MatrixBase< Dest > &  dest 
) const
inlineinherited
286 {
287 eigen_assert(m_factorizationIsOk && "The decomposition is not in a valid state for solving, you must first call either compute() or symbolic()/numeric()");
288 const Index size = m_cholmodFactor->n;
290 eigen_assert(size==b.rows());
291
292 // Cholmod needs column-major stoarge without inner-stride, which corresponds to the default behavior of Ref.
293 Ref<const Matrix<typename Rhs::Scalar,Dynamic,Dynamic,ColMajor> > b_ref(b.derived());
294
295 cholmod_dense b_cd = viewAsCholmod(b_ref);
296 cholmod_dense* x_cd = cholmod_solve(CHOLMOD_A, m_cholmodFactor, &b_cd, &m_cholmod);
297 if(!x_cd)
298 {
299 this->m_info = NumericalIssue;
300 return;
301 }
302 // TODO optimize this copy by swapping when possible (be careful with alignment, etc.)
303 dest = Matrix<Scalar,Dest::RowsAtCompileTime,Dest::ColsAtCompileTime>::Map(reinterpret_cast<Scalar*>(x_cd->x),b.rows(),b.cols());
304 cholmod_free_dense(&x_cd, &m_cholmod);
305 }
#define EIGEN_UNUSED_VARIABLE(var)
Definition Macros.h:618
#define eigen_assert(x)
Definition Macros.h:579
cholmod_factor * m_cholmodFactor
Definition CholmodSupport.h:404
cholmod_common m_cholmod
Definition CholmodSupport.h:403
ComputationInfo m_info
Definition CholmodSupport.h:406
int m_factorizationIsOk
Definition CholmodSupport.h:407
MatrixType::Scalar Scalar
Definition CholmodSupport.h:182
@ NumericalIssue
Definition Constants.h:434
EIGEN_DEFAULT_DENSE_INDEX_TYPE Index
The Index type as used for the API.
Definition Meta.h:33
cholmod_sparse viewAsCholmod(Ref< SparseMatrix< _Scalar, _Options, _StorageIndex > > mat)
Definition CholmodSupport.h:58
constexpr auto size(const C &c) -> decltype(c.size())
Definition span.hpp:183

References Eigen::PlainObjectBase< Matrix< _Scalar, _Rows, _Cols, _Options, _MaxRows, _MaxCols > >::Eigen::Map, eigen_assert, EIGEN_UNUSED_VARIABLE, Eigen::CholmodBase< _MatrixType, _UpLo, Derived >::m_cholmod, Eigen::CholmodBase< _MatrixType, _UpLo, Derived >::m_cholmodFactor, Eigen::CholmodBase< _MatrixType, _UpLo, Derived >::m_factorizationIsOk, Eigen::CholmodBase< _MatrixType, _UpLo, Derived >::m_info, Eigen::NumericalIssue, and Eigen::viewAsCholmod().

+ Here is the call graph for this function:

◆ _solve_impl() [2/3]

template<typename Derived >
template<typename Rhs , typename Dest >
void Eigen::SparseSolverBase< Derived >::_solve_impl ( const SparseMatrixBase< Rhs > &  b,
SparseMatrixBase< Dest > &  dest 
) const
inlineinherited
112 {
113 internal::solve_sparse_through_dense_panels(derived(), b.derived(), dest.derived());
114 }
Derived & derived()
Definition SparseSolverBase.h:79
enable_if< Rhs::ColsAtCompileTime!=1 &&Dest::ColsAtCompileTime!=1 >::type solve_sparse_through_dense_panels(const Decomposition &dec, const Rhs &rhs, Dest &dest)
Definition SparseSolverBase.h:23

References Eigen::SparseSolverBase< Derived >::derived(), Eigen::SparseMatrixBase< Derived >::derived(), and Eigen::internal::solve_sparse_through_dense_panels().

+ Here is the call graph for this function:

◆ _solve_impl() [3/3]

template<typename _MatrixType , int _UpLo, typename Derived >
template<typename RhsDerived , typename DestDerived >
void Eigen::CholmodBase< _MatrixType, _UpLo, Derived >::_solve_impl ( const SparseMatrixBase< RhsDerived > &  b,
SparseMatrixBase< DestDerived > &  dest 
) const
inlineinherited
310 {
311 eigen_assert(m_factorizationIsOk && "The decomposition is not in a valid state for solving, you must first call either compute() or symbolic()/numeric()");
312 const Index size = m_cholmodFactor->n;
314 eigen_assert(size==b.rows());
315
316 // note: cs stands for Cholmod Sparse
317 Ref<SparseMatrix<typename RhsDerived::Scalar,ColMajor,typename RhsDerived::StorageIndex> > b_ref(b.const_cast_derived());
318 cholmod_sparse b_cs = viewAsCholmod(b_ref);
319 cholmod_sparse* x_cs = cholmod_spsolve(CHOLMOD_A, m_cholmodFactor, &b_cs, &m_cholmod);
320 if(!x_cs)
321 {
322 this->m_info = NumericalIssue;
323 return;
324 }
325 // TODO optimize this copy by swapping when possible (be careful with alignment, etc.)
326 dest.derived() = viewAsEigen<typename DestDerived::Scalar,ColMajor,typename DestDerived::StorageIndex>(*x_cs);
327 cholmod_free_sparse(&x_cs, &m_cholmod);
328 }

References Eigen::SparseMatrixBase< Derived >::derived(), eigen_assert, EIGEN_UNUSED_VARIABLE, Eigen::CholmodBase< _MatrixType, _UpLo, Derived >::m_cholmod, Eigen::CholmodBase< _MatrixType, _UpLo, Derived >::m_cholmodFactor, Eigen::CholmodBase< _MatrixType, _UpLo, Derived >::m_factorizationIsOk, Eigen::CholmodBase< _MatrixType, _UpLo, Derived >::m_info, Eigen::NumericalIssue, and Eigen::viewAsCholmod().

+ Here is the call graph for this function:

◆ analyzePattern()

template<typename _MatrixType , int _UpLo, typename Derived >
void Eigen::CholmodBase< _MatrixType, _UpLo, Derived >::analyzePattern ( const MatrixType matrix)
inlineinherited

Performs a symbolic decomposition on the sparsity pattern of matrix.

This function is particularly useful when solving for several problems having the same structure.

See also
factorize()
246 {
248 {
249 cholmod_free_factor(&m_cholmodFactor, &m_cholmod);
250 m_cholmodFactor = 0;
251 }
252 cholmod_sparse A = viewAsCholmod(matrix.template selfadjointView<UpLo>());
253 m_cholmodFactor = cholmod_analyze(&A, &m_cholmod);
254
255 this->m_isInitialized = true;
256 this->m_info = Success;
257 m_analysisIsOk = true;
258 m_factorizationIsOk = false;
259 }
int m_analysisIsOk
Definition CholmodSupport.h:408
bool m_isInitialized
Definition SparseSolverBase.h:119
@ Success
Definition Constants.h:432

References Eigen::CholmodBase< _MatrixType, _UpLo, Derived >::m_analysisIsOk, Eigen::CholmodBase< _MatrixType, _UpLo, Derived >::m_cholmod, Eigen::CholmodBase< _MatrixType, _UpLo, Derived >::m_cholmodFactor, Eigen::CholmodBase< _MatrixType, _UpLo, Derived >::m_factorizationIsOk, Eigen::CholmodBase< _MatrixType, _UpLo, Derived >::m_info, Eigen::CholmodBase< _MatrixType, _UpLo, Derived >::m_isInitialized, Eigen::Success, and Eigen::viewAsCholmod().

Referenced by Eigen::CholmodBase< _MatrixType, _UpLo, Derived >::compute().

+ Here is the call graph for this function:
+ Here is the caller graph for this function:

◆ cholmod()

template<typename _MatrixType , int _UpLo, typename Derived >
cholmod_common & Eigen::CholmodBase< _MatrixType, _UpLo, Derived >::cholmod ( )
inlineinherited

Returns a reference to the Cholmod's configuration structure to get a full control over the performed operations. See the Cholmod user guide for details.

280{ return m_cholmod; }

References Eigen::CholmodBase< _MatrixType, _UpLo, Derived >::m_cholmod.

◆ cols()

template<typename _MatrixType , int _UpLo, typename Derived >
StorageIndex Eigen::CholmodBase< _MatrixType, _UpLo, Derived >::cols ( ) const
inlineinherited
217{ return internal::convert_index<StorageIndex, Index>(m_cholmodFactor->n); }

References Eigen::CholmodBase< _MatrixType, _UpLo, Derived >::m_cholmodFactor.

◆ compute()

template<typename _MatrixType , int _UpLo, typename Derived >
Derived & Eigen::CholmodBase< _MatrixType, _UpLo, Derived >::compute ( const MatrixType matrix)
inlineinherited

Computes the sparse Cholesky decomposition of matrix

233 {
234 analyzePattern(matrix);
235 factorize(matrix);
236 return derived();
237 }
void analyzePattern(const MatrixType &matrix)
Definition CholmodSupport.h:245
void factorize(const MatrixType &matrix)
Definition CholmodSupport.h:267
Derived & derived()
Definition SparseSolverBase.h:79

References Eigen::CholmodBase< _MatrixType, _UpLo, Derived >::analyzePattern(), Eigen::CholmodBase< _MatrixType, _UpLo, Derived >::derived(), and Eigen::CholmodBase< _MatrixType, _UpLo, Derived >::factorize().

Referenced by Eigen::CholmodBase< _MatrixType, _UpLo, Derived >::CholmodBase(), Eigen::CholmodDecomposition< _MatrixType, _UpLo >::CholmodDecomposition(), Eigen::CholmodSimplicialLDLT< _MatrixType, _UpLo >::CholmodSimplicialLDLT(), Eigen::CholmodSimplicialLLT< _MatrixType, _UpLo >::CholmodSimplicialLLT(), Eigen::CholmodSupernodalLLT< _MatrixType, _UpLo >::CholmodSupernodalLLT(), and igl::slim::solve_weighted_arap().

+ Here is the call graph for this function:
+ Here is the caller graph for this function:

◆ derived() [1/2]

template<typename _MatrixType , int _UpLo, typename Derived >
Derived & Eigen::SparseSolverBase< Derived >::derived ( )
inlineprotectedinherited
79{ return *static_cast<Derived*>(this); }

Referenced by Eigen::CholmodBase< _MatrixType, _UpLo, Derived >::compute(), and Eigen::CholmodBase< _MatrixType, _UpLo, Derived >::setShift().

+ Here is the caller graph for this function:

◆ derived() [2/2]

template<typename _MatrixType , int _UpLo, typename Derived >
const Derived & Eigen::SparseSolverBase< Derived >::derived ( ) const
inlineprotectedinherited
80{ return *static_cast<const Derived*>(this); }

◆ determinant()

template<typename _MatrixType , int _UpLo, typename Derived >
Scalar Eigen::CholmodBase< _MatrixType, _UpLo, Derived >::determinant ( ) const
inlineinherited
Returns
the determinant of the underlying matrix from the current factorization
349 {
350 using std::exp;
351 return exp(logDeterminant());
352 }
EIGEN_DEVICE_FUNC const ExpReturnType exp() const
Definition ArrayCwiseUnaryOps.h:88
Scalar logDeterminant() const
Definition CholmodSupport.h:355

References exp(), and Eigen::CholmodBase< _MatrixType, _UpLo, Derived >::logDeterminant().

+ Here is the call graph for this function:

◆ dumpMemory()

template<typename _MatrixType , int _UpLo, typename Derived >
template<typename Stream >
void Eigen::CholmodBase< _MatrixType, _UpLo, Derived >::dumpMemory ( Stream &  )
inlineinherited
400 {}

◆ factorize()

template<typename _MatrixType , int _UpLo, typename Derived >
void Eigen::CholmodBase< _MatrixType, _UpLo, Derived >::factorize ( const MatrixType matrix)
inlineinherited

Performs a numeric decomposition of matrix

The given matrix must have the same sparsity pattern as the matrix on which the symbolic decomposition has been performed.

See also
analyzePattern()
268 {
269 eigen_assert(m_analysisIsOk && "You must first call analyzePattern()");
270 cholmod_sparse A = viewAsCholmod(matrix.template selfadjointView<UpLo>());
271 cholmod_factorize_p(&A, m_shiftOffset, 0, 0, m_cholmodFactor, &m_cholmod);
272
273 // If the factorization failed, minor is the column at which it did. On success minor == n.
275 m_factorizationIsOk = true;
276 }
double m_shiftOffset[2]
Definition CholmodSupport.h:405

References eigen_assert, Eigen::CholmodBase< _MatrixType, _UpLo, Derived >::m_analysisIsOk, Eigen::CholmodBase< _MatrixType, _UpLo, Derived >::m_cholmod, Eigen::CholmodBase< _MatrixType, _UpLo, Derived >::m_cholmodFactor, Eigen::CholmodBase< _MatrixType, _UpLo, Derived >::m_factorizationIsOk, Eigen::CholmodBase< _MatrixType, _UpLo, Derived >::m_info, Eigen::CholmodBase< _MatrixType, _UpLo, Derived >::m_shiftOffset, Eigen::NumericalIssue, Eigen::Success, and Eigen::viewAsCholmod().

Referenced by Eigen::CholmodBase< _MatrixType, _UpLo, Derived >::compute().

+ Here is the call graph for this function:
+ Here is the caller graph for this function:

◆ info()

template<typename _MatrixType , int _UpLo, typename Derived >
ComputationInfo Eigen::CholmodBase< _MatrixType, _UpLo, Derived >::info ( ) const
inlineinherited

Reports whether previous computation was successful.

Returns
Success if computation was succesful, NumericalIssue if the matrix.appears to be negative.
226 {
227 eigen_assert(m_isInitialized && "Decomposition is not initialized.");
228 return m_info;
229 }

References eigen_assert, Eigen::CholmodBase< _MatrixType, _UpLo, Derived >::m_info, and Eigen::CholmodBase< _MatrixType, _UpLo, Derived >::m_isInitialized.

◆ init()

template<typename _MatrixType , int _UpLo = Lower>
void Eigen::CholmodDecomposition< _MatrixType, _UpLo >::init ( )
inlineprotected
631 {
632 m_cholmod.final_asis = 1;
633 m_cholmod.supernodal = CHOLMOD_AUTO;
634 }
cholmod_common m_cholmod
Definition CholmodSupport.h:403

References Eigen::CholmodDecomposition< _MatrixType, _UpLo >::m_cholmod.

Referenced by Eigen::CholmodDecomposition< _MatrixType, _UpLo >::CholmodDecomposition(), and Eigen::CholmodDecomposition< _MatrixType, _UpLo >::CholmodDecomposition().

+ Here is the caller graph for this function:

◆ logDeterminant()

template<typename _MatrixType , int _UpLo, typename Derived >
Scalar Eigen::CholmodBase< _MatrixType, _UpLo, Derived >::logDeterminant ( ) const
inlineinherited
Returns
the log determinant of the underlying matrix from the current factorization
356 {
357 using std::log;
358 using numext::real;
359 eigen_assert(m_factorizationIsOk && "The decomposition is not in a valid state for solving, you must first call either compute() or symbolic()/numeric()");
360
361 RealScalar logDet = 0;
362 Scalar *x = static_cast<Scalar*>(m_cholmodFactor->x);
363 if (m_cholmodFactor->is_super)
364 {
365 // Supernodal factorization stored as a packed list of dense column-major blocs,
366 // as described by the following structure:
367
368 // super[k] == index of the first column of the j-th super node
369 StorageIndex *super = static_cast<StorageIndex*>(m_cholmodFactor->super);
370 // pi[k] == offset to the description of row indices
371 StorageIndex *pi = static_cast<StorageIndex*>(m_cholmodFactor->pi);
372 // px[k] == offset to the respective dense block
373 StorageIndex *px = static_cast<StorageIndex*>(m_cholmodFactor->px);
374
375 Index nb_super_nodes = m_cholmodFactor->nsuper;
376 for (Index k=0; k < nb_super_nodes; ++k)
377 {
378 StorageIndex ncols = super[k + 1] - super[k];
379 StorageIndex nrows = pi[k + 1] - pi[k];
380
381 Map<const Array<Scalar,1,Dynamic>, 0, InnerStride<> > sk(x + px[k], ncols, InnerStride<>(nrows+1));
382 logDet += sk.real().log().sum();
383 }
384 }
385 else
386 {
387 // Simplicial factorization stored as standard CSC matrix.
388 StorageIndex *p = static_cast<StorageIndex*>(m_cholmodFactor->p);
390 for (Index k=0; k<size; ++k)
391 logDet += log(real( x[p[k]] ));
392 }
393 if (m_cholmodFactor->is_ll)
394 logDet *= 2.0;
395 return logDet;
396 };
EIGEN_DEVICE_FUNC const LogReturnType log() const
Definition ArrayCwiseUnaryOps.h:105
EIGEN_DEVICE_FUNC RealReturnType real() const
Definition CommonCwiseUnaryOps.h:86
MatrixType::RealScalar RealScalar
Definition CholmodSupport.h:183
MatrixType::StorageIndex StorageIndex
Definition CholmodSupport.h:185
const double pi
Definition agg_basics.h:268
TCoord< P > x(const P &p)
Definition geometry_traits.hpp:297

References eigen_assert, log(), Eigen::CholmodBase< _MatrixType, _UpLo, Derived >::m_cholmodFactor, Eigen::CholmodBase< _MatrixType, _UpLo, Derived >::m_factorizationIsOk, and real().

Referenced by Eigen::CholmodBase< _MatrixType, _UpLo, Derived >::determinant().

+ Here is the call graph for this function:
+ Here is the caller graph for this function:

◆ rows()

template<typename _MatrixType , int _UpLo, typename Derived >
StorageIndex Eigen::CholmodBase< _MatrixType, _UpLo, Derived >::rows ( ) const
inlineinherited
218{ return internal::convert_index<StorageIndex, Index>(m_cholmodFactor->n); }

References Eigen::CholmodBase< _MatrixType, _UpLo, Derived >::m_cholmodFactor.

◆ setMode()

template<typename _MatrixType , int _UpLo = Lower>
void Eigen::CholmodDecomposition< _MatrixType, _UpLo >::setMode ( CholmodMode  mode)
inline
605 {
606 switch(mode)
607 {
608 case CholmodAuto:
609 m_cholmod.final_asis = 1;
610 m_cholmod.supernodal = CHOLMOD_AUTO;
611 break;
613 m_cholmod.final_asis = 0;
614 m_cholmod.supernodal = CHOLMOD_SIMPLICIAL;
615 m_cholmod.final_ll = 1;
616 break;
618 m_cholmod.final_asis = 1;
619 m_cholmod.supernodal = CHOLMOD_SUPERNODAL;
620 break;
621 case CholmodLDLt:
622 m_cholmod.final_asis = 1;
623 m_cholmod.supernodal = CHOLMOD_SIMPLICIAL;
624 break;
625 default:
626 break;
627 }
628 }
@ CholmodSimplicialLLt
Definition CholmodSupport.h:163
@ CholmodAuto
Definition CholmodSupport.h:163
@ CholmodLDLt
Definition CholmodSupport.h:163
@ CholmodSupernodalLLt
Definition CholmodSupport.h:163

References Eigen::CholmodAuto, Eigen::CholmodLDLt, Eigen::CholmodSimplicialLLt, Eigen::CholmodSupernodalLLt, and Eigen::CholmodDecomposition< _MatrixType, _UpLo >::m_cholmod.

◆ setShift()

template<typename _MatrixType , int _UpLo, typename Derived >
Derived & Eigen::CholmodBase< _MatrixType, _UpLo, Derived >::setShift ( const RealScalar offset)
inlineinherited

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

During the numerical factorization, an offset term is added to the diagonal coefficients:
d_ii = offset + d_ii

The default is offset=0.

Returns
a reference to *this.
342 {
343 m_shiftOffset[0] = double(offset);
344 return derived();
345 }

References Eigen::CholmodBase< _MatrixType, _UpLo, Derived >::derived(), and Eigen::CholmodBase< _MatrixType, _UpLo, Derived >::m_shiftOffset.

+ Here is the call graph for this function:

◆ solve() [1/2]

template<typename Derived >
template<typename Rhs >
const Solve< Derived, Rhs > Eigen::SparseSolverBase< Derived >::solve ( const MatrixBase< Rhs > &  b) const
inlineinherited
Returns
an expression of the solution x of $ A x = b $ using the current decomposition of A.
See also
compute()
89 {
90 eigen_assert(m_isInitialized && "Solver is not initialized.");
91 eigen_assert(derived().rows()==b.rows() && "solve(): invalid number of rows of the right hand side matrix b");
92 return Solve<Derived, Rhs>(derived(), b.derived());
93 }
bool m_isInitialized
Definition SparseSolverBase.h:119
size_t rows(const T &raster)
Definition MarchingSquares.hpp:55

References Eigen::SparseSolverBase< Derived >::derived(), eigen_assert, and Eigen::SparseSolverBase< Derived >::m_isInitialized.

Referenced by igl::embree::bone_heat(), igl::Frame_field_deformer::compute_optimal_positions(), igl::eigs(), igl::copyleft::comiso::FrameInterpolator::interpolateSymmetric(), igl::slim::solve_weighted_arap(), and igl::copyleft::comiso::NRosyField::solveNoRoundings().

+ Here is the call graph for this function:
+ Here is the caller graph for this function:

◆ solve() [2/2]

template<typename Derived >
template<typename Rhs >
const Solve< Derived, Rhs > Eigen::SparseSolverBase< Derived >::solve ( const SparseMatrixBase< Rhs > &  b) const
inlineinherited
Returns
an expression of the solution x of $ A x = b $ using the current decomposition of A.
See also
compute()
102 {
103 eigen_assert(m_isInitialized && "Solver is not initialized.");
104 eigen_assert(derived().rows()==b.rows() && "solve(): invalid number of rows of the right hand side matrix b");
105 return Solve<Derived, Rhs>(derived(), b.derived());
106 }

Member Data Documentation

◆ m_analysisIsOk

template<typename _MatrixType , int _UpLo, typename Derived >
int Eigen::CholmodBase< _MatrixType, _UpLo, Derived >::m_analysisIsOk
protectedinherited

◆ m_cholmod

template<typename _MatrixType , int _UpLo = Lower>
cholmod_common Eigen::CholmodBase< _MatrixType, _UpLo, Derived >::m_cholmod
mutableprivate

◆ m_cholmodFactor

◆ m_factorizationIsOk

◆ m_info

◆ m_isInitialized

template<typename _MatrixType , int _UpLo, typename Derived >
bool Eigen::SparseSolverBase< Derived >::m_isInitialized
mutableprotectedinherited

◆ m_shiftOffset


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