1

这是一个晦涩难懂的问题,我真的不希望有人回答,但我有这个方法,它接受(并返回)一个 Eigen::SparseMatrix。我想把它放到 deal.ii 库中,有没有办法从 deal.ii/Eigen 复制/转换 SparseMatrix?我知道您可以将 deal.ii 复制到 Trilinos SparseMatrix,例如:

  `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);`    

有没有类似的方式 Eigen::SparseMatrix?我猜 Eigen 在 deal.ii 中并没有真正得到那种支持。因此,也许有一些“蛮力”类型的方法,例如这种对代码的尝试显然不起作用:

`

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; 

`

好的,所以我想出了如何从 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; 
  }

不,我只需要弄清楚相反的情况。

4

1 回答 1

2

这个问题很老,但也许我仍然可以提供帮助。我是 deal.II 开发人员之一,我不记得在邮件列表中看到过这个(对于这些类型的问题,它比 SO 更活跃)。

deal.II 中的 ASparseMatrix不存储自己的稀疏模式:相反,它存储指向SparsityPattern对象的指针。您需要对特征矩阵进行两次循环:一次用于设置,SparsityPattern第二次用于复制矩阵值。像下面这样的东西似乎工作:

#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
  }
}

您还必须管理SparsityPattern对象的生命周期。

deal.II 不使用 CSR 或 CSC:它使用自己的类似 CSR 的格式,其中主对角线上的条目首先存储在包含该行的矩阵条目的数组中,因此我们确实需要使用迭代器接口进行复制.

于 2018-01-27T04:54:44.040 回答