Eigen  3.4.90 (git rev a4098ac676528a83cfb73d4d26ce1b42ec05f47c)
Eigen::SparseQR< MatrixType_, OrderingType_ > Class Template Reference

Detailed Description

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

Sparse left-looking QR factorization with numerical column pivoting.

This class implements a left-looking QR decomposition of sparse matrices with numerical column pivoting. When a column has a norm less than a given tolerance it is implicitly permuted to the end. The QR factorization thus obtained is given by A*P = Q*R where R is upper triangular or trapezoidal.

P is the column permutation which is the product of the fill-reducing and the numerical permutations. Use colsPermutation() to get it.

Q is the orthogonal matrix represented as products of Householder reflectors. Use matrixQ() to get an expression and matrixQ().adjoint() to get the adjoint. You can then apply it to a vector.

R is the sparse triangular or trapezoidal matrix. The later occurs when A is rank-deficient. matrixR().topLeftCorner(rank(), rank()) always returns a triangular factor of full rank.

Template Parameters
MatrixType_The type of the sparse matrix A, must be a column-major SparseMatrix<>
OrderingType_The fill-reducing ordering method. See the OrderingMethods module for the list of built-in and external ordering methods.

This class follows the sparse solver concept .

The numerical pivoting strategy and default threshold are the same as in SuiteSparse QR, and detailed in the following paper: Tim Davis, "Algorithm 915, SuiteSparseQR: Multifrontal Multithreaded Rank-Revealing Sparse QR Factorization, ACM Trans. on Math. Soft. 38(1), 2011. </i> Even though it is qualified as "rank-revealing", this strategy might fail for some rank deficient problems. When this class is used to solve linear or least-square problems it is thus strongly recommended to check the accuracy of the computed solution. If it failed, it usually helps to increase the threshold with setPivotThreshold.

Warning
The input sparse matrix A must be in compressed mode (see SparseMatrix::makeCompressed()).
For complex matrices matrixQ().transpose() will actually return the adjoint matrix.

+ Inheritance diagram for Eigen::SparseQR< MatrixType_, OrderingType_ >:

Public Member Functions

void analyzePattern (const MatrixType &mat)
 Preprocessing step of a QR factorization. More...
 
Index cols () const
 
const PermutationTypecolsPermutation () const
 
void compute (const MatrixType &mat)
 
void factorize (const MatrixType &mat)
 Performs the numerical QR factorization of the input matrix. More...
 
ComputationInfo info () const
 Reports whether previous computation was successful. More...
 
std::string lastErrorMessage () const
 
SparseQRMatrixQReturnType< SparseQRmatrixQ () const
 
const QRMatrixTypematrixR () const
 
Index rank () const
 
Index rows () const
 
void setPivotThreshold (const RealScalar &threshold)
 
template<typename Rhs >
const Solve< SparseQR, Rhs > solve (const MatrixBase< Rhs > &B) const
 
 SparseQR (const MatrixType &mat)
 
- Public Member Functions inherited from Eigen::SparseSolverBase< SparseQR< MatrixType_, OrderingType_ > >
const Solve< SparseQR< MatrixType_, OrderingType_ >, Rhs > solve (const MatrixBase< Rhs > &b) const
 
const Solve< SparseQR< MatrixType_, OrderingType_ >, Rhs > solve (const SparseMatrixBase< Rhs > &b) const
 
 SparseSolverBase ()
 

Constructor & Destructor Documentation

◆ SparseQR()

template<typename MatrixType_ , typename OrderingType_ >
Eigen::SparseQR< MatrixType_, OrderingType_ >::SparseQR ( const MatrixType &  mat)
inlineexplicit

Construct a QR factorization of the matrix mat.

Warning
The matrix mat must be in compressed mode (see SparseMatrix::makeCompressed()).
See also
compute()

Member Function Documentation

◆ analyzePattern()

template<typename MatrixType , typename OrderingType >
void Eigen::SparseQR< MatrixType, OrderingType >::analyzePattern ( const MatrixType &  mat)

Preprocessing step of a QR factorization.

Warning
The matrix mat must be in compressed mode (see SparseMatrix::makeCompressed()).

In this step, the fill-reducing permutation is computed and applied to the columns of A and the column elimination tree is computed as well. Only the sparsity pattern of mat is exploited.

Note
In this step it is assumed that there is no empty row in the matrix mat.

◆ cols()

template<typename MatrixType_ , typename OrderingType_ >
Index Eigen::SparseQR< MatrixType_, OrderingType_ >::cols ( void  ) const
inline
Returns
the number of columns of the represented matrix.

◆ colsPermutation()

template<typename MatrixType_ , typename OrderingType_ >
const PermutationType & Eigen::SparseQR< MatrixType_, OrderingType_ >::colsPermutation ( ) const
inline
Returns
a const reference to the column permutation P that was applied to A such that A*P = Q*R It is the combination of the fill-in reducing permutation and numerical column pivoting.

◆ compute()

template<typename MatrixType_ , typename OrderingType_ >
void Eigen::SparseQR< MatrixType_, OrderingType_ >::compute ( const MatrixType &  mat)
inline

Computes the QR factorization of the sparse matrix mat.

Warning
The matrix mat must be in compressed mode (see SparseMatrix::makeCompressed()).
See also
analyzePattern(), factorize()

◆ factorize()

template<typename MatrixType , typename OrderingType >
void Eigen::SparseQR< MatrixType, OrderingType >::factorize ( const MatrixType &  mat)

Performs the numerical QR factorization of the input matrix.

The function SparseQR::analyzePattern(const MatrixType&) must have been called beforehand with a matrix having the same sparsity pattern than mat.

Parameters
matThe sparse column-major matrix

◆ info()

template<typename MatrixType_ , typename OrderingType_ >
ComputationInfo Eigen::SparseQR< MatrixType_, OrderingType_ >::info ( ) const
inline

Reports whether previous computation was successful.

Returns
Success if computation was successful, NumericalIssue if the QR factorization reports a numerical problem InvalidInput if the input matrix is invalid
See also
iparm()

◆ lastErrorMessage()

template<typename MatrixType_ , typename OrderingType_ >
std::string Eigen::SparseQR< MatrixType_, OrderingType_ >::lastErrorMessage ( ) const
inline
Returns
A string describing the type of error. This method is provided to ease debugging, not to handle errors.

◆ matrixQ()

template<typename MatrixType_ , typename OrderingType_ >
SparseQRMatrixQReturnType< SparseQR > Eigen::SparseQR< MatrixType_, OrderingType_ >::matrixQ ( void  ) const
inline
Returns
an expression of the matrix Q as products of sparse Householder reflectors. The common usage of this function is to apply it to a dense matrix or vector
VectorXd B1, B2;
// Initialize B1
B2 = matrixQ() * B1;
SparseQRMatrixQReturnType< SparseQR > matrixQ() const
Definition: SparseQR.h:188

To get a plain SparseMatrix representation of Q:

SparseMatrix<double> Q;
Q = SparseQR<SparseMatrix<double> >(A).matrixQ();

Internally, this call simply performs a sparse product between the matrix Q and a sparse identity matrix. However, due to the fact that the sparse reflectors are stored unsorted, two transpositions are needed to sort them before performing the product.

◆ matrixR()

template<typename MatrixType_ , typename OrderingType_ >
const QRMatrixType & Eigen::SparseQR< MatrixType_, OrderingType_ >::matrixR ( ) const
inline
Returns
a const reference to the sparse upper triangular matrix R of the QR factorization.
Warning
The entries of the returned matrix are not sorted. This means that using it in algorithms expecting sorted entries will fail. This include random coefficient accesses (SpaseMatrix::coeff()), and coefficient-wise operations. Matrix products and triangular solves are fine though.

To sort the entries, you can assign it to a row-major matrix, and if a column-major matrix is required, you can copy it again:

SparseMatrix<double> R = qr.matrixR(); // column-major, not sorted!
SparseMatrix<double,RowMajor> Rr = qr.matrixR(); // row-major, sorted
SparseMatrix<double> Rc = Rr; // column-major, sorted

◆ rank()

template<typename MatrixType_ , typename OrderingType_ >
Index Eigen::SparseQR< MatrixType_, OrderingType_ >::rank ( ) const
inline
Returns
the number of non linearly dependent columns as determined by the pivoting threshold.
See also
setPivotThreshold()

◆ rows()

template<typename MatrixType_ , typename OrderingType_ >
Index Eigen::SparseQR< MatrixType_, OrderingType_ >::rows ( void  ) const
inline
Returns
the number of rows of the represented matrix.

◆ setPivotThreshold()

template<typename MatrixType_ , typename OrderingType_ >
void Eigen::SparseQR< MatrixType_, OrderingType_ >::setPivotThreshold ( const RealScalar &  threshold)
inline

Sets the threshold that is used to determine linearly dependent columns during the factorization.

In practice, if during the factorization the norm of the column that has to be eliminated is below this threshold, then the entire column is treated as zero, and it is moved at the end.

◆ solve()

template<typename MatrixType_ , typename OrderingType_ >
template<typename Rhs >
const Solve< SparseQR, Rhs > Eigen::SparseQR< MatrixType_, OrderingType_ >::solve ( const MatrixBase< Rhs > &  B) const
inline
Returns
the solution X of \( A X = B \) using the current decomposition of A.
See also
compute()

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