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
-
| _MatrixType | the type of the sparse matrix A, it must be a SparseMatrix<> |
| _UpLo | the 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
template<typename _MatrixType , int _UpLo, typename Derived >
template<typename Rhs , typename Dest >
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()");
291
292
293 Ref<const Matrix<typename Rhs::Scalar,Dynamic,Dynamic,ColMajor> > b_ref(
b.derived());
294
297 if(!x_cd)
298 {
300 return;
301 }
302
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
friend class Eigen::Map
Definition PlainObjectBase.h:121
@ 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().
template<typename _MatrixType , int _UpLo, typename Derived >
template<typename RhsDerived , typename DestDerived >
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()");
315
316
317 Ref<SparseMatrix<typename RhsDerived::Scalar,ColMajor,typename RhsDerived::StorageIndex> > b_ref(
b.const_cast_derived());
320 if(!x_cs)
321 {
323 return;
324 }
325
326 dest.derived() = viewAsEigen<typename DestDerived::Scalar,ColMajor,typename DestDerived::StorageIndex>(*x_cs);
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().
template<typename _MatrixType , int _UpLo, typename Derived >
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 {
251 }
252 cholmod_sparse A =
viewAsCholmod(matrix.template selfadjointView<UpLo>());
254
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().
template<typename _MatrixType , int _UpLo, typename Derived >
Computes the sparse Cholesky decomposition of matrix
233 {
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().
template<typename _MatrixType , int _UpLo, typename Derived >
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 {
270 cholmod_sparse A =
viewAsCholmod(matrix.template selfadjointView<UpLo>());
272
273
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().
template<typename _MatrixType , int _UpLo, typename Derived >
- 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
364 {
365
366
367
368
370
372
374
376 for (
Index k=0; k < nb_super_nodes; ++k)
377 {
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
391 logDet +=
log(
real( x[p[k]] ));
392 }
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().
template<typename _MatrixType , int _UpLo, typename Derived >