![]() |
Eigen
3.3.3
|
Sparse supernodal LU factorization for general matrices.
This class implements the supernodal LU factorization for general matrices. It uses the main techniques from the sequential SuperLU package (http://crd-legacy.lbl.gov/~xiaoye/SuperLU/). It handles transparently real and complex arithmetics with single and double precision, depending on the scalar type of your input matrix. The code has been optimized to provide BLAS-3 operations during supernode-panel updates. It benefits directly from the built-in high-performant Eigen BLAS routines. Moreover, when the size of a supernode is very small, the BLAS calls are avoided to enable a better optimization from the compiler. For best performance, you should compile it with NDEBUG flag to avoid the numerous bounds checking on vectors.
An important parameter of this class is the ordering method. It is used to reorder the columns (and eventually the rows) of the matrix to reduce the number of new elements that are created during numerical factorization. The cheapest method available is COLAMD. See the OrderingMethods module for the list of built-in and external ordering methods.
Simple example with key steps
VectorXd x(n), b(n); SparseMatrix<double, ColMajor> A; SparseLU<SparseMatrix<scalar, ColMajor>, COLAMDOrdering<Index> > solver; // fill A and b; // Compute the ordering permutation vector from the structural pattern of A solver.analyzePattern(A); // Compute the numerical factorization solver.factorize(A); //Use the factors to solve the linear system x = solver.solve(b);
_MatrixType | The type of the sparse matrix. It must be a column-major SparseMatrix<> |
_OrderingType | The ordering method to use, either AMD, COLAMD or METIS. Default is COLMAD |
This class follows the sparse solver concept .
Public Member Functions | |
Scalar | absDeterminant () |
void | analyzePattern (const MatrixType &matrix) |
const PermutationType & | colsPermutation () const |
void | compute (const MatrixType &matrix) |
Scalar | determinant () |
void | factorize (const MatrixType &matrix) |
ComputationInfo | info () const |
Reports whether previous computation was successful. | |
void | isSymmetric (bool sym) |
std::string | lastErrorMessage () const |
Scalar | logAbsDeterminant () const |
SparseLUMatrixLReturnType < SCMatrix > | matrixL () const |
SparseLUMatrixUReturnType < SCMatrix, MappedSparseMatrix < Scalar, ColMajor, StorageIndex > > | matrixU () const |
const PermutationType & | rowsPermutation () const |
void | setPivotThreshold (const RealScalar &thresh) |
Scalar | signDeterminant () |
template<typename Rhs > | |
const Solve< SparseLU, Rhs > | solve (const MatrixBase< Rhs > &B) const |
Scalar Eigen::SparseLU< _MatrixType, _OrderingType >::absDeterminant | ( | ) | [inline] |
void Eigen::SparseLU< MatrixType, OrderingType >::analyzePattern | ( | const MatrixType & | mat | ) |
Compute the column permutation to minimize the fill-in
const PermutationType& Eigen::SparseLU< _MatrixType, _OrderingType >::colsPermutation | ( | ) | const [inline] |
void Eigen::SparseLU< _MatrixType, _OrderingType >::compute | ( | const MatrixType & | matrix | ) | [inline] |
Compute the symbolic and numeric factorization of the input sparse matrix. The input matrix should be in column-major storage.
Scalar Eigen::SparseLU< _MatrixType, _OrderingType >::determinant | ( | ) | [inline] |
void Eigen::SparseLU< MatrixType, OrderingType >::factorize | ( | const MatrixType & | matrix | ) |
= 0: successful factorization
> 0: if info = i, and i is
<= A->ncol: U(i,i) is exactly zero. The factorization has been completed, but the factor U is exactly singular, and division by zero will occur if it is used to solve a system of equations.
> A->ncol: number of bytes allocated when memory allocation failure occurred, plus A->ncol. If lwork = -1, it is the estimated amount of space needed, plus A->ncol.
ComputationInfo Eigen::SparseLU< _MatrixType, _OrderingType >::info | ( | ) | const [inline] |
Reports whether previous computation was successful.
Success
if computation was succesful, NumericalIssue
if the LU factorization reports a problem, zero diagonal for instance InvalidInput
if the input matrix is invalidvoid Eigen::SparseLU< _MatrixType, _OrderingType >::isSymmetric | ( | bool | sym | ) | [inline] |
Indicate that the pattern of the input matrix is symmetric
std::string Eigen::SparseLU< _MatrixType, _OrderingType >::lastErrorMessage | ( | ) | const [inline] |
Scalar Eigen::SparseLU< _MatrixType, _OrderingType >::logAbsDeterminant | ( | ) | const [inline] |
SparseLUMatrixLReturnType<SCMatrix> Eigen::SparseLU< _MatrixType, _OrderingType >::matrixL | ( | ) | const [inline] |
y = b; matrixL().solveInPlace(y);
SparseLUMatrixUReturnType<SCMatrix,MappedSparseMatrix<Scalar,ColMajor,StorageIndex> > Eigen::SparseLU< _MatrixType, _OrderingType >::matrixU | ( | ) | const [inline] |
y = b; matrixU().solveInPlace(y);
const PermutationType& Eigen::SparseLU< _MatrixType, _OrderingType >::rowsPermutation | ( | ) | const [inline] |
void Eigen::SparseLU< _MatrixType, _OrderingType >::setPivotThreshold | ( | const RealScalar & | thresh | ) | [inline] |
Set the threshold used for a diagonal entry to be an acceptable pivot.
Scalar Eigen::SparseLU< _MatrixType, _OrderingType >::signDeterminant | ( | ) | [inline] |
const Solve<SparseLU, Rhs> Eigen::SparseLU< _MatrixType, _OrderingType >::solve | ( | const MatrixBase< Rhs > & | B | ) | const [inline] |
Reimplemented from Eigen::SparseSolverBase< SparseLU< _MatrixType, _OrderingType > >.