template<typename _MatrixType>
Eigen::SPQR class

Sparse QR factorization based on SuiteSparseQR library.

Template parameters
_MatrixType The type of the sparse matrix A, must be a column-major SparseMatrix<>

Contents

This class is used to perform a multithreaded and multifrontal rank-revealing QR decomposition of sparse matrices. The result is then used to solve linear leasts_square systems. Clearly, a QR factorization is returned such that A*P = Q*R where :

P is the column permutation. Use colsPermutation() to get it.

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

R is the sparse triangular factor. Use matrixQR() to get it as SparseMatrix. NOTE : The Index type of R is always SuiteSparse_long. You can get it with SPQR::Index

This class follows the sparse solver concept.

Base classes

template<typename Derived>
class SparseSolverBase
A base class for sparse solvers.

Public functions

auto cholmodCommon() const -> cholmod_common*
auto cols() const -> Index
auto colsPermutation() const -> PermutationType
Get the permutation that was applied to columns of A.
auto info() const -> ComputationInfo
Reports whether previous computation was successful.
auto matrixQ() const -> SPQRMatrixQReturnType<SPQR>
Get an expression of the matrix Q.
auto matrixR() const -> const MatrixType
auto rank() const -> Index
auto rows() const -> Index
void setPivotThreshold(const RealScalar& tol)
Set the tolerance tol to treat columns with 2-norm < =tol as zero.
void setSPQROrdering(int ord)
Set the fill-reducing ordering method to be used.

Function documentation

template<typename _MatrixType>
cholmod_common* Eigen::SPQR<_MatrixType>::cholmodCommon() const

Returns a pointer to the SPQR workspace

template<typename _MatrixType>
Index Eigen::SPQR<_MatrixType>::cols() const

Get the number of columns of the input matrix.

template<typename _MatrixType>
ComputationInfo Eigen::SPQR<_MatrixType>::info() const

Reports whether previous computation was successful.

Returns Success if computation was successful, NumericalIssue if the sparse QR can not be computed

template<typename _MatrixType>
const MatrixType Eigen::SPQR<_MatrixType>::matrixR() const

Returns the sparse triangular factor R. It is a sparse matrix

template<typename _MatrixType>
Index Eigen::SPQR<_MatrixType>::rank() const

Gets the rank of the matrix. It should be equal to matrixQR().cols if the matrix is full-rank

template<typename _MatrixType>
Index Eigen::SPQR<_MatrixType>::rows() const

Get the number of rows of the input matrix and the Q matrix