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

A sparse direct LU factorization and solver based on the SuperLU library. More...

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

+ Inheritance diagram for Eigen::SuperLU< _MatrixType >:
+ Collaboration diagram for Eigen::SuperLU< _MatrixType >:

Public Types

typedef SuperLUBase< _MatrixType, SuperLUBase
 
typedef _MatrixType MatrixType
 
typedef Base::Scalar Scalar
 
typedef Base::RealScalar RealScalar
 
typedef Base::StorageIndex StorageIndex
 
typedef Base::IntRowVectorType IntRowVectorType
 
typedef Base::IntColVectorType IntColVectorType
 
typedef Base::PermutationMap PermutationMap
 
typedef Base::LUMatrixType LUMatrixType
 
typedef TriangularView< LUMatrixType, Lower|UnitDiagLMatrixType
 
typedef TriangularView< LUMatrixType, UpperUMatrixType
 
enum  
 
typedef Matrix< Scalar, Dynamic, 1 > Vector
 

Public Member Functions

 SuperLU ()
 
 SuperLU (const MatrixType &matrix)
 
 ~SuperLU ()
 
void analyzePattern (const MatrixType &matrix)
 
void factorize (const MatrixType &matrix)
 
template<typename Rhs , typename Dest >
void _solve_impl (const MatrixBase< Rhs > &b, MatrixBase< Dest > &dest) const
 
const LMatrixTypematrixL () const
 
const UMatrixTypematrixU () const
 
const IntColVectorTypepermutationP () const
 
const IntRowVectorTypepermutationQ () const
 
Scalar determinant () const
 
Index rows () const
 
Index cols () const
 
superlu_options_t & options ()
 
ComputationInfo info () const
 Reports whether previous computation was successful.
 
void compute (const MatrixType &matrix)
 
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
 
template<typename Rhs , typename Dest >
void _solve_impl (const SparseMatrixBase< Rhs > &b, SparseMatrixBase< Dest > &dest) const
 

Protected Member Functions

void init ()
 
void initFactorization (const MatrixType &a)
 
void extractData () const
 
void clearFactors ()
 
SuperLU< _MatrixType > & derived ()
 
const SuperLU< _MatrixType > & derived () const
 

Protected Attributes

LUMatrixType m_matrix
 
superlu_options_t m_sluOptions
 
SluMatrix m_sluA
 
SluMatrix m_sluB
 
SluMatrix m_sluX
 
IntColVectorType m_p
 
IntRowVectorType m_q
 
std::vector< int > m_sluEtree
 
char m_sluEqued
 
Matrix< RealScalar, Dynamic, 1 > m_sluRscale
 
Matrix< RealScalar, Dynamic, 1 > m_sluCscale
 
SuperMatrix m_sluL
 
SuperMatrix m_sluU
 
SuperLUStat_t m_sluStat
 
Matrix< RealScalar, Dynamic, 1 > m_sluFerr
 
Matrix< RealScalar, Dynamic, 1 > m_sluBerr
 
LUMatrixType m_l
 
LUMatrixType m_u
 
int m_analysisIsOk
 
int m_factorizationIsOk
 
bool m_extractedDataAreDirty
 
bool m_isInitialized
 
ComputationInfo m_info
 

Private Member Functions

 SuperLU (SuperLU &)
 

Detailed Description

template<typename _MatrixType>
class Eigen::SuperLU< _MatrixType >

A sparse direct LU factorization and solver based on the SuperLU library.

This class allows to solve for A.X = B sparse linear problems via a direct LU factorization using the SuperLU library. The sparse matrix A must be squared and invertible. The vectors or matrices X and B can be either dense or sparse.

Template Parameters
_MatrixTypethe type of the sparse matrix A, it must be a SparseMatrix<>
Warning
This class is only for the 4.x versions of SuperLU. The 3.x and 5.x versions are not supported.

\implsparsesolverconcept

See also
TutorialSparseSolverConcept, class SparseLU

Member Typedef Documentation

◆ Base

template<typename _MatrixType >
typedef SuperLUBase<_MatrixType,SuperLU> Eigen::SuperLU< _MatrixType >::Base

◆ IntColVectorType

template<typename _MatrixType >
typedef Base::IntColVectorType Eigen::SuperLU< _MatrixType >::IntColVectorType

◆ IntRowVectorType

template<typename _MatrixType >
typedef Base::IntRowVectorType Eigen::SuperLU< _MatrixType >::IntRowVectorType

◆ LMatrixType

template<typename _MatrixType >
typedef TriangularView<LUMatrixType, Lower|UnitDiag> Eigen::SuperLU< _MatrixType >::LMatrixType

◆ LUMatrixType

template<typename _MatrixType >
typedef Base::LUMatrixType Eigen::SuperLU< _MatrixType >::LUMatrixType

◆ MatrixType

template<typename _MatrixType >
typedef _MatrixType Eigen::SuperLU< _MatrixType >::MatrixType

◆ PermutationMap

template<typename _MatrixType >
typedef Base::PermutationMap Eigen::SuperLU< _MatrixType >::PermutationMap

◆ RealScalar

template<typename _MatrixType >
typedef Base::RealScalar Eigen::SuperLU< _MatrixType >::RealScalar

◆ Scalar

template<typename _MatrixType >
typedef Base::Scalar Eigen::SuperLU< _MatrixType >::Scalar

◆ StorageIndex

template<typename _MatrixType >
typedef Base::StorageIndex Eigen::SuperLU< _MatrixType >::StorageIndex

◆ UMatrixType

template<typename _MatrixType >
typedef TriangularView<LUMatrixType, Upper> Eigen::SuperLU< _MatrixType >::UMatrixType

◆ Vector

typedef Matrix<Scalar,Dynamic,1> Eigen::SuperLUBase< _MatrixType, SuperLU< _MatrixType > >::Vector
inherited

Member Enumeration Documentation

◆ anonymous enum

anonymous enum
inherited
333 {
334 ColsAtCompileTime = MatrixType::ColsAtCompileTime,
335 MaxColsAtCompileTime = MatrixType::MaxColsAtCompileTime
336 };
@ MaxColsAtCompileTime
Definition SuperLUSupport.h:335

Constructor & Destructor Documentation

◆ SuperLU() [1/3]

template<typename _MatrixType >
Eigen::SuperLU< _MatrixType >::SuperLU ( )
inline
506: Base() { init(); }
SuperLUBase< _MatrixType, SuperLU > Base
Definition SuperLUSupport.h:491
void init()
Definition SuperLUSupport.h:596

References Eigen::SuperLU< _MatrixType >::init().

+ Here is the call graph for this function:

◆ SuperLU() [2/3]

template<typename _MatrixType >
Eigen::SuperLU< _MatrixType >::SuperLU ( const MatrixType matrix)
inlineexplicit
508 : Base()
509 {
510 init();
511 Base::compute(matrix);
512 }
void compute(const MatrixType &matrix)
Definition SuperLUSupport.h:365

References Eigen::SuperLUBase< _MatrixType, Derived >::compute(), and Eigen::SuperLU< _MatrixType >::init().

+ Here is the call graph for this function:

◆ ~SuperLU()

template<typename _MatrixType >
Eigen::SuperLU< _MatrixType >::~SuperLU ( )
inline
515 {
516 }

◆ SuperLU() [3/3]

template<typename _MatrixType >
Eigen::SuperLU< _MatrixType >::SuperLU ( SuperLU< _MatrixType > &  )
inlineprivate
609{ }

Member Function Documentation

◆ _solve_impl() [1/2]

template<typename MatrixType >
template<typename Rhs , typename Dest >
void Eigen::SuperLU< MatrixType >::_solve_impl ( const MatrixBase< Rhs > &  b,
MatrixBase< Dest > &  dest 
) const
650{
651 eigen_assert(m_factorizationIsOk && "The decomposition is not in a valid state for solving, you must first call either compute() or analyzePattern()/factorize()");
652
653 const Index size = m_matrix.rows();
654 const Index rhsCols = b.cols();
655 eigen_assert(size==b.rows());
656
657 m_sluOptions.Trans = NOTRANS;
658 m_sluOptions.Fact = FACTORED;
659 m_sluOptions.IterRefine = NOREFINE;
660
661
662 m_sluFerr.resize(rhsCols);
663 m_sluBerr.resize(rhsCols);
664
665 Ref<const Matrix<typename Rhs::Scalar,Dynamic,Dynamic,ColMajor> > b_ref(b);
666 Ref<const Matrix<typename Dest::Scalar,Dynamic,Dynamic,ColMajor> > x_ref(x);
667
668 m_sluB = SluMatrix::Map(b_ref.const_cast_derived());
669 m_sluX = SluMatrix::Map(x_ref.const_cast_derived());
670
671 typename Rhs::PlainObject b_cpy;
672 if(m_sluEqued!='N')
673 {
674 b_cpy = b;
675 m_sluB = SluMatrix::Map(b_cpy.const_cast_derived());
676 }
677
678 StatInit(&m_sluStat);
679 int info = 0;
680 RealScalar recip_pivot_growth, rcond;
681 SuperLU_gssvx(&m_sluOptions, &m_sluA,
682 m_q.data(), m_p.data(),
684 &m_sluRscale[0], &m_sluCscale[0],
685 &m_sluL, &m_sluU,
686 NULL, 0,
687 &m_sluB, &m_sluX,
688 &recip_pivot_growth, &rcond,
689 &m_sluFerr[0], &m_sluBerr[0],
690 &m_sluStat, &info, Scalar());
691 StatFree(&m_sluStat);
692
693 if(x.derived().data() != x_ref.data())
694 x = x_ref;
695
697}
#define eigen_assert(x)
Definition Macros.h:579
EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE const Scalar * data() const
Definition PlainObjectBase.h:255
EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void resize(Index rows, Index cols)
Definition PlainObjectBase.h:279
Index rows() const
Definition SparseMatrix.h:136
ComputationInfo info() const
Reports whether previous computation was successful.
Definition SuperLUSupport.h:358
SluMatrix m_sluX
Definition SuperLUSupport.h:453
LUMatrixType m_matrix
Definition SuperLUSupport.h:450
superlu_options_t m_sluOptions
Definition SuperLUSupport.h:455
SluMatrix m_sluB
Definition SuperLUSupport.h:453
Matrix< RealScalar, Dynamic, 1 > m_sluRscale
Definition SuperLUSupport.h:457
Base::RealScalar RealScalar
Definition SuperLUSupport.h:494
Matrix< RealScalar, Dynamic, 1 > m_sluFerr
Definition SuperLUSupport.h:458
ComputationInfo m_info
Definition SuperLUSupport.h:461
std::vector< int > m_sluEtree
Definition SuperLUSupport.h:456
IntRowVectorType m_q
Definition SuperLUSupport.h:448
SuperLUStat_t m_sluStat
Definition SuperLUSupport.h:454
Matrix< RealScalar, Dynamic, 1 > m_sluCscale
Definition SuperLUSupport.h:457
Matrix< RealScalar, Dynamic, 1 > m_sluBerr
Definition SuperLUSupport.h:458
SuperMatrix m_sluL
Definition SuperLUSupport.h:452
SuperMatrix m_sluU
Definition SuperLUSupport.h:452
char m_sluEqued
Definition SuperLUSupport.h:459
SluMatrix m_sluA
Definition SuperLUSupport.h:451
int m_factorizationIsOk
Definition SuperLUSupport.h:462
IntColVectorType m_p
Definition SuperLUSupport.h:447
Base::Scalar Scalar
Definition SuperLUSupport.h:493
@ NumericalIssue
Definition Constants.h:434
@ Success
Definition Constants.h:432
EIGEN_DEFAULT_DENSE_INDEX_TYPE Index
The Index type as used for the API.
Definition Meta.h:33
constexpr auto size(const C &c) -> decltype(c.size())
Definition span.hpp:183
TCoord< P > x(const P &p)
Definition geometry_traits.hpp:297
static SluMatrix Map(MatrixBase< MatrixType > &_mat)
Definition SuperLUSupport.h:175

References eigen_assert, Eigen::SluMatrix::Map(), Eigen::NumericalIssue, and Eigen::Success.

+ Here is the call graph for this function:

◆ _solve_impl() [2/2]

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:

◆ analyzePattern()

template<typename _MatrixType >
void Eigen::SuperLU< _MatrixType >::analyzePattern ( const MatrixType matrix)
inline

Performs a symbolic decomposition on the sparcity of matrix.

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

See also
factorize()
525 {
527 m_isInitialized = false;
528 Base::analyzePattern(matrix);
529 }
void analyzePattern(const MatrixType &)
Definition SuperLUSupport.h:377
bool m_isInitialized
Definition SparseSolverBase.h:119
@ InvalidInput
Definition Constants.h:439

References Eigen::SuperLUBase< _MatrixType, Derived >::analyzePattern(), Eigen::InvalidInput, Eigen::SuperLU< _MatrixType >::m_info, and Eigen::SuperLU< _MatrixType >::m_isInitialized.

+ Here is the call graph for this function:

◆ clearFactors()

void Eigen::SuperLUBase< _MatrixType, SuperLU< _MatrixType > >::clearFactors ( )
inlineprotectedinherited
431 {
432 if(m_sluL.Store)
433 Destroy_SuperNode_Matrix(&m_sluL);
434 if(m_sluU.Store)
435 Destroy_CompCol_Matrix(&m_sluU);
436
437 m_sluL.Store = 0;
438 m_sluU.Store = 0;
439
440 memset(&m_sluL,0,sizeof m_sluL);
441 memset(&m_sluU,0,sizeof m_sluU);
442 }
SuperMatrix m_sluL
Definition SuperLUSupport.h:452
SuperMatrix m_sluU
Definition SuperLUSupport.h:452

◆ cols()

Index Eigen::SuperLUBase< _MatrixType, SuperLU< _MatrixType > >::cols ( ) const
inlineinherited
348{ return m_matrix.cols(); }
Index cols() const
Definition SparseMatrix.h:138
LUMatrixType m_matrix
Definition SuperLUSupport.h:450

◆ compute()

void Eigen::SuperLUBase< _MatrixType, SuperLU< _MatrixType > >::compute ( const MatrixType matrix)
inlineinherited

Computes the sparse Cholesky decomposition of matrix

366 {
367 derived().analyzePattern(matrix);
368 derived().factorize(matrix);
369 }
SuperLU< _MatrixType > & derived()
Definition SparseSolverBase.h:79

◆ derived() [1/2]

SuperLU< _MatrixType > & Eigen::SparseSolverBase< SuperLU< _MatrixType > >::derived ( )
inlineprotectedinherited
79{ return *static_cast<Derived*>(this); }

◆ derived() [2/2]

const SuperLU< _MatrixType > & Eigen::SparseSolverBase< SuperLU< _MatrixType > >::derived ( ) const
inlineprotectedinherited
80{ return *static_cast<const Derived*>(this); }

◆ determinant()

template<typename MatrixType >
SuperLU< MatrixType >::Scalar Eigen::SuperLU< MatrixType >::determinant
794{
795 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
798 this->extractData();
799
800 Scalar det = Scalar(1);
801 for (int j=0; j<m_u.cols(); ++j)
802 {
803 if (m_u.outerIndexPtr()[j+1]-m_u.outerIndexPtr()[j] > 0)
804 {
805 int lastId = m_u.outerIndexPtr()[j+1]-1;
806 eigen_assert(m_u.innerIndexPtr()[lastId]<=j);
807 if (m_u.innerIndexPtr()[lastId]==j)
808 det *= m_u.valuePtr()[lastId];
809 }
810 }
811 if(PermutationMap(m_p.data(),m_p.size()).determinant()*PermutationMap(m_q.data(),m_q.size()).determinant()<0)
812 det = -det;
813 if(m_sluEqued!='N')
814 return det/m_sluRscale.prod()/m_sluCscale.prod();
815 else
816 return det;
817}
const StorageIndex * innerIndexPtr() const
Definition SparseMatrix.h:157
const StorageIndex * outerIndexPtr() const
Definition SparseMatrix.h:166
const Scalar * valuePtr() const
Definition SparseMatrix.h:148
void extractData() const
Definition SuperLUSupport.h:707
bool m_extractedDataAreDirty
Definition SuperLUSupport.h:464
Base::PermutationMap PermutationMap
Definition SuperLUSupport.h:498
LUMatrixType m_u
Definition SuperLUSupport.h:446

References eigen_assert.

◆ dumpMemory()

void Eigen::SuperLUBase< _MatrixType, SuperLU< _MatrixType > >::dumpMemory ( Stream &  )
inlineinherited
387 {}

◆ extractData()

void Eigen::SuperLUBase< MatrixType, SuperLU< _MatrixType > >::extractData
protectedinherited
708{
709 eigen_assert(m_factorizationIsOk && "The decomposition is not in a valid state for extracting factors, you must first call either compute() or analyzePattern()/factorize()");
710 if (m_extractedDataAreDirty)
711 {
712 int upper;
713 int fsupc, istart, nsupr;
714 int lastl = 0, lastu = 0;
715 SCformat *Lstore = static_cast<SCformat*>(m_sluL.Store);
716 NCformat *Ustore = static_cast<NCformat*>(m_sluU.Store);
717 Scalar *SNptr;
718
719 const Index size = m_matrix.rows();
720 m_l.resize(size,size);
721 m_l.resizeNonZeros(Lstore->nnz);
722 m_u.resize(size,size);
723 m_u.resizeNonZeros(Ustore->nnz);
724
725 int* Lcol = m_l.outerIndexPtr();
726 int* Lrow = m_l.innerIndexPtr();
727 Scalar* Lval = m_l.valuePtr();
728
729 int* Ucol = m_u.outerIndexPtr();
730 int* Urow = m_u.innerIndexPtr();
731 Scalar* Uval = m_u.valuePtr();
732
733 Ucol[0] = 0;
734 Ucol[0] = 0;
735
736 /* for each supernode */
737 for (int k = 0; k <= Lstore->nsuper; ++k)
738 {
739 fsupc = L_FST_SUPC(k);
740 istart = L_SUB_START(fsupc);
741 nsupr = L_SUB_START(fsupc+1) - istart;
742 upper = 1;
743
744 /* for each column in the supernode */
745 for (int j = fsupc; j < L_FST_SUPC(k+1); ++j)
746 {
747 SNptr = &((Scalar*)Lstore->nzval)[L_NZ_START(j)];
748
749 /* Extract U */
750 for (int i = U_NZ_START(j); i < U_NZ_START(j+1); ++i)
751 {
752 Uval[lastu] = ((Scalar*)Ustore->nzval)[i];
753 /* Matlab doesn't like explicit zero. */
754 if (Uval[lastu] != 0.0)
755 Urow[lastu++] = U_SUB(i);
756 }
757 for (int i = 0; i < upper; ++i)
758 {
759 /* upper triangle in the supernode */
760 Uval[lastu] = SNptr[i];
761 /* Matlab doesn't like explicit zero. */
762 if (Uval[lastu] != 0.0)
763 Urow[lastu++] = L_SUB(istart+i);
764 }
765 Ucol[j+1] = lastu;
766
767 /* Extract L */
768 Lval[lastl] = 1.0; /* unit diagonal */
769 Lrow[lastl++] = L_SUB(istart + upper - 1);
770 for (int i = upper; i < nsupr; ++i)
771 {
772 Lval[lastl] = SNptr[i];
773 /* Matlab doesn't like explicit zero. */
774 if (Lval[lastl] != 0.0)
775 Lrow[lastl++] = L_SUB(istart+i);
776 }
777 Lcol[j+1] = lastl;
778
779 ++upper;
780 } /* for j ... */
781
782 } /* for k ... */
783
784 // squeeze the matrices :
785 m_l.resizeNonZeros(lastl);
786 m_u.resizeNonZeros(lastu);
787
789 }
790}
void resizeNonZeros(Index size)
Definition SparseMatrix.h:644
void resize(Index rows, Index cols)
Definition SparseMatrix.h:621
bool m_extractedDataAreDirty
Definition SuperLUSupport.h:464
MatrixType::Scalar Scalar
Definition SuperLUSupport.h:325
LUMatrixType m_l
Definition SuperLUSupport.h:445
LUMatrixType m_u
Definition SuperLUSupport.h:446

◆ factorize()

template<typename MatrixType >
void Eigen::SuperLU< MatrixType >::factorize ( const MatrixType matrix)

Performs a numeric decomposition of matrix

The given matrix must has the same sparcity than the matrix on which the symbolic decomposition has been performed.

See also
analyzePattern()
614{
615 eigen_assert(m_analysisIsOk && "You must first call analyzePattern()");
616 if(!m_analysisIsOk)
617 {
619 return;
620 }
621
622 this->initFactorization(a);
623
624 m_sluOptions.ColPerm = COLAMD;
625 int info = 0;
626 RealScalar recip_pivot_growth, rcond;
627 RealScalar ferr, berr;
628
629 StatInit(&m_sluStat);
630 SuperLU_gssvx(&m_sluOptions, &m_sluA, m_q.data(), m_p.data(), &m_sluEtree[0],
632 &m_sluL, &m_sluU,
633 NULL, 0,
634 &m_sluB, &m_sluX,
635 &recip_pivot_growth, &rcond,
636 &ferr, &berr,
637 &m_sluStat, &info, Scalar());
638 StatFree(&m_sluStat);
639
641
642 // FIXME how to better check for errors ???
644 m_factorizationIsOk = true;
645}
void initFactorization(const MatrixType &a)
Definition SuperLUSupport.h:391
int m_analysisIsOk
Definition SuperLUSupport.h:463

References eigen_assert, Eigen::InvalidInput, Eigen::NumericalIssue, and Eigen::Success.

◆ info()

ComputationInfo Eigen::SuperLUBase< _MatrixType, SuperLU< _MatrixType > >::info ( ) const
inlineinherited

Reports whether previous computation was successful.

Returns
Success if computation was succesful, NumericalIssue if the matrix.appears to be negative.
359 {
360 eigen_assert(m_isInitialized && "Decomposition is not initialized.");
361 return m_info;
362 }
ComputationInfo m_info
Definition SuperLUSupport.h:461

◆ init()

template<typename _MatrixType >
void Eigen::SuperLU< _MatrixType >::init ( )
inlineprotected
597 {
598 Base::init();
599
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;
605 }
void init()
Definition SuperLUSupport.h:420

References Eigen::SuperLUBase< _MatrixType, Derived >::init(), and Eigen::SuperLU< _MatrixType >::m_sluOptions.

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

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

◆ initFactorization()

void Eigen::SuperLUBase< _MatrixType, SuperLU< _MatrixType > >::initFactorization ( const MatrixType a)
inlineprotectedinherited
392 {
393 set_default_options(&this->m_sluOptions);
394
395 const Index size = a.rows();
396 m_matrix = a;
397
398 m_sluA = internal::asSluMatrix(m_matrix);
399 clearFactors();
400
401 m_p.resize(size);
402 m_q.resize(size);
403 m_sluRscale.resize(size);
404 m_sluCscale.resize(size);
405 m_sluEtree.resize(size);
406
407 // set empty B and X
408 m_sluB.setStorageType(SLU_DN);
410 m_sluB.Mtype = SLU_GE;
411 m_sluB.storage.values = 0;
412 m_sluB.nrow = 0;
413 m_sluB.ncol = 0;
414 m_sluB.storage.lda = internal::convert_index<int>(size);
415 m_sluX = m_sluB;
416
418 }
SluMatrix m_sluX
Definition SuperLUSupport.h:453
SluMatrix m_sluB
Definition SuperLUSupport.h:453
Matrix< RealScalar, Dynamic, 1 > m_sluRscale
Definition SuperLUSupport.h:457
std::vector< int > m_sluEtree
Definition SuperLUSupport.h:456
IntRowVectorType m_q
Definition SuperLUSupport.h:448
void clearFactors()
Definition SuperLUSupport.h:430
Matrix< RealScalar, Dynamic, 1 > m_sluCscale
Definition SuperLUSupport.h:457
SluMatrix m_sluA
Definition SuperLUSupport.h:451
IntColVectorType m_p
Definition SuperLUSupport.h:447
struct Eigen::SluMatrix::@636 storage
void setScalarType()
Definition SuperLUSupport.h:158
void setStorageType(Stype_t t)
Definition SuperLUSupport.h:145

◆ matrixL()

template<typename _MatrixType >
const LMatrixType & Eigen::SuperLU< _MatrixType >::matrixL ( ) const
inline
544 {
546 return m_l;
547 }
LUMatrixType m_l
Definition SuperLUSupport.h:445

References Eigen::SuperLUBase< _MatrixType, SuperLU< _MatrixType > >::extractData(), Eigen::SuperLU< _MatrixType >::m_extractedDataAreDirty, and Eigen::SuperLU< _MatrixType >::m_l.

+ Here is the call graph for this function:

◆ matrixU()

template<typename _MatrixType >
const UMatrixType & Eigen::SuperLU< _MatrixType >::matrixU ( ) const
inline
550 {
552 return m_u;
553 }

References Eigen::SuperLUBase< _MatrixType, SuperLU< _MatrixType > >::extractData(), Eigen::SuperLU< _MatrixType >::m_extractedDataAreDirty, and Eigen::SuperLU< _MatrixType >::m_u.

+ Here is the call graph for this function:

◆ options()

superlu_options_t & Eigen::SuperLUBase< _MatrixType, SuperLU< _MatrixType > >::options ( )
inlineinherited
Returns
a reference to the Super LU option object to configure the Super LU algorithms.
351{ return m_sluOptions; }
superlu_options_t m_sluOptions
Definition SuperLUSupport.h:455

◆ permutationP()

template<typename _MatrixType >
const IntColVectorType & Eigen::SuperLU< _MatrixType >::permutationP ( ) const
inline
556 {
558 return m_p;
559 }

References Eigen::SuperLUBase< _MatrixType, SuperLU< _MatrixType > >::extractData(), Eigen::SuperLU< _MatrixType >::m_extractedDataAreDirty, and Eigen::SuperLU< _MatrixType >::m_p.

+ Here is the call graph for this function:

◆ permutationQ()

template<typename _MatrixType >
const IntRowVectorType & Eigen::SuperLU< _MatrixType >::permutationQ ( ) const
inline
562 {
564 return m_q;
565 }

References Eigen::SuperLUBase< _MatrixType, SuperLU< _MatrixType > >::extractData(), Eigen::SuperLU< _MatrixType >::m_extractedDataAreDirty, and Eigen::SuperLU< _MatrixType >::m_q.

+ Here is the call graph for this function:

◆ rows()

Index Eigen::SuperLUBase< _MatrixType, SuperLU< _MatrixType > >::rows ( ) const
inlineinherited
347{ return m_matrix.rows(); }

◆ 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 Eigen::SuperLUBase< _MatrixType, Derived >::m_analysisIsOk
protected

◆ m_extractedDataAreDirty

template<typename _MatrixType >
bool Eigen::SuperLUBase< _MatrixType, Derived >::m_extractedDataAreDirty
mutableprotected

◆ m_factorizationIsOk

template<typename _MatrixType >
int Eigen::SuperLUBase< _MatrixType, Derived >::m_factorizationIsOk
protected

◆ m_info

template<typename _MatrixType >
ComputationInfo Eigen::SuperLUBase< _MatrixType, Derived >::m_info
mutableprotected

◆ m_isInitialized

template<typename _MatrixType >
bool Eigen::SparseSolverBase< Derived >::m_isInitialized
mutableprotected

◆ m_l

template<typename _MatrixType >
LUMatrixType Eigen::SuperLUBase< _MatrixType, Derived >::m_l
mutableprotected

◆ m_matrix

template<typename _MatrixType >
LUMatrixType Eigen::SuperLUBase< _MatrixType, Derived >::m_matrix
mutableprotected

◆ m_p

template<typename _MatrixType >
IntColVectorType Eigen::SuperLUBase< _MatrixType, Derived >::m_p
mutableprotected

◆ m_q

template<typename _MatrixType >
IntRowVectorType Eigen::SuperLUBase< _MatrixType, Derived >::m_q
mutableprotected

◆ m_sluA

template<typename _MatrixType >
SluMatrix Eigen::SuperLUBase< _MatrixType, Derived >::m_sluA
mutableprotected

◆ m_sluB

template<typename _MatrixType >
SluMatrix Eigen::SuperLUBase< _MatrixType, Derived >::m_sluB
mutableprotected

◆ m_sluBerr

template<typename _MatrixType >
Matrix<RealScalar,Dynamic,1> Eigen::SuperLUBase< _MatrixType, Derived >::m_sluBerr
protected

◆ m_sluCscale

template<typename _MatrixType >
Matrix<RealScalar,Dynamic,1> Eigen::SuperLUBase< _MatrixType, Derived >::m_sluCscale
protected

◆ m_sluEqued

template<typename _MatrixType >
char Eigen::SuperLUBase< _MatrixType, Derived >::m_sluEqued
mutableprotected

◆ m_sluEtree

template<typename _MatrixType >
std::vector<int> Eigen::SuperLUBase< _MatrixType, Derived >::m_sluEtree
mutableprotected

◆ m_sluFerr

template<typename _MatrixType >
Matrix<RealScalar,Dynamic,1> Eigen::SuperLUBase< _MatrixType, Derived >::m_sluFerr
mutableprotected

◆ m_sluL

template<typename _MatrixType >
SuperMatrix Eigen::SuperLUBase< _MatrixType, Derived >::m_sluL
mutableprotected

◆ m_sluOptions

template<typename _MatrixType >
superlu_options_t Eigen::SuperLUBase< _MatrixType, Derived >::m_sluOptions
mutableprotected

◆ m_sluRscale

template<typename _MatrixType >
Matrix<RealScalar,Dynamic,1> Eigen::SuperLUBase< _MatrixType, Derived >::m_sluRscale
mutableprotected

◆ m_sluStat

template<typename _MatrixType >
SuperLUStat_t Eigen::SuperLUBase< _MatrixType, Derived >::m_sluStat
mutableprotected

◆ m_sluU

template<typename _MatrixType >
SuperMatrix Eigen::SuperLUBase< _MatrixType, Derived >::m_sluU
protected

◆ m_sluX

template<typename _MatrixType >
SluMatrix Eigen::SuperLUBase< _MatrixType, Derived >::m_sluX
protected

◆ m_u

template<typename _MatrixType >
LUMatrixType Eigen::SuperLUBase< _MatrixType, Derived >::m_u
mutableprotected

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