AMDiS  0.3
The Adaptive Multi-Dimensional Simulation Toolbox
MatrixBackend.hpp
1 #pragma once
2 
3 #include <list>
4 #include <memory>
5 #include <string>
6 #include <vector>
7 
8 #include <Eigen/SparseCore>
9 
10 #include <dune/common/timer.hh>
11 
12 #include <amdis/Output.hpp>
13 
14 namespace AMDiS
15 {
18  template <class T, int Orientation = Eigen::RowMajor>
20  {
21  public:
23  using BaseMatrix = Eigen::SparseMatrix<T, Orientation>;
24 
26  using value_type = typename BaseMatrix::Scalar;
27 
29  using size_type = typename BaseMatrix::Index;
30 
31  public:
33  template <class Basis>
34  explicit EigenSparseMatrix(Basis const&, Basis const&) {}
35 
38  {
39  return matrix_;
40  }
41 
43  BaseMatrix const& matrix() const
44  {
45  return matrix_;
46  }
47 
48 
52  void insert(size_type r, size_type c, value_type const& value)
53  {
54  test_exit_dbg(r < matrix_.rows() && c < matrix_.cols(),
55  "Indices out of range [0,{})x[0,{})", matrix_.rows(), matrix_.cols());
56  tripletList_.emplace_back(r, c, value);
57  }
58 
59  template <class Ind, class LocalMat>
60  void scatter(Ind const& idx, LocalMat const& mat)
61  {
62  scatter(idx, idx, mat);
63  }
64 
65  template <class RowInd, class ColInd, class LocalMat>
66  void scatter(RowInd const& rows, ColInd const& cols, LocalMat const& mat)
67  {
68  //tripletList_.reserve(tripletList_.size() + rows.size()*cols.size());
69  for (size_type i = 0; i < size_type(rows.size()); ++i)
70  for (size_type j = 0; j < size_type(cols.size()); ++j)
71  tripletList_.emplace_back(rows[i], cols[j], mat[i][j]);
72  }
73 
74 
76  template <class Pattern>
77  void init(Pattern const& pattern)
78  {
79  matrix_.resize(pattern.rows(), pattern.cols());
80  matrix_.setZero();
81  }
82 
84  void init()
85  {
86  matrix_.setZero();
87  }
88 
90  void finish()
91  {
92  matrix_.setFromTriplets(tripletList_.begin(), tripletList_.end());
93  matrix_.makeCompressed();
94 
95  tripletList_.clear(); // NOTE: this does not free the memory
96  }
97 
99  std::size_t nnz() const
100  {
101  return matrix_.nonZeros();
102  }
103 
104  private:
106  BaseMatrix matrix_;
107 
109  std::vector<Eigen::Triplet<value_type, Eigen::Index>> tripletList_;
110  };
111 
112 } // end namespace AMDiS
void init()
Set all matrix entries to zero while keeping the size unchanged.
Definition: MatrixBackend.hpp:84
void init(Pattern const &pattern)
Resize the matrix according to the pattern provided and set all entries to zero.
Definition: MatrixBackend.hpp:77
std::size_t nnz() const
Return the number of nonzeros in the matrix.
Definition: MatrixBackend.hpp:99
typename BaseMatrix::Index size_type
The index/size - type.
Definition: MatrixBackend.hpp:29
void test_exit_dbg(bool condition, Args &&... args)
call assert_msg, in debug mode only
Definition: Output.hpp:205
Contains all classes needed for solving linear and non linear equation systems.
Definition: AdaptBase.hpp:6
The basic container that stores a base matrix and a corresponding row/column feSpace.
Definition: MatrixBackend.hpp:19
BaseMatrix & matrix()
Return a reference to the data-matrix matrix.
Definition: MatrixBackend.hpp:37
BaseMatrix const & matrix() const
Return a reference to the data-matrix matrix.
Definition: MatrixBackend.hpp:43
EigenSparseMatrix(Basis const &, Basis const &)
Constructor. Constructs new BaseMatrix.
Definition: MatrixBackend.hpp:34
void finish()
Set the matrix entries from the triplet list and compress the matrix afterwards.
Definition: MatrixBackend.hpp:90
typename BaseMatrix::Scalar value_type
The type of the elements of the DOFMatrix.
Definition: MatrixBackend.hpp:26
void insert(size_type r, size_type c, value_type const &value)
Returns an update-proxy of the inserter, to inserte/update a value at position (r, c) in the matrix. Need an insertionMode inserter, that can be created by init and must be closed by finish after insertion.
Definition: MatrixBackend.hpp:52
Eigen::SparseMatrix< T, Orientation > BaseMatrix
The matrix type of the underlying base matrix.
Definition: MatrixBackend.hpp:23