10 #ifndef EIGEN_SUPERLUSUPPORT_H 11 #define EIGEN_SUPERLUSUPPORT_H 15 #if defined(SUPERLU_MAJOR_VERSION) && (SUPERLU_MAJOR_VERSION >= 5) 16 #define DECL_GSSVX(PREFIX,FLOATTYPE,KEYTYPE) \ 18 extern void PREFIX##gssvx(superlu_options_t *, SuperMatrix *, int *, int *, int *, \ 19 char *, FLOATTYPE *, FLOATTYPE *, SuperMatrix *, SuperMatrix *, \ 20 void *, int, SuperMatrix *, SuperMatrix *, \ 21 FLOATTYPE *, FLOATTYPE *, FLOATTYPE *, FLOATTYPE *, \ 22 GlobalLU_t *, mem_usage_t *, SuperLUStat_t *, int *); \ 24 inline float SuperLU_gssvx(superlu_options_t *options, SuperMatrix *A, \ 25 int *perm_c, int *perm_r, int *etree, char *equed, \ 26 FLOATTYPE *R, FLOATTYPE *C, SuperMatrix *L, \ 27 SuperMatrix *U, void *work, int lwork, \ 28 SuperMatrix *B, SuperMatrix *X, \ 29 FLOATTYPE *recip_pivot_growth, \ 30 FLOATTYPE *rcond, FLOATTYPE *ferr, FLOATTYPE *berr, \ 31 SuperLUStat_t *stats, int *info, KEYTYPE) { \ 32 mem_usage_t mem_usage; \ 34 PREFIX##gssvx(options, A, perm_c, perm_r, etree, equed, R, C, L, \ 35 U, work, lwork, B, X, recip_pivot_growth, rcond, \ 36 ferr, berr, &gLU, &mem_usage, stats, info); \ 37 return mem_usage.for_lu; \ 39 #else // version < 5.0 40 #define DECL_GSSVX(PREFIX,FLOATTYPE,KEYTYPE) \ 42 extern void PREFIX##gssvx(superlu_options_t *, SuperMatrix *, int *, int *, int *, \ 43 char *, FLOATTYPE *, FLOATTYPE *, SuperMatrix *, SuperMatrix *, \ 44 void *, int, SuperMatrix *, SuperMatrix *, \ 45 FLOATTYPE *, FLOATTYPE *, FLOATTYPE *, FLOATTYPE *, \ 46 mem_usage_t *, SuperLUStat_t *, int *); \ 48 inline float SuperLU_gssvx(superlu_options_t *options, SuperMatrix *A, \ 49 int *perm_c, int *perm_r, int *etree, char *equed, \ 50 FLOATTYPE *R, FLOATTYPE *C, SuperMatrix *L, \ 51 SuperMatrix *U, void *work, int lwork, \ 52 SuperMatrix *B, SuperMatrix *X, \ 53 FLOATTYPE *recip_pivot_growth, \ 54 FLOATTYPE *rcond, FLOATTYPE *ferr, FLOATTYPE *berr, \ 55 SuperLUStat_t *stats, int *info, KEYTYPE) { \ 56 mem_usage_t mem_usage; \ 57 PREFIX##gssvx(options, A, perm_c, perm_r, etree, equed, R, C, L, \ 58 U, work, lwork, B, X, recip_pivot_growth, rcond, \ 59 ferr, berr, &mem_usage, stats, info); \ 60 return mem_usage.for_lu; \ 64 DECL_GSSVX(s,
float,
float)
65 DECL_GSSVX(c,
float,std::complex<float>)
66 DECL_GSSVX(d,
double,
double)
67 DECL_GSSVX(z,
double,std::complex<double>)
70 #define EIGEN_SUPERLU_HAS_ILU 73 #ifdef EIGEN_SUPERLU_HAS_ILU 76 #define DECL_GSISX(PREFIX,FLOATTYPE,KEYTYPE) \ 78 extern void PREFIX##gsisx(superlu_options_t *, SuperMatrix *, int *, int *, int *, \ 79 char *, FLOATTYPE *, FLOATTYPE *, SuperMatrix *, SuperMatrix *, \ 80 void *, int, SuperMatrix *, SuperMatrix *, FLOATTYPE *, FLOATTYPE *, \ 81 mem_usage_t *, SuperLUStat_t *, int *); \ 83 inline float SuperLU_gsisx(superlu_options_t *options, SuperMatrix *A, \ 84 int *perm_c, int *perm_r, int *etree, char *equed, \ 85 FLOATTYPE *R, FLOATTYPE *C, SuperMatrix *L, \ 86 SuperMatrix *U, void *work, int lwork, \ 87 SuperMatrix *B, SuperMatrix *X, \ 88 FLOATTYPE *recip_pivot_growth, \ 90 SuperLUStat_t *stats, int *info, KEYTYPE) { \ 91 mem_usage_t mem_usage; \ 92 PREFIX##gsisx(options, A, perm_c, perm_r, etree, equed, R, C, L, \ 93 U, work, lwork, B, X, recip_pivot_growth, rcond, \ 94 &mem_usage, stats, info); \ 95 return mem_usage.for_lu; \ 98 DECL_GSISX(s,
float,
float)
99 DECL_GSISX(c,
float,std::complex<float>)
100 DECL_GSISX(d,
double,
double)
101 DECL_GSISX(z,
double,std::complex<double>)
105 template<
typename MatrixType>
106 struct SluMatrixMapHelper;
115 struct SluMatrix : SuperMatrix
122 SluMatrix(
const SluMatrix& other)
126 storage = other.storage;
129 SluMatrix& operator=(
const SluMatrix& other)
131 SuperMatrix::operator=(static_cast<const SuperMatrix&>(other));
133 storage = other.storage;
139 union {
int nnz;
int lda;};
145 void setStorageType(Stype_t t)
148 if (t==SLU_NC || t==SLU_NR || t==SLU_DN)
152 eigen_assert(
false &&
"storage type not supported");
157 template<
typename Scalar>
160 if (internal::is_same<Scalar,float>::value)
162 else if (internal::is_same<Scalar,double>::value)
164 else if (internal::is_same<Scalar,std::complex<float> >::value)
166 else if (internal::is_same<Scalar,std::complex<double> >::value)
170 eigen_assert(
false &&
"Scalar type not supported by SuperLU");
174 template<
typename MatrixType>
175 static SluMatrix Map(MatrixBase<MatrixType>& _mat)
177 MatrixType& mat(_mat.derived());
178 eigen_assert( ((MatrixType::Flags&
RowMajorBit)!=RowMajorBit) &&
"row-major dense matrices are not supported by SuperLU");
180 res.setStorageType(SLU_DN);
181 res.setScalarType<
typename MatrixType::Scalar>();
184 res.nrow = internal::convert_index<int>(mat.rows());
185 res.ncol = internal::convert_index<int>(mat.cols());
187 res.storage.lda = internal::convert_index<int>(MatrixType::IsVectorAtCompileTime ? mat.size() : mat.outerStride());
188 res.storage.values = (
void*)(mat.data());
192 template<
typename MatrixType>
193 static SluMatrix Map(SparseMatrixBase<MatrixType>& a_mat)
195 MatrixType &mat(a_mat.derived());
199 res.setStorageType(SLU_NR);
200 res.nrow = internal::convert_index<int>(mat.cols());
201 res.ncol = internal::convert_index<int>(mat.rows());
205 res.setStorageType(SLU_NC);
206 res.nrow = internal::convert_index<int>(mat.rows());
207 res.ncol = internal::convert_index<int>(mat.cols());
212 res.storage.nnz = internal::convert_index<int>(mat.nonZeros());
213 res.storage.values = mat.valuePtr();
214 res.storage.innerInd = mat.innerIndexPtr();
215 res.storage.outerInd = mat.outerIndexPtr();
217 res.setScalarType<
typename MatrixType::Scalar>();
220 if (
int(MatrixType::Flags) & int(
Upper))
222 if (
int(MatrixType::Flags) & int(
Lower))
225 eigen_assert(((
int(MatrixType::Flags) &
int(
SelfAdjoint))==0) &&
"SelfAdjoint matrix shape not supported by SuperLU");
231 template<
typename Scalar,
int Rows,
int Cols,
int Options,
int MRows,
int MCols>
232 struct SluMatrixMapHelper<Matrix<Scalar,Rows,Cols,Options,MRows,MCols> >
234 typedef Matrix<Scalar,Rows,Cols,Options,MRows,MCols> MatrixType;
235 static void run(MatrixType& mat, SluMatrix& res)
237 eigen_assert( ((Options&
RowMajor)!=RowMajor) &&
"row-major dense matrices is not supported by SuperLU");
238 res.setStorageType(SLU_DN);
239 res.setScalarType<Scalar>();
242 res.nrow = mat.rows();
243 res.ncol = mat.cols();
245 res.storage.lda = mat.outerStride();
246 res.storage.values = mat.data();
250 template<
typename Derived>
251 struct SluMatrixMapHelper<SparseMatrixBase<Derived> >
253 typedef Derived MatrixType;
254 static void run(MatrixType& mat, SluMatrix& res)
256 if ((MatrixType::Flags&RowMajorBit)==RowMajorBit)
258 res.setStorageType(SLU_NR);
259 res.nrow = mat.cols();
260 res.ncol = mat.rows();
264 res.setStorageType(SLU_NC);
265 res.nrow = mat.rows();
266 res.ncol = mat.cols();
271 res.storage.nnz = mat.nonZeros();
272 res.storage.values = mat.valuePtr();
273 res.storage.innerInd = mat.innerIndexPtr();
274 res.storage.outerInd = mat.outerIndexPtr();
276 res.setScalarType<
typename MatrixType::Scalar>();
279 if (MatrixType::Flags &
Upper)
281 if (MatrixType::Flags &
Lower)
284 eigen_assert(((MatrixType::Flags &
SelfAdjoint)==0) &&
"SelfAdjoint matrix shape not supported by SuperLU");
290 template<
typename MatrixType>
291 SluMatrix asSluMatrix(MatrixType& mat)
293 return SluMatrix::Map(mat);
297 template<
typename Scalar,
int Flags,
typename Index>
298 MappedSparseMatrix<Scalar,Flags,Index> map_superlu(SluMatrix& sluMat)
300 eigen_assert(((Flags&
RowMajor)==RowMajor && sluMat.Stype == SLU_NR)
301 || ((Flags&
ColMajor)==ColMajor && sluMat.Stype == SLU_NC));
305 return MappedSparseMatrix<Scalar,Flags,Index>(
306 sluMat.nrow, sluMat.ncol, sluMat.storage.outerInd[outerSize],
307 sluMat.storage.outerInd, sluMat.storage.innerInd, reinterpret_cast<Scalar*>(sluMat.storage.values) );
316 template<
typename _MatrixType,
typename Derived>
322 using Base::m_isInitialized;
324 typedef _MatrixType MatrixType;
325 typedef typename MatrixType::Scalar Scalar;
326 typedef typename MatrixType::RealScalar RealScalar;
327 typedef typename MatrixType::StorageIndex StorageIndex;
334 ColsAtCompileTime = MatrixType::ColsAtCompileTime,
335 MaxColsAtCompileTime = MatrixType::MaxColsAtCompileTime
347 inline Index rows()
const {
return m_matrix.rows(); }
348 inline Index cols()
const {
return m_matrix.cols(); }
351 inline superlu_options_t&
options() {
return m_sluOptions; }
360 eigen_assert(m_isInitialized &&
"Decomposition is not initialized.");
367 derived().analyzePattern(matrix);
368 derived().factorize(matrix);
379 m_isInitialized =
true;
381 m_analysisIsOk =
true;
382 m_factorizationIsOk =
false;
385 template<
typename Stream>
386 void dumpMemory(Stream& )
391 void initFactorization(
const MatrixType& a)
393 set_default_options(&this->m_sluOptions);
395 const Index size = a.rows();
398 m_sluA = internal::asSluMatrix(m_matrix);
403 m_sluRscale.resize(size);
404 m_sluCscale.resize(size);
405 m_sluEtree.resize(size);
408 m_sluB.setStorageType(SLU_DN);
409 m_sluB.setScalarType<Scalar>();
410 m_sluB.Mtype = SLU_GE;
411 m_sluB.storage.values = 0;
414 m_sluB.storage.lda = internal::convert_index<int>(size);
417 m_extractedDataAreDirty =
true;
423 m_isInitialized =
false;
428 void extractData()
const;
433 Destroy_SuperNode_Matrix(&m_sluL);
435 Destroy_CompCol_Matrix(&m_sluU);
440 memset(&m_sluL,0,
sizeof m_sluL);
441 memset(&m_sluU,0,
sizeof m_sluU);
445 mutable LUMatrixType m_l;
446 mutable LUMatrixType m_u;
447 mutable IntColVectorType m_p;
448 mutable IntRowVectorType m_q;
450 mutable LUMatrixType m_matrix;
451 mutable SluMatrix m_sluA;
452 mutable SuperMatrix m_sluL, m_sluU;
453 mutable SluMatrix m_sluB, m_sluX;
454 mutable SuperLUStat_t m_sluStat;
455 mutable superlu_options_t m_sluOptions;
456 mutable std::vector<int> m_sluEtree;
459 mutable char m_sluEqued;
462 int m_factorizationIsOk;
464 mutable bool m_extractedDataAreDirty;
487 template<
typename _MatrixType>
492 typedef _MatrixType MatrixType;
493 typedef typename Base::Scalar Scalar;
494 typedef typename Base::RealScalar RealScalar;
495 typedef typename Base::StorageIndex StorageIndex;
504 using Base::_solve_impl;
508 explicit SuperLU(
const MatrixType& matrix) : Base()
511 Base::compute(matrix);
527 m_isInitialized =
false;
528 Base::analyzePattern(matrix);
537 void factorize(
const MatrixType& matrix);
540 template<
typename Rhs,
typename Dest>
543 inline const LMatrixType& matrixL()
const 545 if (m_extractedDataAreDirty) this->extractData();
549 inline const UMatrixType& matrixU()
const 551 if (m_extractedDataAreDirty) this->extractData();
555 inline const IntColVectorType& permutationP()
const 557 if (m_extractedDataAreDirty) this->extractData();
561 inline const IntRowVectorType& permutationQ()
const 563 if (m_extractedDataAreDirty) this->extractData();
567 Scalar determinant()
const;
571 using Base::m_matrix;
572 using Base::m_sluOptions;
578 using Base::m_sluEtree;
579 using Base::m_sluEqued;
580 using Base::m_sluRscale;
581 using Base::m_sluCscale;
584 using Base::m_sluStat;
585 using Base::m_sluFerr;
586 using Base::m_sluBerr;
590 using Base::m_analysisIsOk;
591 using Base::m_factorizationIsOk;
592 using Base::m_extractedDataAreDirty;
593 using Base::m_isInitialized;
600 set_default_options(&this->m_sluOptions);
601 m_sluOptions.PrintStat = NO;
602 m_sluOptions.ConditionNumber = NO;
603 m_sluOptions.Trans = NOTRANS;
604 m_sluOptions.ColPerm = COLAMD;
612 template<
typename MatrixType>
615 eigen_assert(m_analysisIsOk &&
"You must first call analyzePattern()");
622 this->initFactorization(a);
624 m_sluOptions.ColPerm = COLAMD;
626 RealScalar recip_pivot_growth, rcond;
627 RealScalar ferr, berr;
629 StatInit(&m_sluStat);
630 SuperLU_gssvx(&m_sluOptions, &m_sluA, m_q.data(), m_p.data(), &m_sluEtree[0],
631 &m_sluEqued, &m_sluRscale[0], &m_sluCscale[0],
635 &recip_pivot_growth, &rcond,
637 &m_sluStat, &info, Scalar());
638 StatFree(&m_sluStat);
640 m_extractedDataAreDirty =
true;
644 m_factorizationIsOk =
true;
647 template<
typename MatrixType>
648 template<
typename Rhs,
typename Dest>
651 eigen_assert(m_factorizationIsOk &&
"The decomposition is not in a valid state for solving, you must first call either compute() or analyzePattern()/factorize()");
654 eigen_assert(m_matrix.rows()==b.
rows());
656 m_sluOptions.Trans = NOTRANS;
657 m_sluOptions.Fact = FACTORED;
658 m_sluOptions.IterRefine = NOREFINE;
661 m_sluFerr.resize(rhsCols);
662 m_sluBerr.resize(rhsCols);
667 m_sluB = SluMatrix::Map(b_ref.const_cast_derived());
668 m_sluX = SluMatrix::Map(x_ref.const_cast_derived());
670 typename Rhs::PlainObject b_cpy;
674 m_sluB = SluMatrix::Map(b_cpy.const_cast_derived());
677 StatInit(&m_sluStat);
679 RealScalar recip_pivot_growth, rcond;
680 SuperLU_gssvx(&m_sluOptions, &m_sluA,
681 m_q.data(), m_p.data(),
682 &m_sluEtree[0], &m_sluEqued,
683 &m_sluRscale[0], &m_sluCscale[0],
687 &recip_pivot_growth, &rcond,
688 &m_sluFerr[0], &m_sluBerr[0],
689 &m_sluStat, &info, Scalar());
690 StatFree(&m_sluStat);
692 if(x.
derived().data() != x_ref.data())
705 template<
typename MatrixType,
typename Derived>
708 eigen_assert(m_factorizationIsOk &&
"The decomposition is not in a valid state for extracting factors, you must first call either compute() or analyzePattern()/factorize()");
709 if (m_extractedDataAreDirty)
712 int fsupc, istart, nsupr;
713 int lastl = 0, lastu = 0;
714 SCformat *Lstore =
static_cast<SCformat*
>(m_sluL.Store);
715 NCformat *Ustore =
static_cast<NCformat*
>(m_sluU.Store);
718 const Index size = m_matrix.rows();
719 m_l.resize(size,size);
720 m_l.resizeNonZeros(Lstore->nnz);
721 m_u.resize(size,size);
722 m_u.resizeNonZeros(Ustore->nnz);
724 int* Lcol = m_l.outerIndexPtr();
725 int* Lrow = m_l.innerIndexPtr();
726 Scalar* Lval = m_l.valuePtr();
728 int* Ucol = m_u.outerIndexPtr();
729 int* Urow = m_u.innerIndexPtr();
730 Scalar* Uval = m_u.valuePtr();
736 for (
int k = 0; k <= Lstore->nsuper; ++k)
738 fsupc = L_FST_SUPC(k);
739 istart = L_SUB_START(fsupc);
740 nsupr = L_SUB_START(fsupc+1) - istart;
744 for (
int j = fsupc; j < L_FST_SUPC(k+1); ++j)
746 SNptr = &((Scalar*)Lstore->nzval)[L_NZ_START(j)];
749 for (
int i = U_NZ_START(j); i < U_NZ_START(j+1); ++i)
751 Uval[lastu] = ((Scalar*)Ustore->nzval)[i];
753 if (Uval[lastu] != 0.0)
754 Urow[lastu++] = U_SUB(i);
756 for (
int i = 0; i < upper; ++i)
759 Uval[lastu] = SNptr[i];
761 if (Uval[lastu] != 0.0)
762 Urow[lastu++] = L_SUB(istart+i);
768 Lrow[lastl++] = L_SUB(istart + upper - 1);
769 for (
int i = upper; i < nsupr; ++i)
771 Lval[lastl] = SNptr[i];
773 if (Lval[lastl] != 0.0)
774 Lrow[lastl++] = L_SUB(istart+i);
784 m_l.resizeNonZeros(lastl);
785 m_u.resizeNonZeros(lastu);
787 m_extractedDataAreDirty =
false;
791 template<
typename MatrixType>
794 eigen_assert(m_factorizationIsOk &&
"The decomposition is not in a valid state for computing the determinant, you must first call either compute() or analyzePattern()/factorize()");
796 if (m_extractedDataAreDirty)
799 Scalar det = Scalar(1);
800 for (
int j=0; j<m_u.cols(); ++j)
802 if (m_u.outerIndexPtr()[j+1]-m_u.outerIndexPtr()[j] > 0)
804 int lastId = m_u.outerIndexPtr()[j+1]-1;
805 eigen_assert(m_u.innerIndexPtr()[lastId]<=j);
806 if (m_u.innerIndexPtr()[lastId]==j)
807 det *= m_u.valuePtr()[lastId];
813 return det/m_sluRscale.prod()/m_sluCscale.prod();
818 #ifdef EIGEN_PARSED_BY_DOXYGEN 819 #define EIGEN_SUPERLU_HAS_ILU 822 #ifdef EIGEN_SUPERLU_HAS_ILU 840 template<
typename _MatrixType>
845 typedef _MatrixType MatrixType;
846 typedef typename Base::Scalar Scalar;
847 typedef typename Base::RealScalar RealScalar;
850 using Base::_solve_impl;
854 SuperILU(
const MatrixType& matrix) : Base()
857 Base::compute(matrix);
872 Base::analyzePattern(matrix);
881 void factorize(
const MatrixType& matrix);
883 #ifndef EIGEN_PARSED_BY_DOXYGEN 885 template<
typename Rhs,
typename Dest>
887 #endif // EIGEN_PARSED_BY_DOXYGEN 891 using Base::m_matrix;
892 using Base::m_sluOptions;
898 using Base::m_sluEtree;
899 using Base::m_sluEqued;
900 using Base::m_sluRscale;
901 using Base::m_sluCscale;
904 using Base::m_sluStat;
905 using Base::m_sluFerr;
906 using Base::m_sluBerr;
910 using Base::m_analysisIsOk;
911 using Base::m_factorizationIsOk;
912 using Base::m_extractedDataAreDirty;
913 using Base::m_isInitialized;
920 ilu_set_default_options(&m_sluOptions);
921 m_sluOptions.PrintStat = NO;
922 m_sluOptions.ConditionNumber = NO;
923 m_sluOptions.Trans = NOTRANS;
924 m_sluOptions.ColPerm = MMD_AT_PLUS_A;
927 m_sluOptions.ILU_MILU = SILU;
931 m_sluOptions.ILU_DropRule = DROP_BASIC;
939 template<
typename MatrixType>
942 eigen_assert(m_analysisIsOk &&
"You must first call analyzePattern()");
949 this->initFactorization(a);
952 RealScalar recip_pivot_growth, rcond;
954 StatInit(&m_sluStat);
955 SuperLU_gsisx(&m_sluOptions, &m_sluA, m_q.data(), m_p.data(), &m_sluEtree[0],
956 &m_sluEqued, &m_sluRscale[0], &m_sluCscale[0],
960 &recip_pivot_growth, &rcond,
961 &m_sluStat, &info, Scalar());
962 StatFree(&m_sluStat);
966 m_factorizationIsOk =
true;
969 #ifndef EIGEN_PARSED_BY_DOXYGEN 970 template<
typename MatrixType>
971 template<
typename Rhs,
typename Dest>
974 eigen_assert(m_factorizationIsOk &&
"The decomposition is not in a valid state for solving, you must first call either compute() or analyzePattern()/factorize()");
976 const int rhsCols = b.
cols();
977 eigen_assert(m_matrix.rows()==b.
rows());
979 m_sluOptions.Trans = NOTRANS;
980 m_sluOptions.Fact = FACTORED;
981 m_sluOptions.IterRefine = NOREFINE;
983 m_sluFerr.resize(rhsCols);
984 m_sluBerr.resize(rhsCols);
989 m_sluB = SluMatrix::Map(b_ref.const_cast_derived());
990 m_sluX = SluMatrix::Map(x_ref.const_cast_derived());
992 typename Rhs::PlainObject b_cpy;
996 m_sluB = SluMatrix::Map(b_cpy.const_cast_derived());
1000 RealScalar recip_pivot_growth, rcond;
1002 StatInit(&m_sluStat);
1003 SuperLU_gsisx(&m_sluOptions, &m_sluA,
1004 m_q.data(), m_p.data(),
1005 &m_sluEtree[0], &m_sluEqued,
1006 &m_sluRscale[0], &m_sluCscale[0],
1010 &recip_pivot_growth, &rcond,
1011 &m_sluStat, &info, Scalar());
1012 StatFree(&m_sluStat);
1014 if(x.
derived().data() != x_ref.data())
1025 #endif // EIGEN_SUPERLUSUPPORT_H Definition: Constants.h:319
void compute(const MatrixType &matrix)
Definition: SuperLUSupport.h:365
void analyzePattern(const MatrixType &matrix)
Definition: SuperLUSupport.h:870
A sparse direct LU factorization and solver based on the SuperLU library.
Definition: SuperLUSupport.h:488
A matrix or vector expression mapping an existing array of data.
Definition: Map.h:94
A sparse direct incomplete LU factorization and solver based on the SuperLU library.
Definition: SuperLUSupport.h:841
A base class for sparse solvers.
Definition: SparseSolverBase.h:67
Namespace containing all symbols from the Eigen library.
Definition: Core:141
Holds information about the various numeric (i.e. scalar) types allowed by Eigen. ...
Definition: NumTraits.h:232
Derived & derived()
Definition: EigenBase.h:46
const unsigned int RowMajorBit
Definition: Constants.h:66
void analyzePattern(const MatrixType &)
Definition: SuperLUSupport.h:377
Matrix< Type, Size, 1 > Vector
[c++11]
Definition: Matrix.h:551
Definition: Constants.h:209
void factorize(const MatrixType &matrix)
Definition: SuperLUSupport.h:940
Definition: Constants.h:444
EIGEN_DEFAULT_DENSE_INDEX_TYPE Index
The Index type as used for the API.
Definition: Meta.h:74
EIGEN_CONSTEXPR Index cols() const EIGEN_NOEXCEPT
Definition: EigenBase.h:63
The base class for the direct and incomplete LU factorization of SuperLU.
Definition: SuperLUSupport.h:317
Definition: Constants.h:449
Definition: Constants.h:211
Definition: Constants.h:442
A matrix or vector expression mapping an existing expression.
Definition: Ref.h:281
Definition: Eigen_Colamd.h:50
Definition: Constants.h:225
Definition: Constants.h:321
Expression of a triangular part in a matrix.
Definition: TriangularMatrix.h:187
EIGEN_CONSTEXPR Index rows() const EIGEN_NOEXCEPT
Definition: EigenBase.h:60
void factorize(const MatrixType &matrix)
Definition: SuperLUSupport.h:613
superlu_options_t & options()
Definition: SuperLUSupport.h:351
ComputationInfo
Definition: Constants.h:440
Base class for all dense matrices, vectors, and expressions.
Definition: MatrixBase.h:48
void analyzePattern(const MatrixType &matrix)
Definition: SuperLUSupport.h:524
ComputationInfo info() const
Reports whether previous computation was successful.
Definition: SuperLUSupport.h:358