Medial Code Documentation
Loading...
Searching...
No Matches
Public Types | Public Member Functions | Protected Types | Protected Attributes | Friends
Eigen::SPQR< _MatrixType > Class Template Reference

Sparse QR factorization based on SuiteSparseQR library. More...

#include <SuiteSparseQRSupport.h>

Inheritance diagram for Eigen::SPQR< _MatrixType >:
Eigen::SparseSolverBase< SPQR< _MatrixType > > Eigen::internal::noncopyable

Public Types

enum  { ColsAtCompileTime = Dynamic , MaxColsAtCompileTime = Dynamic }
 
typedef _MatrixType::Scalar Scalar
 
typedef _MatrixType::RealScalar RealScalar
 
typedef SuiteSparse_long StorageIndex
 
typedef SparseMatrix< Scalar, ColMajor, StorageIndexMatrixType
 
typedef Map< PermutationMatrix< Dynamic, Dynamic, StorageIndex > > PermutationType
 

Public Member Functions

 SPQR (const _MatrixType &matrix)
 
void SPQR_free ()
 
void compute (const _MatrixType &matrix)
 
Index rows () const
 Get the number of rows of the input matrix and the Q matrix.
 
Index cols () const
 Get the number of columns of the input matrix.
 
template<typename Rhs , typename Dest >
void _solve_impl (const MatrixBase< Rhs > &b, MatrixBase< Dest > &dest) const
 
const MatrixType matrixR () const
 
SPQRMatrixQReturnType< SPQRmatrixQ () const
 Get an expression of the matrix Q.
 
PermutationType colsPermutation () const
 Get the permutation that was applied to columns of A.
 
Index rank () const
 Gets the rank of the matrix.
 
void setSPQROrdering (int ord)
 Set the fill-reducing ordering method to be used.
 
void setPivotThreshold (const RealScalar &tol)
 Set the tolerance tol to treat columns with 2-norm < =tol as zero.
 
cholmod_commoncholmodCommon () const
 
ComputationInfo info () const
 Reports whether previous computation was successful.
 
- Public Member Functions inherited from Eigen::SparseSolverBase< SPQR< _MatrixType > >
 SparseSolverBase ()
 Default constructor.
 
SPQR< _MatrixType > & derived ()
 
const SPQR< _MatrixType > & derived () const
 
const Solve< SPQR< _MatrixType >, Rhs > solve (const MatrixBase< Rhs > &b) const
 
const Solve< SPQR< _MatrixType >, Rhs > solve (const SparseMatrixBase< Rhs > &b) const
 
void _solve_impl (const SparseMatrixBase< Rhs > &b, SparseMatrixBase< Dest > &dest) const
 

Protected Types

typedef SparseSolverBase< SPQR< _MatrixType > > Base
 

Protected Attributes

bool m_analysisIsOk
 
bool m_factorizationIsOk
 
bool m_isRUpToDate
 
ComputationInfo m_info
 
int m_ordering
 
int m_allow_tol
 
RealScalar m_tolerance
 
cholmod_sparsem_cR
 
MatrixType m_R
 
StorageIndexm_E
 
cholmod_sparsem_H
 
StorageIndexm_HPinv
 
cholmod_densem_HTau
 
Index m_rank
 
cholmod_common m_cc
 
bool m_useDefaultThreshold
 
- Protected Attributes inherited from Eigen::SparseSolverBase< SPQR< _MatrixType > >
bool m_isInitialized
 

Friends

template<typename , typename >
struct SPQR_QProduct
 

Detailed Description

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

Sparse QR factorization based on SuiteSparseQR library.

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

Template Parameters
_MatrixTypeThe type of the sparse matrix A, must be a column-major SparseMatrix<>

\implsparsesolverconcept

Member Function Documentation

◆ cholmodCommon()

template<typename _MatrixType >
cholmod_common * Eigen::SPQR< _MatrixType >::cholmodCommon ( ) const
inline
Returns
a pointer to the SPQR workspace

◆ info()

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

Reports whether previous computation was successful.

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

◆ matrixR()

template<typename _MatrixType >
const MatrixType Eigen::SPQR< _MatrixType >::matrixR ( ) const
inline
Returns
the sparse triangular factor R. It is a sparse matrix

◆ rank()

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

Gets the rank of the matrix.

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


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