Passing argument as Ref<const _MatrixType> instead of const _MatrixType&

[dariomangoni]

I'm a first-time user of Eigen and Spectra, so what I'm experiencing may not be critical, even not correct, anyway...

I need to solve a Generalized Eigenvalue Problem on CSR-format matrices.
I already have CSR arrays so I would like to wrap them in a Eigen::SparseMatrix in order to feed the Spectra::SymGEigsSolver by mean of Spectra::SparseSymMatProd and Spectra::SparseCholesky.

To wrap CSR C-style arrays I would use Map<SparseMatrix<double>>(...CSRmatrixdata...),
but then I came to an empasse. How can I call the SparseSymMatProd constructor? There are few ways...

1. pass the Map<> itself
2. pass a Ref<> to the Map<>
3. pass an Eigen::SparseMatrix& referring to Map<>
4. pass an Eigen::SparseMatrix& referring to Ref<> that refers to the Map<>

What happens?
Fact 1) inside the SparseSymMatProd constructor, a const Eigen::SparseMatrix& is initialized with the argument of the constructor.
i.e. const Eigen::SparseMatrix& m_mat = _argument_

Fact 2a) these passages do not call copy*
Map<>(CSRarrays)
Ref<> = Map<>
and obviously Eigen::SparseMatrix& = Eigen::SparseMatrix&

Fact 2b) these passages do call copy*
Eigen::SparseMatrix& = Ref<>
Eigen::SparseMatrix& = Map<>

*in order to prove that, I compared the valuePtr() address with the original CSR value address

Now:
in [1] and [2]: the reference-initialization inside the constructor forces a copy, but suddenly the destructors of SparseMatrix first and CompressedStorage after are called and they destroy the arrays of m_mat that have just been created! The matrix is invalidated an no product can be actually performed.

in [3] and [4]: everything works fine, the reference-initialization called inside the constructor is of the last case of Fact 2a, so there is no copy and avoid this strange behaviour above, but there do is a copy in the reference initialization that is made before calling the constructor.

So:
I would recommend to change the first lines of Spectra::SparseSymMatProd
from

private:
  const SparseMatrix& m_mat;
public:
    SparseSymMatProd(const SparseMatrix& mat_) :
        m_mat(mat_) { }

to

private:
  Eigen::Ref<const SparseMatrix> m_mat;
public:
    SparseSymMatProd(Eigen::Ref<const SparseMatrix> mat_) :
        m_mat(mat_) { }

This is retro-compatible: it works with Eigen::SparseMatrix
and it should works also with Refs
and it doesn't force a copy.

I'm not aware of all the implications and if this is actually correct.
If you, instead, have a better way to pass some arrays to the eigen solver without this mess and that doesn't force copies, please tell me!

[yixuan]

Hi @dariomangoni , thanks for the comments!
This is indeed a good point. Map<SparseMatrix<double>> and Ref<SparseMatrix<double>> are new features introduced in the recent Eigen 3.3, and they didn't exist before. I need to think of a way to make it compatible with older versions of Eigen, but at least for your own use, it's perfectly a good idea to make that change.

[yixuan]

This has been solved in dec0d1f. Thank you for the suggestion!