Eigen::SparseInverse< Scalar > Class Template Reference

calculate sparse subset of inverse of sparse matrix More...

Public Types
typedef SparseMatrix< Scalar, ColMajor >	MatrixType
typedef SparseMatrix< Scalar, RowMajor >	RowMatrixType

Public Member Functions
SparseInverse &	compute (const SparseMatrix< Scalar > &A)
	Calculate the sparse inverse from a given sparse input. More...
const MatrixType &	inverse () const
	return the already-calculated sparse inverse, or a 0x0 matrix if it could not be computed More...
	SparseInverse ()
	SparseInverse (const SparseLU< MatrixType > &slu)
	This Constructor is for if you already have a factored SparseLU and would like to use it to calculate a sparse inverse. More...

Static Public Member Functions
static MatrixType	computeInverse (const RowMatrixType &Upper, const Matrix< Scalar, Dynamic, 1 > &inverseDiagonal, const MatrixType &Lower)
	Internal function to calculate the inverse from strictly upper, diagonal and strictly lower components. More...
static MatrixType	computeInverse (const SparseLU< MatrixType > &slu)
	Internal function to calculate the sparse inverse in a functional way. More...

Private Attributes
MatrixType	_result

Detailed Description

template<typename Scalar>
class Eigen::SparseInverse< Scalar >

calculate sparse subset of inverse of sparse matrix

This class returns a sparse subset of the inverse of the input matrix. The nonzeros correspond to the nonzeros of the input, plus any additional elements required due to fill-in of the internal LU factorization. This is is minimized via a applying a fill-reducing permutation as part of the LU factorization.

If there are specific entries of the input matrix which you need inverse values for, which are zero for the input, you need to insert entries into the input sparse matrix for them to be calculated.

Due to the sensitive nature of matrix inversion, particularly on large matrices which are made possible via sparsity, high accuracy dot products based on Kahan summation are used to reduce numerical error. If you still encounter numerical errors you may with to equilibrate your matrix before calculating the inverse, as well as making sure it is actually full rank.

Definition at line 127 of file SparseInverse.h.

Member Typedef Documentation

MatrixType

template<typename Scalar >

typedef SparseMatrix<Scalar, ColMajor> Eigen::SparseInverse< Scalar >::MatrixType

Definition at line 129 of file SparseInverse.h.

RowMatrixType

template<typename Scalar >

typedef SparseMatrix<Scalar, RowMajor> Eigen::SparseInverse< Scalar >::RowMatrixType

Definition at line 130 of file SparseInverse.h.

Constructor & Destructor Documentation

SparseInverse() [1/2]

template<typename Scalar >

Eigen::SparseInverse< Scalar >::SparseInverse ( )

inline

Definition at line 132 of file SparseInverse.h.

{}

SparseInverse() [2/2]

template<typename Scalar >

Eigen::SparseInverse< Scalar >::SparseInverse ( const SparseLU< MatrixType > &  slu )

inline

This Constructor is for if you already have a factored SparseLU and would like to use it to calculate a sparse inverse.

Just call this constructor with your already factored SparseLU class and you can directly call the .inverse() method to get the result.

Definition at line 141 of file SparseInverse.h.

{ _result = computeInverse(slu); }

Member Function Documentation

compute()

template<typename Scalar >

SparseInverse& Eigen::SparseInverse< Scalar >::compute ( const SparseMatrix< Scalar > &  A )

inline

Calculate the sparse inverse from a given sparse input.

Definition at line 146 of file SparseInverse.h.

                                                        {
    SparseLU<MatrixType> slu;
    slu.compute(A);
    _result = computeInverse(slu);
    return *this;
  }

computeInverse() [1/2]

template<typename Scalar >

static MatrixType Eigen::SparseInverse< Scalar >::computeInverse ( const RowMatrixType &  Upper, const Matrix< Scalar, Dynamic, 1 > &  inverseDiagonal, const MatrixType &  Lower  )

inlinestatic

Internal function to calculate the inverse from strictly upper, diagonal and strictly lower components.

Definition at line 184 of file SparseInverse.h.

                                                            {
    // Calculate the 'minimal set', which is the nonzeros of (L+U).transpose()
    // It could be zeroed, but we will overwrite all non-zeros anyways.
    MatrixType colInv = Lower.transpose().template triangularView<UnitUpper>();
    colInv += Upper.transpose();
 
    // We also need rowmajor representation in order to do efficient row-wise dot products
    RowMatrixType rowInv = Upper.transpose().template triangularView<UnitLower>();
    rowInv += Lower.transpose();
 
    // Use the Takahashi algorithm to build the supporting elements of the inverse
    // upwards and to the left, from the bottom right element, 1 col/row at a time
    for (Index recurseLevel = Upper.cols() - 1; recurseLevel >= 0; recurseLevel--) {
      const auto& col = Lower.col(recurseLevel);
      const auto& row = Upper.row(recurseLevel);
 
      // Calculate the inverse values for the nonzeros in this column
      typename MatrixType::ReverseInnerIterator colIter(colInv, recurseLevel);
      for (; recurseLevel < colIter.index(); --colIter) {
        const Scalar element = -accurateDot(col, rowInv.row(colIter.index()));
        colIter.valueRef() = element;
        rowInv.coeffRef(colIter.index(), recurseLevel) = element;
      }
 
      // Calculate the inverse values for the nonzeros in this row
      typename RowMatrixType::ReverseInnerIterator rowIter(rowInv, recurseLevel);
      for (; recurseLevel < rowIter.index(); --rowIter) {
        const Scalar element = -accurateDot(row, colInv.col(rowIter.index()));
        rowIter.valueRef() = element;
        colInv.coeffRef(recurseLevel, rowIter.index()) = element;
      }
 
      // And finally the diagonal, which corresponds to both row and col iterator now
      const Scalar diag = inverseDiagonal(recurseLevel) - accurateDot(row, colInv.col(recurseLevel));
      rowIter.valueRef() = diag;
      colIter.valueRef() = diag;
    }
 
    return colInv;
  }

computeInverse() [2/2]

template<typename Scalar >

static MatrixType Eigen::SparseInverse< Scalar >::computeInverse ( const SparseLU< MatrixType > &  slu )

inlinestatic

Internal function to calculate the sparse inverse in a functional way.

Returns

A sparse inverse representation, or, if the decomposition didn't complete, a 0x0 matrix.

Definition at line 162 of file SparseInverse.h.

                                                                    {
    if (slu.info() != Success) {
      return MatrixType(0, 0);
    }
 
    // Extract from SparseLU and decompose into L, inverse D and U terms
    Matrix<Scalar, Dynamic, 1> invD;
    RowMatrixType Upper;
    {
      RowMatrixType DU = slu.matrixU().toSparse();
      invD = DU.diagonal().cwiseInverse();
      Upper = (invD.asDiagonal() * DU).template triangularView<StrictlyUpper>();
    }
    MatrixType Lower = slu.matrixL().toSparse().template triangularView<StrictlyLower>();
 
    // Compute the inverse and reapply the permutation matrix from the LU decomposition
    return slu.colsPermutation().transpose() * computeInverse(Upper, invD, Lower) * slu.rowsPermutation();
  }

inverse()

template<typename Scalar >

const MatrixType& Eigen::SparseInverse< Scalar >::inverse ( ) const

inline

return the already-calculated sparse inverse, or a 0x0 matrix if it could not be computed

Definition at line 156 of file SparseInverse.h.

{ return _result; }

Member Data Documentation

_result

template<typename Scalar >

MatrixType Eigen::SparseInverse< Scalar >::_result

private

Definition at line 227 of file SparseInverse.h.

---

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

• SparseInverse.h