#include <Eigen/src/SparseLU/SparseLU.h>

Detailed Description

template<typename MatrixType_, typename OrderingType_>
class Eigen::SparseLU< MatrixType_, OrderingType_ >

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 arithmetic 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> A;
SparseLU<SparseMatrix<double>, COLAMDOrdering<int> > 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);

We can directly call compute() instead of analyzePattern() and factorize()

VectorXd x(n), b(n);
SparseMatrix<double> A;
SparseLU<SparseMatrix<double>, COLAMDOrdering<int> > solver;
// Fill A and b.
solver.compute(A);
// Use the factors to solve the linear system.
x = solver.solve(b);

Or give the matrix to the constructor SparseLU(const MatrixType& matrix)

VectorXd x(n), b(n);
SparseMatrix<double> A;
// Fill A and b.
SparseLU<SparseMatrix<double>, COLAMDOrdering<int> > solver(A);
// Use the factors to solve the linear system.
x = solver.solve(b);

Warning: The input matrix A should be in a compressed and column-major form. Otherwise an expensive copy will be made. You can call the inexpensive makeCompressed() to get a compressed matrix.

Note: Unlike the initial SuperLU implementation, there is no step to equilibrate the matrix. For badly scaled matrices, this step can be useful to reduce the pivoting during factorization. If this is the case for your matrices, you can try the basic scaling method at "unsupported/Eigen/src/IterativeSolvers/Scaling.h"

Template Parameters

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 .

See also: Sparse solver concept; OrderingMethods module

Inheritance diagram for Eigen::SparseLU< MatrixType_, OrderingType_ >:

Public Member Functions
Scalar	absDeterminant ()
	Give the absolute value of the determinant. More...

const SparseLUTransposeView< true, SparseLU< MatrixType_, OrderingType_ > >	adjoint ()
	Return a solver for the adjointed matrix. More...

void	analyzePattern (const MatrixType &matrix)
	Compute the column permutation. More...

Index	cols () const
	Give the number of columns.

const PermutationType &	colsPermutation () const
	Give the column matrix permutation. More...

void	compute (const MatrixType &matrix)
	Analyze and factorize the matrix so the solver is ready to solve. More...

Scalar	determinant ()
	Give the determinant. More...

void	factorize (const MatrixType &matrix)
	Factorize the matrix to get the solver ready. More...

ComputationInfo	info () const
	Reports whether previous computation was successful. More...

void	isSymmetric (bool sym)
	Let you set that the pattern of the input matrix is symmetric.

std::string	lastErrorMessage () const
	Give a human readable error. More...

Scalar	logAbsDeterminant () const
	Give the natural log of the absolute determinant. More...

SparseLUMatrixLReturnType< SCMatrix >	matrixL () const
	Give the matrixL. More...

SparseLUMatrixUReturnType< SCMatrix, Map< SparseMatrix< Scalar, ColMajor, StorageIndex > > >	matrixU () const
	Give the MatrixU. More...

Index	nnzL () const
	Give the number of non zero in matrix L.

Index	nnzU () const
	Give the number of non zero in matrix U.

Index	rows () const
	Give the number of rows.

const PermutationType &	rowsPermutation () const
	Give the row matrix permutation. More...

void	setPivotThreshold (const RealScalar &thresh)

Scalar	signDeterminant ()
	Give the sign of the determinant. More...

template<typename Rhs >
const Solve< SparseLU, Rhs >	solve (const MatrixBase< Rhs > &B) const
	Solve a system \( A X = B \). More...

	SparseLU ()
	Basic constructor of the solver. More...

	SparseLU (const MatrixType &matrix)
	Constructor of the solver already based on a specific matrix. More...

const SparseLUTransposeView< false, SparseLU< MatrixType_, OrderingType_ > >	transpose ()
	Return a solver for the transposed matrix. More...

Public Member Functions inherited from Eigen::SparseSolverBase< SparseLU< MatrixType_, OrderingType_ > >
const Solve< SparseLU< MatrixType_, OrderingType_ >, Rhs >	solve (const MatrixBase< Rhs > &b) const

const Solve< SparseLU< MatrixType_, OrderingType_ >, Rhs >	solve (const SparseMatrixBase< Rhs > &b) const

	SparseSolverBase ()

Additional Inherited Members
Protected Member Functions inherited from Eigen::internal::SparseLUImpl< MatrixType_::Scalar, MatrixType_::StorageIndex >
Index	column_bmod (const Index jcol, const Index nseg, BlockScalarVector dense, ScalarVector &tempv, BlockIndexVector segrep, BlockIndexVector repfnz, Index fpanelc, GlobalLU_t &glu)
	Performs numeric block updates (sup-col) in topological order. More...

Index	column_dfs (const Index m, const Index jcol, IndexVector &perm_r, Index maxsuper, Index &nseg, BlockIndexVector lsub_col, IndexVector &segrep, BlockIndexVector repfnz, IndexVector &xprune, IndexVector &marker, IndexVector &parent, IndexVector &xplore, GlobalLU_t &glu)
	Performs a symbolic factorization on column jcol and decide the supernode boundary. More...

Index	copy_to_ucol (const Index jcol, const Index nseg, IndexVector &segrep, BlockIndexVector repfnz, IndexVector &perm_r, BlockScalarVector dense, GlobalLU_t &glu)
	Performs numeric block updates (sup-col) in topological order. More...

void	countnz (const Index n, Index &nnzL, Index &nnzU, GlobalLU_t &glu)
	Count Nonzero elements in the factors.

Index	expand (VectorType &vec, Index &length, Index nbElts, Index keep_prev, Index &num_expansions)

void	fixupL (const Index n, const IndexVector &perm_r, GlobalLU_t &glu)
	Fix up the data storage lsub for L-subscripts. More...

void	heap_relax_snode (const Index n, IndexVector &et, const Index relax_columns, IndexVector &descendants, IndexVector &relax_end)
	Identify the initial relaxed supernodes. More...

Index	memInit (Index m, Index n, Index annz, Index lwork, Index fillratio, Index panel_size, GlobalLU_t &glu)
	Allocate various working space for the numerical factorization phase. More...

Index	memXpand (VectorType &vec, Index &maxlen, Index nbElts, MemType memtype, Index &num_expansions)
	Expand the existing storage. More...

void	panel_bmod (const Index m, const Index w, const Index jcol, const Index nseg, ScalarVector &dense, ScalarVector &tempv, IndexVector &segrep, IndexVector &repfnz, GlobalLU_t &glu)
	Performs numeric block updates (sup-panel) in topological order. More...

void	panel_dfs (const Index m, const Index w, const Index jcol, MatrixType &A, IndexVector &perm_r, Index &nseg, ScalarVector &dense, IndexVector &panel_lsub, IndexVector &segrep, IndexVector &repfnz, IndexVector &xprune, IndexVector &marker, IndexVector &parent, IndexVector &xplore, GlobalLU_t &glu)
	Performs a symbolic factorization on a panel of columns [jcol, jcol+w) More...

Index	pivotL (const Index jcol, const RealScalar &diagpivotthresh, IndexVector &perm_r, IndexVector &iperm_c, Index &pivrow, GlobalLU_t &glu)
	Performs the numerical pivoting on the current column of L, and the CDIV operation. More...

void	pruneL (const Index jcol, const IndexVector &perm_r, const Index pivrow, const Index nseg, const IndexVector &segrep, BlockIndexVector repfnz, IndexVector &xprune, GlobalLU_t &glu)
	Prunes the L-structure. More...

void	relax_snode (const Index n, IndexVector &et, const Index relax_columns, IndexVector &descendants, IndexVector &relax_end)
	Identify the initial relaxed supernodes. More...

Constructor & Destructor Documentation

◆ SparseLU() [1/2]

template<typename MatrixType_ , typename OrderingType_ >

Eigen::SparseLU< MatrixType_, OrderingType_ >::SparseLU ( )

inline

Basic constructor of the solver.

Construct a SparseLU. As no matrix is given as argument, compute() should be called afterward with a matrix.

◆ SparseLU() [2/2]

template<typename MatrixType_ , typename OrderingType_ >

Eigen::SparseLU< MatrixType_, OrderingType_ >::SparseLU ( const MatrixType & matrix )

inlineexplicit

Constructor of the solver already based on a specific matrix.

Construct a SparseLU. compute() is already called with the given matrix.

Member Function Documentation

◆ absDeterminant()

template<typename MatrixType_ , typename OrderingType_ >

Scalar Eigen::SparseLU< MatrixType_, OrderingType_ >::absDeterminant ( )

inline

Give the absolute value of the determinant.

Returns: the absolute value of the determinant of the matrix of which *this is the QR decomposition.

Warning: a determinant can be very big or small, so for matrices of large enough dimension, there is a risk of overflow/underflow. One way to work around that is to use logAbsDeterminant() instead.

See also: logAbsDeterminant(), signDeterminant()

◆ adjoint()

template<typename MatrixType_ , typename OrderingType_ >

const SparseLUTransposeView<true, SparseLU<MatrixType_, OrderingType_> > Eigen::SparseLU< MatrixType_, OrderingType_ >::adjoint ( )

inline

Return a solver for the adjointed matrix.

Returns: an expression of the adjoint of the factored matrix

A typical usage is to solve for the adjoint problem A' x = b:

solver.compute(A);

x = solver.adjoint().solve(b);

For real scalar types, this function is equivalent to transpose().

See also: transpose(), solve()

◆ analyzePattern()

template<typename MatrixType , typename OrderingType >

void Eigen::SparseLU< MatrixType, OrderingType >::analyzePattern ( const MatrixType & mat )

Compute the column permutation.

Compute the column permutation to minimize the fill-in

Apply this permutation to the input matrix -
Compute the column elimination tree on the permuted matrix
Postorder the elimination tree and the column permutation

It is possible to call compute() instead of analyzePattern() + factorize().

If the matrix is row-major this function will do an heavy copy.

See also: factorize(), compute()

◆ colsPermutation()

template<typename MatrixType_ , typename OrderingType_ >

const PermutationType& Eigen::SparseLU< MatrixType_, OrderingType_ >::colsPermutation ( ) const

inline

Give the column matrix permutation.

Returns: a reference to the column matrix permutation \( P_c^T \) such that \(P_r A P_c^T = L U\)

See also: rowsPermutation()

◆ compute()

template<typename MatrixType_ , typename OrderingType_ >

void Eigen::SparseLU< MatrixType_, OrderingType_ >::compute ( const MatrixType & matrix )

inline

Analyze and factorize the matrix so the solver is ready to solve.

Compute the symbolic and numeric factorization of the input sparse matrix. The input matrix should be in column-major storage, otherwise analyzePattern() will do a heavy copy.

Call analyzePattern() followed by factorize()

See also: analyzePattern(), factorize()

◆ determinant()

template<typename MatrixType_ , typename OrderingType_ >

Scalar Eigen::SparseLU< MatrixType_, OrderingType_ >::determinant ( )

inline

Give the determinant.

Returns: The determinant of the matrix.

See also: absDeterminant(), logAbsDeterminant()

◆ factorize()

template<typename MatrixType , typename OrderingType >

void Eigen::SparseLU< MatrixType, OrderingType >::factorize ( const MatrixType & matrix )

Factorize the matrix to get the solver ready.

Numerical factorization
Interleaved with the symbolic factorization

To get error of this function you should check info(), you can get more info of errors with lastErrorMessage().

In the past (before 2012 (git history is not older)), this function was returning an integer. This exit was 0 if successful factorization.

0 if info = i, and i is been completed, but the factor U is exactly singular,

and division by zero will occur if it is used to solve a system of equation.

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.

It seems that A was the name of the matrix in the past.

See also: analyzePattern(), compute(), SparseLU(), info(), lastErrorMessage()

◆ info()

template<typename MatrixType_ , typename OrderingType_ >

ComputationInfo Eigen::SparseLU< MatrixType_, OrderingType_ >::info ( ) const

inline

Reports whether previous computation was successful.

Returns: Success if computation was successful, NumericalIssue if the LU factorization reports a problem, zero diagonal for instance InvalidInput if the input matrix is invalid

You can get a readable error message with lastErrorMessage().

See also: lastErrorMessage()

◆ lastErrorMessage()

template<typename MatrixType_ , typename OrderingType_ >

std::string Eigen::SparseLU< MatrixType_, OrderingType_ >::lastErrorMessage ( ) const

inline

Give a human readable error.

Returns: A string describing the type of error

◆ logAbsDeterminant()

template<typename MatrixType_ , typename OrderingType_ >

Scalar Eigen::SparseLU< MatrixType_, OrderingType_ >::logAbsDeterminant ( ) const

inline

Give the natural log of the absolute determinant.

Returns: the natural log of the absolute value of the determinant of the matrix of which **this is the QR decomposition

Note: This method is useful to work around the risk of overflow/underflow that's inherent to the determinant computation.

See also: absDeterminant(), signDeterminant()

◆ matrixL()

template<typename MatrixType_ , typename OrderingType_ >

SparseLUMatrixLReturnType<SCMatrix> Eigen::SparseLU< MatrixType_, OrderingType_ >::matrixL ( ) const

inline

Give the matrixL.

Returns: an expression of the matrix L, internally stored as supernodes The only operation available with this expression is the triangular solve
y = b; matrixL().solveInPlace(y);

◆ matrixU()

template<typename MatrixType_ , typename OrderingType_ >

SparseLUMatrixUReturnType<SCMatrix, Map<SparseMatrix<Scalar, ColMajor, StorageIndex> > > Eigen::SparseLU< MatrixType_, OrderingType_ >::matrixU ( ) const

inline

Give the MatrixU.

Returns: an expression of the matrix U, The only operation available with this expression is the triangular solve
y = b; matrixU().solveInPlace(y);

◆ rowsPermutation()

template<typename MatrixType_ , typename OrderingType_ >

const PermutationType& Eigen::SparseLU< MatrixType_, OrderingType_ >::rowsPermutation ( ) const

inline

Give the row matrix permutation.

Returns: a reference to the row matrix permutation \( P_r \) such that \(P_r A P_c^T = L U\)

See also: colsPermutation()

◆ setPivotThreshold()

template<typename MatrixType_ , typename OrderingType_ >

void Eigen::SparseLU< MatrixType_, OrderingType_ >::setPivotThreshold ( const RealScalar & thresh )

inline

Set the threshold used for a diagonal entry to be an acceptable pivot.

◆ signDeterminant()

template<typename MatrixType_ , typename OrderingType_ >

Scalar Eigen::SparseLU< MatrixType_, OrderingType_ >::signDeterminant ( )

inline

Give the sign of the determinant.

Returns: A number representing the sign of the determinant

See also: absDeterminant(), logAbsDeterminant()

◆ solve()

template<typename MatrixType_ , typename OrderingType_ >

template<typename Rhs >

const Solve<SparseLU, Rhs> Eigen::SparseLU< MatrixType_, OrderingType_ >::solve ( const MatrixBase< Rhs > & B ) const

inline

Solve a system \( A X = B \).

Returns: the solution X of \( A X = B \) using the current decomposition of A.

Warning: the destination matrix X in X = this->solve(B) must be colmun-major.

See also: compute()

◆ transpose()

template<typename MatrixType_ , typename OrderingType_ >

const SparseLUTransposeView<false, SparseLU<MatrixType_, OrderingType_> > Eigen::SparseLU< MatrixType_, OrderingType_ >::transpose ( )

inline

Return a solver for the transposed matrix.

Returns: an expression of the transposed of the factored matrix.

A typical usage is to solve for the transposed problem A^T x = b:

solver.compute(A);

x = solver.transpose().solve(b);

See also: adjoint(), solve()

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

/builddir/build/BUILD/eigen-5.0.1/Eigen/src/SparseLU/SparseLU.h

	Eigen 5.0.1-dev

Detailed Description

template<typename MatrixType_, typename OrderingType_> class Eigen::SparseLU< MatrixType_, OrderingType_ >

Public Member Functions

Additional Inherited Members

Constructor & Destructor Documentation

◆ SparseLU() [1/2]

◆ SparseLU() [2/2]

Member Function Documentation

◆ absDeterminant()

◆ adjoint()

◆ analyzePattern()

◆ colsPermutation()

◆ compute()

◆ determinant()

◆ factorize()

◆ info()

◆ lastErrorMessage()

◆ logAbsDeterminant()

◆ matrixL()

◆ matrixU()

◆ rowsPermutation()

◆ setPivotThreshold()

◆ signDeterminant()

◆ solve()

◆ transpose()

template<typename MatrixType_, typename OrderingType_>
class Eigen::SparseLU< MatrixType_, OrderingType_ >