g2o::LinearSolverEigen< MatrixType > Class Template Reference

linear solver which uses the sparse Cholesky solver from Eigen More...

#include <linear_solver_eigen.h>

Inheritance diagram for g2o::LinearSolverEigen< MatrixType >:

Collaboration diagram for g2o::LinearSolverEigen< MatrixType >:

Classes
class	CholeskyDecomposition
	Sub-classing Eigen's SimplicialLDLT to perform ordering with a given ordering. More...

Public Types
typedef Eigen::SparseMatrix< double, Eigen::ColMajor >	SparseMatrix
typedef Eigen::Triplet< double >	Triplet
typedef Eigen::PermutationMatrix< Eigen::Dynamic, Eigen::Dynamic >	PermutationMatrix

Public Member Functions
	LinearSolverEigen ()
virtual	~LinearSolverEigen ()
virtual bool	init ()
bool	solve (const SparseBlockMatrix< MatrixType > &A, double *x, double *b)
bool	blockOrdering () const
	do the AMD ordering on the blocks or on the scalar matrix More...
void	setBlockOrdering (bool blockOrdering)
virtual bool	writeDebug () const
	write a debug dump of the system matrix if it is not SPD in solve More...
virtual void	setWriteDebug (bool b)
Public Member Functions inherited from g2o::LinearSolver< MatrixType >
	LinearSolver ()
virtual	~LinearSolver ()
virtual bool	solveBlocks (double **&blocks, const SparseBlockMatrix< MatrixType > &A)
virtual bool	solvePattern (SparseBlockMatrix< MatrixXD > &spinv, const std::vector< std::pair< int, int > > &blockIndices, const SparseBlockMatrix< MatrixType > &A)

Protected Member Functions
void	computeSymbolicDecomposition (const SparseBlockMatrix< MatrixType > &A)
void	fillSparseMatrix (const SparseBlockMatrix< MatrixType > &A, bool onlyValues)

Protected Attributes
bool	_init
bool	_blockOrdering
bool	_writeDebug
SparseMatrix	_sparseMatrix
CholeskyDecomposition	_cholesky

Detailed Description

template<typename MatrixType>
class g2o::LinearSolverEigen< MatrixType >

linear solver which uses the sparse Cholesky solver from Eigen

Has no dependencies except Eigen. Hence, should compile almost everywhere without to much issues. Performance should be similar to CSparse, I guess.

Definition at line 49 of file linear_solver_eigen.h.

Member Typedef Documentation

template<typename MatrixType >

typedef Eigen::PermutationMatrix<Eigen::Dynamic, Eigen::Dynamic> g2o::LinearSolverEigen< MatrixType >::PermutationMatrix

Definition at line 54 of file linear_solver_eigen.h.

template<typename MatrixType >

typedef Eigen::SparseMatrix<double, Eigen::ColMajor> g2o::LinearSolverEigen< MatrixType >::SparseMatrix

Definition at line 52 of file linear_solver_eigen.h.

template<typename MatrixType >

typedef Eigen::Triplet<double> g2o::LinearSolverEigen< MatrixType >::Triplet

Definition at line 53 of file linear_solver_eigen.h.

Constructor & Destructor Documentation

template<typename MatrixType >

g2o::LinearSolverEigen< MatrixType >::LinearSolverEigen ( )

inline

Definition at line 76 of file linear_solver_eigen.h.

                        :
      LinearSolver<MatrixType>(),
      _init(true), _blockOrdering(false), _writeDebug(false)
    {
    }

template<typename MatrixType >

virtual g2o::LinearSolverEigen< MatrixType >::~LinearSolverEigen ( )

inlinevirtual

Definition at line 82 of file linear_solver_eigen.h.

    {
    }

Member Function Documentation

template<typename MatrixType >

bool g2o::LinearSolverEigen< MatrixType >::blockOrdering ( ) const

inline

do the AMD ordering on the blocks or on the scalar matrix

Definition at line 125 of file linear_solver_eigen.h.

References g2o::LinearSolverEigen< MatrixType >::_blockOrdering.

Referenced by g2o::LinearSolverEigen< MatrixType >::setBlockOrdering().

{ return _blockOrdering;}

template<typename MatrixType >

void g2o::LinearSolverEigen< MatrixType >::computeSymbolicDecomposition ( const SparseBlockMatrix< MatrixType > &  A )

inlineprotected

compute the symbolic decompostion of the matrix only once. Since A has the same pattern in all the iterations, we only compute the fill-in reducing ordering once and re-use for all the following iterations.

Definition at line 145 of file linear_solver_eigen.h.

References g2o::LinearSolverEigen< MatrixType >::CholeskyDecomposition::analyzePatternWithPermutation(), g2o::SparseBlockMatrix< MatrixType >::blockCols(), g2o::SparseBlockMatrix< MatrixType >::colBaseOfBlock(), g2o::SparseBlockMatrix< MatrixType >::cols(), g2o::SparseBlockMatrix< MatrixType >::colsOfBlock(), g2o::get_monotonic_time(), g2o::G2OBatchStatistics::globalStats(), g2o::SparseBlockMatrix< MatrixType >::rows(), and g2o::G2OBatchStatistics::timeSymbolicDecomposition.

Referenced by g2o::LinearSolverEigen< MatrixType >::solve().

    {
      double t=get_monotonic_time();
      if (! _blockOrdering) {
        _cholesky.analyzePattern(_sparseMatrix);
      } else {
        // block ordering with the Eigen Interface
        // This is really ugly currently, as it calls internal functions from Eigen
        // and modifies the SparseMatrix class
        Eigen::PermutationMatrix<Eigen::Dynamic,Eigen::Dynamic> blockP;
        {
          // prepare a block structure matrix for calling AMD
          std::vector<Triplet> triplets;
          for (size_t c = 0; c < A.blockCols().size(); ++c){
            const typename SparseBlockMatrix<MatrixType>::IntBlockMap& column = A.blockCols()[c];
            for (typename SparseBlockMatrix<MatrixType>::IntBlockMap::const_iterator it = column.begin(); it != column.end(); ++it) {
              const int& r = it->first;
              if (r > static_cast<int>(c)) // only upper triangle
                break;
              triplets.push_back(Triplet(r, c, 0.));
            }
          }

          // call the AMD ordering on the block matrix.
          // Relies on Eigen's internal stuff, probably bad idea
          SparseMatrix auxBlockMatrix(A.blockCols().size(), A.blockCols().size());
          auxBlockMatrix.setFromTriplets(triplets.begin(), triplets.end());
          typename CholeskyDecomposition::CholMatrixType C;
          C = auxBlockMatrix.selfadjointView<Eigen::Upper>();
          Eigen::internal::minimum_degree_ordering(C, blockP);
        }

        int rows = A.rows();
        assert(rows == A.cols() && "Matrix A is not square");

        // Adapt the block permutation to the scalar matrix
        PermutationMatrix scalarP;
        scalarP.resize(rows);
        int scalarIdx = 0;
        for (int i = 0; i < blockP.size(); ++i) {
          const int& p = blockP.indices()(i);
          int base  = A.colBaseOfBlock(p);
          int nCols = A.colsOfBlock(p);
          for (int j = 0; j < nCols; ++j)
            scalarP.indices()(scalarIdx++) = base++;
        }
        assert(scalarIdx == rows && "did not completely fill the permutation matrix");
        // analyze with the scalar permutation
        _cholesky.analyzePatternWithPermutation(_sparseMatrix, scalarP);

      }
      G2OBatchStatistics* globalStats = G2OBatchStatistics::globalStats();
      if (globalStats)
        globalStats->timeSymbolicDecomposition = get_monotonic_time() - t;
    }

template<typename MatrixType >

void g2o::LinearSolverEigen< MatrixType >::fillSparseMatrix ( const SparseBlockMatrix< MatrixType > &  A, bool  onlyValues  )

inlineprotected

Definition at line 201 of file linear_solver_eigen.h.

References g2o::SparseBlockMatrix< MatrixType >::blockCols(), g2o::SparseBlockMatrix< MatrixType >::colBaseOfBlock(), g2o::SparseBlockMatrix< MatrixType >::fillCCS(), g2o::SparseBlockMatrix< MatrixType >::nonZeros(), and g2o::SparseBlockMatrix< MatrixType >::rowBaseOfBlock().

Referenced by g2o::LinearSolverEigen< MatrixType >::solve().

    {
      if (onlyValues) {
        A.fillCCS(_sparseMatrix.valuePtr(), true);
      } else {

        // create from triplet structure
        std::vector<Triplet> triplets;
        triplets.reserve(A.nonZeros());
        for (size_t c = 0; c < A.blockCols().size(); ++c) {
          int colBaseOfBlock = A.colBaseOfBlock(c);
          const typename SparseBlockMatrix<MatrixType>::IntBlockMap& column = A.blockCols()[c];
          for (typename SparseBlockMatrix<MatrixType>::IntBlockMap::const_iterator it = column.begin(); it != column.end(); ++it) {
            int rowBaseOfBlock = A.rowBaseOfBlock(it->first);
            const MatrixType& m = *(it->second);
            for (int cc = 0; cc < m.cols(); ++cc) {
              int aux_c = colBaseOfBlock + cc;
              for (int rr = 0; rr < m.rows(); ++rr) {
                int aux_r = rowBaseOfBlock + rr;
                if (aux_r > aux_c)
                  break;
                triplets.push_back(Triplet(aux_r, aux_c, m(rr, cc)));
              }
            }
          }
        }
        _sparseMatrix.setFromTriplets(triplets.begin(), triplets.end());

      }
    }

template<typename MatrixType >

virtual bool g2o::LinearSolverEigen< MatrixType >::init ( )

inlinevirtual

init for operating on matrices with a different non-zero pattern like before

Implements g2o::LinearSolver< MatrixType >.

Definition at line 86 of file linear_solver_eigen.h.

References g2o::LinearSolverEigen< MatrixType >::_init.

    {
      _init = true;
      return true;
    }

template<typename MatrixType >

void g2o::LinearSolverEigen< MatrixType >::setBlockOrdering ( bool  blockOrdering )

inline

Definition at line 126 of file linear_solver_eigen.h.

References g2o::LinearSolverEigen< MatrixType >::_blockOrdering, and g2o::LinearSolverEigen< MatrixType >::blockOrdering().

{ _blockOrdering = blockOrdering;}

template<typename MatrixType >

virtual void g2o::LinearSolverEigen< MatrixType >::setWriteDebug ( bool  b )

inlinevirtual

Reimplemented from g2o::LinearSolver< MatrixType >.

Definition at line 130 of file linear_solver_eigen.h.

References g2o::LinearSolverEigen< MatrixType >::_writeDebug.

{ _writeDebug = b;}

template<typename MatrixType >

bool g2o::LinearSolverEigen< MatrixType >::solve ( const SparseBlockMatrix< MatrixType > &  A, double *  x, double *  b  )

inlinevirtual

Assumes that A is the same matrix for several calls. Among other assumptions, the non-zero pattern does not change! If the matrix changes call init() before. solve system Ax = b, x and b have to allocated beforehand!!

Implements g2o::LinearSolver< MatrixType >.

Definition at line 92 of file linear_solver_eigen.h.

References g2o::LinearSolverEigen< MatrixType >::_cholesky, g2o::LinearSolverEigen< MatrixType >::_init, g2o::LinearSolverEigen< MatrixType >::_sparseMatrix, g2o::LinearSolverEigen< MatrixType >::_writeDebug, g2o::G2OBatchStatistics::choleskyNNZ, g2o::SparseBlockMatrix< MatrixType >::cols(), g2o::LinearSolverEigen< MatrixType >::computeSymbolicDecomposition(), g2o::LinearSolverEigen< MatrixType >::fillSparseMatrix(), g2o::get_monotonic_time(), g2o::G2OBatchStatistics::globalStats(), g2o::SparseBlockMatrix< MatrixType >::rows(), g2o::G2OBatchStatistics::timeNumericDecomposition, and g2o::SparseBlockMatrix< MatrixType >::writeOctave().

    {
      if (_init)
        _sparseMatrix.resize(A.rows(), A.cols());
      fillSparseMatrix(A, !_init);
      if (_init) // compute the symbolic composition once
        computeSymbolicDecomposition(A);
      _init = false;

      double t=get_monotonic_time();
      _cholesky.factorize(_sparseMatrix);
      if (_cholesky.info() != Eigen::Success) { // the matrix is not positive definite
        if (_writeDebug) {
          std::cerr << "Cholesky failure, writing debug.txt (Hessian loadable by Octave)" << std::endl;
          A.writeOctave("debug.txt");
        }
        return false;
      }

      // Solving the system
      VectorXD::MapType xx(x, _sparseMatrix.cols());
      VectorXD::ConstMapType bb(b, _sparseMatrix.cols());
      xx = _cholesky.solve(bb);
      G2OBatchStatistics* globalStats = G2OBatchStatistics::globalStats();
      if (globalStats) {
        globalStats->timeNumericDecomposition = get_monotonic_time() - t;
        globalStats->choleskyNNZ = _cholesky.matrixL().nestedExpression().nonZeros() + _sparseMatrix.cols(); // the elements of D
      }

      return true;
    }

template<typename MatrixType >

virtual bool g2o::LinearSolverEigen< MatrixType >::writeDebug ( ) const

inlinevirtual

write a debug dump of the system matrix if it is not SPD in solve

Reimplemented from g2o::LinearSolver< MatrixType >.

Definition at line 129 of file linear_solver_eigen.h.

References g2o::LinearSolverEigen< MatrixType >::_writeDebug.

{ return _writeDebug;}

Member Data Documentation

template<typename MatrixType >

bool g2o::LinearSolverEigen< MatrixType >::_blockOrdering

protected

Definition at line 134 of file linear_solver_eigen.h.

Referenced by g2o::LinearSolverEigen< MatrixType >::blockOrdering(), and g2o::LinearSolverEigen< MatrixType >::setBlockOrdering().

template<typename MatrixType >

CholeskyDecomposition g2o::LinearSolverEigen< MatrixType >::_cholesky

protected

Definition at line 137 of file linear_solver_eigen.h.

Referenced by g2o::LinearSolverEigen< MatrixType >::solve().

template<typename MatrixType >

bool g2o::LinearSolverEigen< MatrixType >::_init

protected

Definition at line 133 of file linear_solver_eigen.h.

Referenced by g2o::LinearSolverEigen< MatrixType >::init(), and g2o::LinearSolverEigen< MatrixType >::solve().

template<typename MatrixType >

SparseMatrix g2o::LinearSolverEigen< MatrixType >::_sparseMatrix

protected

Definition at line 136 of file linear_solver_eigen.h.

Referenced by g2o::LinearSolverEigen< MatrixType >::solve().

template<typename MatrixType >

bool g2o::LinearSolverEigen< MatrixType >::_writeDebug

protected

Definition at line 135 of file linear_solver_eigen.h.

Referenced by g2o::LinearSolverEigen< MatrixType >::setWriteDebug(), g2o::LinearSolverEigen< MatrixType >::solve(), and g2o::LinearSolverEigen< MatrixType >::writeDebug().

---

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

• /home/xuezhisd/CLionProjects/g2o/g2o/solvers/eigen/linear_solver_eigen.h