10 #ifndef EIGEN_CHOLMODSUPPORT_H 11 #define EIGEN_CHOLMODSUPPORT_H 17 template<
typename Scalar,
typename CholmodType>
18 void cholmod_configure_matrix(CholmodType& mat)
20 if (internal::is_same<Scalar,float>::value)
22 mat.xtype = CHOLMOD_REAL;
23 mat.dtype = CHOLMOD_SINGLE;
25 else if (internal::is_same<Scalar,double>::value)
27 mat.xtype = CHOLMOD_REAL;
28 mat.dtype = CHOLMOD_DOUBLE;
30 else if (internal::is_same<Scalar,std::complex<float> >::value)
32 mat.xtype = CHOLMOD_COMPLEX;
33 mat.dtype = CHOLMOD_SINGLE;
35 else if (internal::is_same<Scalar,std::complex<double> >::value)
37 mat.xtype = CHOLMOD_COMPLEX;
38 mat.dtype = CHOLMOD_DOUBLE;
42 eigen_assert(
false &&
"Scalar type not supported by CHOLMOD");
51 template<
typename _Scalar,
int _Options,
typename _Index>
52 cholmod_sparse viewAsCholmod(SparseMatrix<_Scalar,_Options,_Index>& mat)
55 res.nzmax = mat.nonZeros();
56 res.nrow = mat.rows();;
57 res.ncol = mat.cols();
58 res.p = mat.outerIndexPtr();
59 res.i = mat.innerIndexPtr();
60 res.x = mat.valuePtr();
63 if(mat.isCompressed())
71 res.nz = mat.innerNonZeroPtr();
77 if (internal::is_same<_Index,int>::value)
79 res.itype = CHOLMOD_INT;
81 else if (internal::is_same<_Index,SuiteSparse_long>::value)
83 res.itype = CHOLMOD_LONG;
87 eigen_assert(
false &&
"Index type not supported yet");
91 internal::cholmod_configure_matrix<_Scalar>(res);
98 template<
typename _Scalar,
int _Options,
typename _Index>
99 const cholmod_sparse viewAsCholmod(
const SparseMatrix<_Scalar,_Options,_Index>& mat)
101 cholmod_sparse res = viewAsCholmod(mat.const_cast_derived());
107 template<
typename _Scalar,
int _Options,
typename _Index,
unsigned int UpLo>
108 cholmod_sparse viewAsCholmod(
const SparseSelfAdjointView<SparseMatrix<_Scalar,_Options,_Index>, UpLo>& mat)
110 cholmod_sparse res = viewAsCholmod(mat.matrix().const_cast_derived());
112 if(UpLo==
Upper) res.stype = 1;
113 if(UpLo==
Lower) res.stype = -1;
120 template<
typename Derived>
121 cholmod_dense viewAsCholmod(MatrixBase<Derived>& mat)
123 EIGEN_STATIC_ASSERT((internal::traits<Derived>::Flags&
RowMajorBit)==0,THIS_METHOD_IS_ONLY_FOR_COLUMN_MAJOR_MATRICES);
124 typedef typename Derived::Scalar Scalar;
127 res.nrow = mat.rows();
128 res.ncol = mat.cols();
129 res.nzmax = res.nrow * res.ncol;
130 res.d = Derived::IsVectorAtCompileTime ? mat.derived().size() : mat.derived().outerStride();
131 res.x = (
void*)(mat.derived().data());
134 internal::cholmod_configure_matrix<Scalar>(res);
141 template<
typename Scalar,
int Flags,
typename Index>
142 MappedSparseMatrix<Scalar,Flags,Index> viewAsEigen(cholmod_sparse& cm)
144 return MappedSparseMatrix<Scalar,Flags,Index>
145 (cm.nrow, cm.ncol,
static_cast<Index*
>(cm.p)[cm.ncol],
146 static_cast<Index*>(cm.p),
static_cast<Index*
>(cm.i),static_cast<Scalar*>(cm.x) );
150 CholmodAuto, CholmodSimplicialLLt, CholmodSupernodalLLt, CholmodLDLt
159 template<
typename _MatrixType,
int _UpLo,
typename Derived>
163 typedef _MatrixType MatrixType;
164 enum { UpLo = _UpLo };
165 typedef typename MatrixType::Scalar Scalar;
166 typedef typename MatrixType::RealScalar RealScalar;
167 typedef MatrixType CholMatrixType;
168 typedef typename MatrixType::Index Index;
173 : m_cholmodFactor(0), m_info(
Success), m_isInitialized(
false)
175 m_shiftOffset[0] = m_shiftOffset[1] = RealScalar(0.0);
176 cholmod_start(&m_cholmod);
180 : m_cholmodFactor(0), m_info(
Success), m_isInitialized(
false)
182 m_shiftOffset[0] = m_shiftOffset[1] = RealScalar(0.0);
183 cholmod_start(&m_cholmod);
190 cholmod_free_factor(&m_cholmodFactor, &m_cholmod);
191 cholmod_finish(&m_cholmod);
194 inline Index cols()
const {
return m_cholmodFactor->n; }
195 inline Index rows()
const {
return m_cholmodFactor->n; }
197 Derived& derived() {
return *
static_cast<Derived*
>(
this); }
198 const Derived& derived()
const {
return *
static_cast<const Derived*
>(
this); }
207 eigen_assert(m_isInitialized &&
"Decomposition is not initialized.");
214 analyzePattern(matrix);
223 template<
typename Rhs>
224 inline const internal::solve_retval<CholmodBase, Rhs>
227 eigen_assert(m_isInitialized &&
"LLT is not initialized.");
228 eigen_assert(rows()==b.rows()
229 &&
"CholmodDecomposition::solve(): invalid number of rows of the right hand side matrix b");
230 return internal::solve_retval<CholmodBase, Rhs>(*
this, b.derived());
237 template<
typename Rhs>
238 inline const internal::sparse_solve_retval<CholmodBase, Rhs>
241 eigen_assert(m_isInitialized &&
"LLT is not initialized.");
242 eigen_assert(rows()==b.
rows()
243 &&
"CholmodDecomposition::solve(): invalid number of rows of the right hand side matrix b");
244 return internal::sparse_solve_retval<CholmodBase, Rhs>(*
this, b.
derived());
257 cholmod_free_factor(&m_cholmodFactor, &m_cholmod);
260 cholmod_sparse A = viewAsCholmod(matrix.template selfadjointView<UpLo>());
261 m_cholmodFactor = cholmod_analyze(&A, &m_cholmod);
263 this->m_isInitialized =
true;
265 m_analysisIsOk =
true;
266 m_factorizationIsOk =
false;
277 eigen_assert(m_analysisIsOk &&
"You must first call analyzePattern()");
278 cholmod_sparse A = viewAsCholmod(matrix.template selfadjointView<UpLo>());
279 cholmod_factorize_p(&A, m_shiftOffset, 0, 0, m_cholmodFactor, &m_cholmod);
283 m_factorizationIsOk =
true;
288 cholmod_common&
cholmod() {
return m_cholmod; }
290 #ifndef EIGEN_PARSED_BY_DOXYGEN 292 template<
typename Rhs,
typename Dest>
295 eigen_assert(m_factorizationIsOk &&
"The decomposition is not in a valid state for solving, you must first call either compute() or symbolic()/numeric()");
296 const Index size = m_cholmodFactor->n;
297 EIGEN_UNUSED_VARIABLE(size);
298 eigen_assert(size==b.rows());
301 Rhs& b_ref(b.const_cast_derived());
302 cholmod_dense b_cd = viewAsCholmod(b_ref);
303 cholmod_dense* x_cd = cholmod_solve(CHOLMOD_A, m_cholmodFactor, &b_cd, &m_cholmod);
310 cholmod_free_dense(&x_cd, &m_cholmod);
314 template<
typename RhsScalar,
int RhsOptions,
typename RhsIndex,
typename DestScalar,
int DestOptions,
typename DestIndex>
317 eigen_assert(m_factorizationIsOk &&
"The decomposition is not in a valid state for solving, you must first call either compute() or symbolic()/numeric()");
318 const Index size = m_cholmodFactor->n;
319 EIGEN_UNUSED_VARIABLE(size);
320 eigen_assert(size==b.
rows());
323 cholmod_sparse b_cs = viewAsCholmod(b);
324 cholmod_sparse* x_cs = cholmod_spsolve(CHOLMOD_A, m_cholmodFactor, &b_cs, &m_cholmod);
330 dest = viewAsEigen<DestScalar,DestOptions,DestIndex>(*x_cs);
331 cholmod_free_sparse(&x_cs, &m_cholmod);
333 #endif // EIGEN_PARSED_BY_DOXYGEN 347 m_shiftOffset[0] = offset;
351 template<
typename Stream>
352 void dumpMemory(Stream& )
356 mutable cholmod_common m_cholmod;
357 cholmod_factor* m_cholmodFactor;
358 RealScalar m_shiftOffset[2];
360 bool m_isInitialized;
361 int m_factorizationIsOk;
383 template<
typename _MatrixType,
int _UpLo = Lower>
387 using Base::m_cholmod;
391 typedef _MatrixType MatrixType;
398 Base::compute(matrix);
405 m_cholmod.final_asis = 0;
406 m_cholmod.supernodal = CHOLMOD_SIMPLICIAL;
407 m_cholmod.final_ll = 1;
430 template<
typename _MatrixType,
int _UpLo = Lower>
434 using Base::m_cholmod;
438 typedef _MatrixType MatrixType;
445 Base::compute(matrix);
452 m_cholmod.final_asis = 1;
453 m_cholmod.supernodal = CHOLMOD_SIMPLICIAL;
475 template<
typename _MatrixType,
int _UpLo = Lower>
479 using Base::m_cholmod;
483 typedef _MatrixType MatrixType;
490 Base::compute(matrix);
497 m_cholmod.final_asis = 1;
498 m_cholmod.supernodal = CHOLMOD_SUPERNODAL;
522 template<
typename _MatrixType,
int _UpLo = Lower>
526 using Base::m_cholmod;
530 typedef _MatrixType MatrixType;
537 Base::compute(matrix);
542 void setMode(CholmodMode mode)
547 m_cholmod.final_asis = 1;
548 m_cholmod.supernodal = CHOLMOD_AUTO;
550 case CholmodSimplicialLLt:
551 m_cholmod.final_asis = 0;
552 m_cholmod.supernodal = CHOLMOD_SIMPLICIAL;
553 m_cholmod.final_ll = 1;
555 case CholmodSupernodalLLt:
556 m_cholmod.final_asis = 1;
557 m_cholmod.supernodal = CHOLMOD_SUPERNODAL;
560 m_cholmod.final_asis = 1;
561 m_cholmod.supernodal = CHOLMOD_SIMPLICIAL;
570 m_cholmod.final_asis = 1;
571 m_cholmod.supernodal = CHOLMOD_AUTO;
577 template<
typename _MatrixType,
int _UpLo,
typename Derived,
typename Rhs>
578 struct solve_retval<CholmodBase<_MatrixType,_UpLo,Derived>, Rhs>
579 : solve_retval_base<CholmodBase<_MatrixType,_UpLo,Derived>, Rhs>
582 EIGEN_MAKE_SOLVE_HELPERS(Dec,Rhs)
584 template<
typename Dest>
void evalTo(Dest& dst)
const 586 dec()._solve(rhs(),dst);
590 template<
typename _MatrixType,
int _UpLo,
typename Derived,
typename Rhs>
591 struct sparse_solve_retval<CholmodBase<_MatrixType,_UpLo,Derived>, Rhs>
592 : sparse_solve_retval_base<CholmodBase<_MatrixType,_UpLo,Derived>, Rhs>
595 EIGEN_MAKE_SPARSE_SOLVE_HELPERS(Dec,Rhs)
597 template<
typename Dest>
void evalTo(Dest& dst)
const 599 dec()._solve(rhs(),dst);
607 #endif // EIGEN_CHOLMODSUPPORT_H void factorize(const MatrixType &matrix)
Definition: CholmodSupport.h:275
Index rows() const
Definition: SparseMatrix.h:119
void analyzePattern(const MatrixType &matrix)
Definition: CholmodSupport.h:253
Definition: Constants.h:167
A versatible sparse matrix representation.
Definition: SparseMatrix.h:85
ComputationInfo info() const
Reports whether previous computation was successful.
Definition: CholmodSupport.h:205
Definition: Constants.h:378
const internal::sparse_solve_retval< CholmodBase, Rhs > solve(const SparseMatrixBase< Rhs > &b) const
Definition: CholmodSupport.h:239
A supernodal Cholesky (LLT) factorization and solver based on Cholmod.
Definition: CholmodSupport.h:476
Definition: Constants.h:169
Base class of any sparse matrices or sparse expressions.
Definition: ForwardDeclarations.h:239
Derived & setShift(const RealScalar &offset)
Definition: CholmodSupport.h:345
Derived & derived()
Definition: EigenBase.h:34
A general Cholesky factorization and solver based on Cholmod.
Definition: CholmodSupport.h:523
const internal::solve_retval< CholmodBase, Rhs > solve(const MatrixBase< Rhs > &b) const
Definition: CholmodSupport.h:225
The base class for the direct Cholesky factorization of Cholmod.
Definition: CholmodSupport.h:160
A simplicial direct Cholesky (LLT) factorization and solver based on Cholmod.
Definition: CholmodSupport.h:384
Definition: Eigen_Colamd.h:50
A simplicial direct Cholesky (LDLT) factorization and solver based on Cholmod.
Definition: CholmodSupport.h:431
Definition: Constants.h:376
const unsigned int RowMajorBit
Definition: Constants.h:53
The matrix class, also used for vectors and row-vectors.
Definition: Matrix.h:127
Index rows() const
Definition: SparseMatrixBase.h:159
ComputationInfo
Definition: Constants.h:374
Derived & compute(const MatrixType &matrix)
Definition: CholmodSupport.h:212
Base class for all dense matrices, vectors, and expressions.
Definition: MatrixBase.h:48
cholmod_common & cholmod()
Definition: CholmodSupport.h:288