Converting Eigen::SparseMatrix to deal.ii ::SparseMatrix?

Question

This is kind of an obscure question and I don't really expect anyone to answer, but I have this method that takes (and returns) an Eigen::SparseMatrix. I want to put it into the deal.ii library, is there a way to copy/convert a SparseMatrix from deal.ii/Eigen? I know you can copy deal.ii to Trilinos SparseMatrix something like:

  `SparseMatrix<double> matrix(sparsity);
...//fill matrix

  Epetra_Map map(TrilinosWrappers::types::int_type(5),
                 TrilinosWrappers::types::int_type(5),
                 0,
                 Utilities::Trilinos::comm_world());

  TrilinosWrappers::SparseMatrix tmatrix;
  tmatrix.reinit (map, map, matrix, 0, false);`

Is there a similar way Eigen::SparseMatrix? I guess Eigen don't really have that kind of support in deal.ii. So perhaps there is some 'brute force' type method, like this attempt at code which obviously doesn't work:

`

Eigen::SparseMatrix<double> ConvertToEigenMatrix(SparseMatrix<double> data)
{
    Eigen::SparseMatrix<double> eMatrix(data.m(), data.n());
    for (int i = 0; i < data.m(); ++i)
        eMatrix.row(i) =  Eigen::SparseMatrix<double> ::Map(&data[i][0], data.n());
    return eMatrix;

`

Ok, so I figured out how to convert from dealii::SparseMatrix -> Eigen::SparseMatrix.

  SparseMatrix<double>::iterator smi = matrix.begin();
  SparseMatrix<double>::iterator smi_end = matrix.end();

  unsigned int row,col;
  double val;
  for (; smi!=smi_end; ++smi)
  {
       row = smi->row();
       col = smi->column();
       val = smi->value();

       spMat.insert(row, col) = val; 
       std::cout << val << std::endl; 
  }

No, I just need to figure out the reverse.

I don't know about `deal.ii` but if there is a way to pass to it raw buffer in CRS or CCS format, then you can pass a `Eigen::SparseMatrix` too by accessing the buffers with `outerIndexPtr()`, `innerIndexPtr()`, and `valuesPtr()` methods. — ggael, Dec 08 '17 at 13:21

David Wells · Accepted Answer · 2018-01-27T14:49:29.467

This question is old but maybe I can still help. I am one of the deal.II developers and I don't remember seeing this on the mailing list (which is much more active for these types of questions than SO).

A SparseMatrix in deal.II does not store its own sparsity pattern: instead, it stores a pointer to a SparsityPattern object. You'll need to loop over the eigen matrix twice: once to set up the SparsityPattern and a second time to copy matrix values. Something like the following seems to work:

#include <deal.II/lac/dynamic_sparsity_pattern.h>
#include <deal.II/lac/sparsity_pattern.h>
#include <deal.II/lac/sparse_matrix.h>

#include <eigen3/Eigen/Sparse>

#include <iostream>

int main()
{
  const std::size_t shape = 3;
  Eigen::SparseMatrix<double> matrix(shape, shape);
  matrix.insert(0, 0) = 1.0;
  matrix.insert(0, 1) = 2.0;
  matrix.insert(0, 2) = 1.0;
  matrix.insert(2, 2) = 2.0;
  matrix.makeCompressed();

  {
    dealii::SparsityPattern sparsity_pattern(matrix.rows(), matrix.cols());
    dealii::DynamicSparsityPattern dynamic_sparsity_pattern(matrix.rows(), matrix.cols());

    for (decltype(matrix.outerSize()) row_n = 0; row_n < matrix.outerSize(); ++row_n)
      for (Eigen::SparseMatrix<double>::InnerIterator it(matrix, row_n); it; ++it)
        dynamic_sparsity_pattern.add(it.row(), it.col());

    sparsity_pattern.copy_from(dynamic_sparsity_pattern);
    dealii::SparseMatrix<double> matrix2(sparsity_pattern);

    for (decltype(matrix.outerSize()) row_n = 0; row_n < matrix.outerSize(); ++row_n)
      for (Eigen::SparseMatrix<double>::InnerIterator it(matrix, row_n); it; ++it)
        matrix2.set(it.row(), it.col(), it.value());

    matrix2.print(std::cout); // prints the right matrix
  }
}

You will have to manage the lifetime of the SparsityPattern object too.

deal.II does not use CSR or CSC: it uses its own CSR-like format where the entry on the main diagonal is stored first in the array containing the matrix entries for that row, so we really do need to copy with the iterator interfaces.

Converting Eigen::SparseMatrix to deal.ii ::SparseMatrix?

1 Answers1