Eigen::PartialPivLU< MatrixType_, PermutationIndex_ > Class Template Reference

LU decomposition of a matrix with partial pivoting, and related features. More...

#include <PartialPivLU.h>

+ Inheritance diagram for Eigen::PartialPivLU< MatrixType_, PermutationIndex_ >:

Public Types

enum  { MaxRowsAtCompileTime = MatrixType::MaxRowsAtCompileTime , MaxColsAtCompileTime = MatrixType::MaxColsAtCompileTime }
 
typedef MatrixType_ MatrixType
 
typedef SolverBase< PartialPivLUBase
 
using PermutationIndex = PermutationIndex_
 
typedef PermutationMatrix< RowsAtCompileTime, MaxRowsAtCompileTime, PermutationIndexPermutationType
 
typedef Transpositions< RowsAtCompileTime, MaxRowsAtCompileTime, PermutationIndexTranspositionType
 
typedef MatrixType::PlainObject PlainObject
 
- Public Types inherited from Eigen::SolverBase< PartialPivLU< MatrixType_, PermutationIndex_ > >
enum  
 
typedef EigenBase< PartialPivLU< MatrixType_, PermutationIndex_ > > Base
 
typedef internal::traits< PartialPivLU< MatrixType_, PermutationIndex_ > >::Scalar Scalar
 
typedef Scalar CoeffReturnType
 
typedef Transpose< const PartialPivLU< MatrixType_, PermutationIndex_ > > ConstTransposeReturnType
 
typedef std::conditional_t< NumTraits< Scalar >::IsComplex, CwiseUnaryOp< internal::scalar_conjugate_op< Scalar >, const ConstTransposeReturnType >, const ConstTransposeReturnTypeAdjointReturnType
 
- Public Types inherited from Eigen::EigenBase< Derived >
typedef Eigen::Index Index
 The interface type of indices. More...
 
typedef internal::traits< Derived >::StorageKind StorageKind
 

Public Member Functions

 PartialPivLU ()
 Default Constructor. More...
 
 PartialPivLU (Index size)
 Default Constructor with memory preallocation. More...
 
template<typename InputType >
 PartialPivLU (const EigenBase< InputType > &matrix)
 
template<typename InputType >
 PartialPivLU (EigenBase< InputType > &matrix)
 
template<typename InputType >
PartialPivLUcompute (const EigenBase< InputType > &matrix)
 
const MatrixTypematrixLU () const
 
const PermutationTypepermutationP () const
 
RealScalar rcond () const
 
const Inverse< PartialPivLUinverse () const
 
Scalar determinant () const
 
MatrixType reconstructedMatrix () const
 
EIGEN_CONSTEXPR Index rows () const EIGEN_NOEXCEPT
 
EIGEN_CONSTEXPR Index cols () const EIGEN_NOEXCEPT
 
template<typename RhsType , typename DstType >
EIGEN_DEVICE_FUNC void _solve_impl (const RhsType &rhs, DstType &dst) const
 
template<bool Conjugate, typename RhsType , typename DstType >
EIGEN_DEVICE_FUNC void _solve_impl_transposed (const RhsType &rhs, DstType &dst) const
 
- Public Member Functions inherited from Eigen::SolverBase< PartialPivLU< MatrixType_, PermutationIndex_ > >
 SolverBase ()
 
 ~SolverBase ()
 
const Solve< PartialPivLU< MatrixType_, PermutationIndex_ >, Rhs > solve (const MatrixBase< Rhs > &b) const
 
const ConstTransposeReturnType transpose () const
 
const AdjointReturnType adjoint () const
 
constexpr EIGEN_DEVICE_FUNC PartialPivLU< MatrixType_, PermutationIndex_ > & derived ()
 
constexpr EIGEN_DEVICE_FUNC const PartialPivLU< MatrixType_, PermutationIndex_ > & derived () const
 
- Public Member Functions inherited from Eigen::EigenBase< Derived >
constexpr EIGEN_DEVICE_FUNC Derived & derived ()
 
constexpr EIGEN_DEVICE_FUNC const Derived & derived () const
 
EIGEN_DEVICE_FUNC Derived & const_cast_derived () const
 
EIGEN_DEVICE_FUNC const Derived & const_derived () const
 
EIGEN_DEVICE_FUNC EIGEN_CONSTEXPR Index rows () const EIGEN_NOEXCEPT
 
EIGEN_DEVICE_FUNC EIGEN_CONSTEXPR Index cols () const EIGEN_NOEXCEPT
 
EIGEN_DEVICE_FUNC EIGEN_CONSTEXPR Index size () const EIGEN_NOEXCEPT
 
template<typename Dest >
EIGEN_DEVICE_FUNC void evalTo (Dest &dst) const
 
template<typename Dest >
EIGEN_DEVICE_FUNC void addTo (Dest &dst) const
 
template<typename Dest >
EIGEN_DEVICE_FUNC void subTo (Dest &dst) const
 
template<typename Dest >
EIGEN_DEVICE_FUNC void applyThisOnTheRight (Dest &dst) const
 
template<typename Dest >
EIGEN_DEVICE_FUNC void applyThisOnTheLeft (Dest &dst) const
 
template<typename Device >
EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE DeviceWrapper< Derived, Device > device (Device &device)
 
template<typename Device >
EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE DeviceWrapper< const Derived, Device > device (Device &device) const
 

Protected Member Functions

void compute ()
 
- Protected Member Functions inherited from Eigen::SolverBase< PartialPivLU< MatrixType_, PermutationIndex_ > >
void _check_solve_assertion (const Rhs &b) const
 

Protected Attributes

MatrixType m_lu
 
PermutationType m_p
 
TranspositionType m_rowsTranspositions
 
RealScalar m_l1_norm
 
signed char m_det_p
 
bool m_isInitialized
 

Friends

class SolverBase< PartialPivLU >
 

Detailed Description

template<typename MatrixType_, typename PermutationIndex_>
class Eigen::PartialPivLU< MatrixType_, PermutationIndex_ >

LU decomposition of a matrix with partial pivoting, and related features.

Template Parameters
MatrixType_the type of the matrix of which we are computing the LU decomposition

This class represents a LU decomposition of a square invertible matrix, with partial pivoting: the matrix A is decomposed as A = PLU where L is unit-lower-triangular, U is upper-triangular, and P is a permutation matrix.

Typically, partial pivoting LU decomposition is only considered numerically stable for square invertible matrices. Thus LAPACK's dgesv and dgesvx require the matrix to be square and invertible. The present class does the same. It will assert that the matrix is square, but it won't (actually it can't) check that the matrix is invertible: it is your task to check that you only use this decomposition on invertible matrices.

The guaranteed safe alternative, working for all matrices, is the full pivoting LU decomposition, provided by class FullPivLU.

This is not a rank-revealing LU decomposition. Many features are intentionally absent from this class, such as rank computation. If you need these features, use class FullPivLU.

This LU decomposition is suitable to invert invertible matrices. It is what MatrixBase::inverse() uses in the general case. On the other hand, it is not suitable to determine whether a given matrix is invertible.

The data of the LU decomposition can be directly accessed through the methods matrixLU(), permutationP().

This class supports the inplace decomposition mechanism.

See also
MatrixBase::partialPivLu(), MatrixBase::determinant(), MatrixBase::inverse(), MatrixBase::computeInverse(), class FullPivLU

Member Typedef Documentation

◆ Base

template<typename MatrixType_ , typename PermutationIndex_ >
typedef SolverBase<PartialPivLU> Eigen::PartialPivLU< MatrixType_, PermutationIndex_ >::Base

◆ MatrixType

template<typename MatrixType_ , typename PermutationIndex_ >
typedef MatrixType_ Eigen::PartialPivLU< MatrixType_, PermutationIndex_ >::MatrixType

◆ PermutationIndex

template<typename MatrixType_ , typename PermutationIndex_ >
using Eigen::PartialPivLU< MatrixType_, PermutationIndex_ >::PermutationIndex = PermutationIndex_

◆ PermutationType

template<typename MatrixType_ , typename PermutationIndex_ >
typedef PermutationMatrix<RowsAtCompileTime, MaxRowsAtCompileTime, PermutationIndex> Eigen::PartialPivLU< MatrixType_, PermutationIndex_ >::PermutationType

◆ PlainObject

template<typename MatrixType_ , typename PermutationIndex_ >
typedef MatrixType::PlainObject Eigen::PartialPivLU< MatrixType_, PermutationIndex_ >::PlainObject

◆ TranspositionType

template<typename MatrixType_ , typename PermutationIndex_ >
typedef Transpositions<RowsAtCompileTime, MaxRowsAtCompileTime, PermutationIndex> Eigen::PartialPivLU< MatrixType_, PermutationIndex_ >::TranspositionType

Member Enumeration Documentation

◆ anonymous enum

template<typename MatrixType_ , typename PermutationIndex_ >
anonymous enum
Enumerator
MaxRowsAtCompileTime 
MaxColsAtCompileTime 
84  {
85  MaxRowsAtCompileTime = MatrixType::MaxRowsAtCompileTime,
86  MaxColsAtCompileTime = MatrixType::MaxColsAtCompileTime
87  };
@ MaxRowsAtCompileTime
Definition: PartialPivLU.h:85
@ MaxColsAtCompileTime
Definition: PartialPivLU.h:86

Constructor & Destructor Documentation

◆ PartialPivLU() [1/4]

template<typename MatrixType , typename PermutationIndex >
Eigen::PartialPivLU< MatrixType, PermutationIndex >::PartialPivLU

Default Constructor.

The default constructor is useful in cases in which the user intends to perform decompositions via PartialPivLU::compute(const MatrixType&).

271  : m_lu(), m_p(), m_rowsTranspositions(), m_l1_norm(0), m_det_p(0), m_isInitialized(false) {}
bool m_isInitialized
Definition: PartialPivLU.h:266
RealScalar m_l1_norm
Definition: PartialPivLU.h:264
TranspositionType m_rowsTranspositions
Definition: PartialPivLU.h:263
PermutationType m_p
Definition: PartialPivLU.h:262
signed char m_det_p
Definition: PartialPivLU.h:265
MatrixType m_lu
Definition: PartialPivLU.h:261

◆ PartialPivLU() [2/4]

template<typename MatrixType , typename PermutationIndex >
Eigen::PartialPivLU< MatrixType, PermutationIndex >::PartialPivLU ( Index  size)
explicit

Default Constructor with memory preallocation.

Like the default constructor but with preallocation of the internal data according to the specified problem size.

See also
PartialPivLU()
EIGEN_DEVICE_FUNC EIGEN_CONSTEXPR Index size() const EIGEN_NOEXCEPT
Definition: EigenBase.h:64

◆ PartialPivLU() [3/4]

template<typename MatrixType , typename PermutationIndex >
template<typename InputType >
Eigen::PartialPivLU< MatrixType, PermutationIndex >::PartialPivLU ( const EigenBase< InputType > &  matrix)
explicit

Constructor.

Parameters
matrixthe matrix of which to compute the LU decomposition.
Warning
The matrix should have full rank (e.g. if it's square, it should be invertible). If you need to deal with non-full rank, use class FullPivLU instead.
280  : m_lu(matrix.rows(), matrix.cols()),
281  m_p(matrix.rows()),
283  m_l1_norm(0),
284  m_det_p(0),
285  m_isInitialized(false) {
286  compute(matrix.derived());
287 }
void compute()
Definition: PartialPivLU.h:481
Eigen::Map< Eigen::Matrix< T, Eigen::Dynamic, Eigen::Dynamic, Eigen::ColMajor >, 0, Eigen::OuterStride<> > matrix(T *data, int rows, int cols, int stride)
Definition: common.h:85

References Eigen::PartialPivLU< MatrixType_, PermutationIndex_ >::compute(), and matrix().

◆ PartialPivLU() [4/4]

template<typename MatrixType , typename PermutationIndex >
template<typename InputType >
Eigen::PartialPivLU< MatrixType, PermutationIndex >::PartialPivLU ( EigenBase< InputType > &  matrix)
explicit

Constructor for inplace decomposition

Parameters
matrixthe matrix of which to compute the LU decomposition.
Warning
The matrix should have full rank (e.g. if it's square, it should be invertible). If you need to deal with non-full rank, use class FullPivLU instead.
292  : m_lu(matrix.derived()),
293  m_p(matrix.rows()),
295  m_l1_norm(0),
296  m_det_p(0),
297  m_isInitialized(false) {
298  compute();
299 }

References Eigen::PartialPivLU< MatrixType_, PermutationIndex_ >::compute().

Member Function Documentation

◆ _solve_impl()

template<typename MatrixType_ , typename PermutationIndex_ >
template<typename RhsType , typename DstType >
EIGEN_DEVICE_FUNC void Eigen::PartialPivLU< MatrixType_, PermutationIndex_ >::_solve_impl ( const RhsType &  rhs,
DstType &  dst 
) const
inline
218  {
219  /* The decomposition PA = LU can be rewritten as A = P^{-1} L U.
220  * So we proceed as follows:
221  * Step 1: compute c = Pb.
222  * Step 2: replace c by the solution x to Lx = c.
223  * Step 3: replace c by the solution x to Ux = c.
224  */
225 
226  // Step 1
227  dst = permutationP() * rhs;
228 
229  // Step 2
230  m_lu.template triangularView<UnitLower>().solveInPlace(dst);
231 
232  // Step 3
233  m_lu.template triangularView<Upper>().solveInPlace(dst);
234  }
const PermutationType & permutationP() const
Definition: PartialPivLU.h:149

References Eigen::PartialPivLU< MatrixType_, PermutationIndex_ >::m_lu, and Eigen::PartialPivLU< MatrixType_, PermutationIndex_ >::permutationP().

◆ _solve_impl_transposed()

template<typename MatrixType_ , typename PermutationIndex_ >
template<bool Conjugate, typename RhsType , typename DstType >
EIGEN_DEVICE_FUNC void Eigen::PartialPivLU< MatrixType_, PermutationIndex_ >::_solve_impl_transposed ( const RhsType &  rhs,
DstType &  dst 
) const
inline
237  {
238  /* The decomposition PA = LU can be rewritten as A^T = U^T L^T P.
239  * So we proceed as follows:
240  * Step 1: compute c as the solution to L^T c = b
241  * Step 2: replace c by the solution x to U^T x = c.
242  * Step 3: update c = P^-1 c.
243  */
244 
245  eigen_assert(rhs.rows() == m_lu.cols());
246 
247  // Step 1
248  dst = m_lu.template triangularView<Upper>().transpose().template conjugateIf<Conjugate>().solve(rhs);
249  // Step 2
250  m_lu.template triangularView<UnitLower>().transpose().template conjugateIf<Conjugate>().solveInPlace(dst);
251  // Step 3
252  dst = permutationP().transpose() * dst;
253  }
#define eigen_assert(x)
Definition: Macros.h:910
InverseReturnType transpose() const
Definition: PermutationMatrix.h:177

References eigen_assert, Eigen::PartialPivLU< MatrixType_, PermutationIndex_ >::m_lu, Eigen::PartialPivLU< MatrixType_, PermutationIndex_ >::permutationP(), and Eigen::PermutationBase< Derived >::transpose().

◆ cols()

◆ compute() [1/2]

template<typename MatrixType , typename PermutationIndex >
void Eigen::PartialPivLU< MatrixType, PermutationIndex >::compute
protected
481  {
482  eigen_assert(m_lu.rows() < NumTraits<PermutationIndex>::highest());
483 
484  if (m_lu.cols() > 0)
485  m_l1_norm = m_lu.cwiseAbs().colwise().sum().maxCoeff();
486  else
487  m_l1_norm = RealScalar(0);
488 
489  eigen_assert(m_lu.rows() == m_lu.cols() && "PartialPivLU is only for square (and moreover invertible) matrices");
490  const Index size = m_lu.rows();
491 
493 
494  typename TranspositionType::StorageIndex nb_transpositions;
496  m_det_p = (nb_transpositions % 2) ? -1 : 1;
497 
499 
500  m_isInitialized = true;
501 }
NumTraits< Scalar >::Real RealScalar
Definition: bench_gemm.cpp:46
void resize(Index newSize)
Definition: Transpositions.h:63
IndicesType::Scalar StorageIndex
Definition: Transpositions.h:147
void partial_lu_inplace(MatrixType &lu, TranspositionType &row_transpositions, typename TranspositionType::StorageIndex &nb_transpositions)
Definition: PartialPivLU.h:460
Eigen::Index Index
The interface type of indices.
Definition: EigenBase.h:43

References eigen_assert, Eigen::internal::partial_lu_inplace(), and size.

Referenced by Eigen::PartialPivLU< MatrixType_, PermutationIndex_ >::compute(), and Eigen::PartialPivLU< MatrixType_, PermutationIndex_ >::PartialPivLU().

◆ compute() [2/2]

template<typename MatrixType_ , typename PermutationIndex_ >
template<typename InputType >
PartialPivLU& Eigen::PartialPivLU< MatrixType_, PermutationIndex_ >::compute ( const EigenBase< InputType > &  matrix)
inline

◆ determinant()

template<typename MatrixType , typename PermutationIndex >
PartialPivLU< MatrixType, PermutationIndex >::Scalar Eigen::PartialPivLU< MatrixType, PermutationIndex >::determinant
Returns
the determinant of the matrix of which *this is the LU decomposition. It has only linear complexity (that is, O(n) where n is the dimension of the square matrix) as the LU decomposition has already been computed.
Note
For fixed-size matrices of size up to 4, MatrixBase::determinant() offers optimized paths.
Warning
a determinant can be very big or small, so for matrices of large enough dimension, there is a risk of overflow/underflow.
See also
MatrixBase::determinant()
505  {
506  eigen_assert(m_isInitialized && "PartialPivLU is not initialized.");
507  return Scalar(m_det_p) * m_lu.diagonal().prod();
508 }
internal::traits< PartialPivLU< MatrixType_, PermutationIndex_ > >::Scalar Scalar
Definition: SolverBase.h:75

References eigen_assert.

Referenced by lu_verify_assert().

◆ inverse()

template<typename MatrixType_ , typename PermutationIndex_ >
const Inverse<PartialPivLU> Eigen::PartialPivLU< MatrixType_, PermutationIndex_ >::inverse ( ) const
inline
Returns
the inverse of the matrix of which *this is the LU decomposition.
Warning
The matrix being decomposed here is assumed to be invertible. If you need to check for invertibility, use class FullPivLU instead.
See also
MatrixBase::inverse(), LU::inverse()
191  {
192  eigen_assert(m_isInitialized && "PartialPivLU is not initialized.");
193  return Inverse<PartialPivLU>(*this);
194  }

References eigen_assert, and Eigen::PartialPivLU< MatrixType_, PermutationIndex_ >::m_isInitialized.

Referenced by lu_partial_piv(), and lu_verify_assert().

◆ matrixLU()

template<typename MatrixType_ , typename PermutationIndex_ >
const MatrixType& Eigen::PartialPivLU< MatrixType_, PermutationIndex_ >::matrixLU ( ) const
inline
Returns
the LU decomposition matrix: the upper-triangular part is U, the unit-lower-triangular part is L (at least for square matrices; in the non-square case, special care is needed, see the documentation of class FullPivLU).
See also
matrixL(), matrixU()
142  {
143  eigen_assert(m_isInitialized && "PartialPivLU is not initialized.");
144  return m_lu;
145  }

References eigen_assert, Eigen::PartialPivLU< MatrixType_, PermutationIndex_ >::m_isInitialized, and Eigen::PartialPivLU< MatrixType_, PermutationIndex_ >::m_lu.

Referenced by lu_verify_assert().

◆ permutationP()

template<typename MatrixType_ , typename PermutationIndex_ >
const PermutationType& Eigen::PartialPivLU< MatrixType_, PermutationIndex_ >::permutationP ( ) const
inline

◆ rcond()

template<typename MatrixType_ , typename PermutationIndex_ >
RealScalar Eigen::PartialPivLU< MatrixType_, PermutationIndex_ >::rcond ( ) const
inline
Returns
an estimate of the reciprocal condition number of the matrix of which *this is the LU decomposition.
179  {
180  eigen_assert(m_isInitialized && "PartialPivLU is not initialized.");
182  }
Decomposition::RealScalar rcond_estimate_helper(typename Decomposition::RealScalar matrix_norm, const Decomposition &dec)
Reciprocal condition number estimator.
Definition: ConditionEstimator.h:157

References eigen_assert, Eigen::PartialPivLU< MatrixType_, PermutationIndex_ >::m_isInitialized, Eigen::PartialPivLU< MatrixType_, PermutationIndex_ >::m_l1_norm, and Eigen::internal::rcond_estimate_helper().

Referenced by lu_partial_piv().

◆ reconstructedMatrix()

template<typename MatrixType , typename PermutationIndex >
MatrixType Eigen::PartialPivLU< MatrixType, PermutationIndex >::reconstructedMatrix
Returns
the matrix represented by the decomposition, i.e., it returns the product: P^{-1} L U. This function is provided for debug purpose.
514  {
515  eigen_assert(m_isInitialized && "LU is not initialized.");
516  // LU
517  MatrixType res = m_lu.template triangularView<UnitLower>().toDenseMatrix() * m_lu.template triangularView<Upper>();
518 
519  // P^{-1}(LU)
520  res = m_p.inverse() * res;
521 
522  return res;
523 }
cout<< "Here is the matrix m:"<< endl<< m<< endl;Matrix< ptrdiff_t, 3, 1 > res
Definition: PartialRedux_count.cpp:3
MatrixXf MatrixType
Definition: benchmark-blocking-sizes.cpp:52
InverseReturnType inverse() const
Definition: PermutationMatrix.h:172

References eigen_assert, and res.

Referenced by lu_partial_piv().

◆ rows()

Friends And Related Function Documentation

◆ SolverBase< PartialPivLU >

template<typename MatrixType_ , typename PermutationIndex_ >
friend class SolverBase< PartialPivLU >
friend

Member Data Documentation

◆ m_det_p

template<typename MatrixType_ , typename PermutationIndex_ >
signed char Eigen::PartialPivLU< MatrixType_, PermutationIndex_ >::m_det_p
protected

◆ m_isInitialized

◆ m_l1_norm

template<typename MatrixType_ , typename PermutationIndex_ >
RealScalar Eigen::PartialPivLU< MatrixType_, PermutationIndex_ >::m_l1_norm
protected

◆ m_lu

◆ m_p

template<typename MatrixType_ , typename PermutationIndex_ >
PermutationType Eigen::PartialPivLU< MatrixType_, PermutationIndex_ >::m_p
protected

◆ m_rowsTranspositions

template<typename MatrixType_ , typename PermutationIndex_ >
TranspositionType Eigen::PartialPivLU< MatrixType_, PermutationIndex_ >::m_rowsTranspositions
protected

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