Eigen: Compare two sparse matrix with probably different sparsity pattern

Question

I want to compare two Eigen::SparseMatrix

There are exists res.isApprox(ans) method but unfortunately it fail with assertion in case of different sparsity pattern, as far as i know

AffineInvariantDeformerTest01: Eigen/src/SparseCore/SparseMatrix.h:934: void Eigen::internal::set_from_triplets(const InputIterator&, const InputIterator&, SparseMatrixType&, DupFunctor) [with InputIterator = __gnu_cxx::__normal_iterator<Eigen::Triplet<double, int>*, std::vector<Eigen::Triplet<double, int>, std::allocator<Eigen::Triplet<double, int> > > >; SparseMatrixType = Eigen::SparseMatrix<double, 0, int>; DupFunctor = Eigen::internal::scalar_sum_op<double, double>]: Assertion `it->row()>=0 && it->row()<mat.rows() && it->col()>=0 && it->col()<mat.cols()' failed.

I want it for unit tests so it's ok if it will be not so fast. Probably i can transform my sparse matrices to dense, but i hope that there are exist more elegant solution

EDIT: Just a method which can compare sparsity patterns of two matrices is ok for me too

Answer 1

First of all, isApprox() works fine, even with different sparsity patterns. It looks like your error is elsewhere (likely already when you set your matrix using setFromTriplets())

When given two matrices with different sparsity patterns they will be considered approximately equal, if the differing entries are (almost) zero. The following should calculate true twice (increase the 1e-19 to see a difference):

#include <Eigen/SparseCore>
#include <iostream>
#include <array>

int main() {
    Eigen::SparseMatrix<double> Mat1(2,2), Mat2(2,2);

    std::array<Eigen::Triplet<double,int>, 2> data1  {{{0,0,1.0}, {1,1, 1e-19}}};
    std::array<Eigen::Triplet<double,int>, 2> data2  {{{0,0,1.0}, {1,0, 1e-19}}};

    Mat1.setFromTriplets(data1.begin(), data1.end());
    Mat2.setFromTriplets(data2.begin(), data2.end());

    std::cout << "Mat1.isApprox(Mat1) == " << Mat1.isApprox(Mat1) << "\nMat1.isApprox(Mat2) == " << Mat1.isApprox(Mat2) << "\n";
}

If you want to compare the sparsity patterns, you can check the data starting at the inner and outer index pointers (using a Map for both). The following works if both matrices have the same type and are compressed:

template<class Derived>
bool hasSamePattern(Eigen::SparseCompressedBase<Derived> const& A, Eigen::SparseCompressedBase<Derived> const& B)
{
    assert(A.isCompressed() && B.isCompressed());
    if(A.rows() != B.rows() || A.cols() != B.cols() || A.nonZeros() != B.nonZeros())
        return false;
    typedef Eigen::Matrix<typename Derived::StorageIndex, Eigen::Dynamic, 1> IndexVector;
    Eigen::Index outerSize = A.outerSize(), nnz = A.nonZeros();
    if(IndexVector::Map(A.outerIndexPtr(), outerSize) != IndexVector::Map(B.outerIndexPtr(), outerSize))
        return false;
    if(IndexVector::Map(A.innerIndexPtr(), nnz) != IndexVector::Map(B.innerIndexPtr(), nnz))
        return false;

    return true;
}