AMDiS 2.10
The Adaptive Multi-Dimensional Simulation Toolbox
SparsityPattern.hpp
1#pragma once
2
3#include <memory>
4
5#include <dune/istl/matrixindexset.hh>
6
7#include <amdis/common/Index.hpp>
8#include <amdis/linearalgebra/SymmetryStructure.hpp>
9
10namespace AMDiS
11{
15 {
16 public:
17 template <class Basis>
18 SparsityPattern(Basis const& basis, SymmetryStructure symmetry = SymmetryStructure::unknown)
19 : rows_(basis.dimension())
20 , cols_(rows_)
21 , pattern_(rows_, cols_)
22 {
23 init(basis);
24 }
25
26 template <class RowBasis, class ColBasis>
27 SparsityPattern(RowBasis const& rowBasis, ColBasis const& colBasis,
28 SymmetryStructure /*symmetry*/ = SymmetryStructure::unknown)
29 : rows_(rowBasis.dimension())
30 , cols_(colBasis.dimension())
31 , pattern_(rows_, cols_)
32 {
33 if (uintptr_t(&rowBasis) == uintptr_t(&colBasis))
34 init(rowBasis);
35 else
36 init(rowBasis, colBasis);
37 }
38
40 std::size_t rows() const
41 {
42 return rows_;
43 }
44
46 std::size_t cols() const
47 {
48 return cols_;
49 }
50
52 std::size_t avgRowSize() const
53 {
54 assert(rows_ == pattern_.rows());
55 std::size_t sum = 0;
56 for (std::size_t r = 0; r < rows_; ++r)
57 sum += rowSize(r);
58 return (sum + rows_ - 1) / rows_;
59 }
60
62 std::size_t rowSize(std::size_t r) const
63 {
64 return pattern_.rowsize(r);
65 }
66
68 std::size_t nnz() const
69 {
70 return pattern_.size();
71 }
72
74 template <class Matrix>
75 void applyTo(Matrix& matrix) const
76 {
77 pattern_.exportIdx(matrix);
78 }
79
80
81 protected:
82 // Create pattern from basis. This method is called if rowBasis == colBasis.
83 template <class Basis>
84 void init(Basis const& basis)
85 {
86 auto localView = basis.localView();
87 for (const auto& element : elements(basis.gridView())) {
88 localView.bind(element);
89
90 for (std::size_t i = 0, size = localView.size(); i < size; ++i) {
91 // The global index of the i-th vertex of the element
92 auto row = localView.index(i);
93 for (std::size_t j = 0; j < size; ++j) {
94 // The global index of the j-th vertex of the element
95 auto col = localView.index(j);
96 pattern_.add(row, col);
97 }
98 }
99 }
100 }
101
102 // Create pattern from basis. This method is called if rowBasis != colBasis.
103 template <class RowBasis, class ColBasis>
104 void init(RowBasis const& rowBasis, ColBasis const& colBasis)
105 {
106 auto rowLocalView = rowBasis.localView();
107 auto colLocalView = colBasis.localView();
108 for (const auto& element : elements(rowBasis.gridView())) {
109 rowLocalView.bind(element);
110 colLocalView.bind(element);
111
112 for (std::size_t i = 0; i < rowLocalView.size(); ++i) {
113 // The global index of the i-th vertex of the element
114 auto row = rowLocalView.index(i);
115 for (std::size_t j = 0; j < colLocalView.size(); ++j) {
116 // The global index of the j-th vertex of the element
117 auto col = colLocalView.index(j);
118 pattern_.add(row, col);
119 }
120 }
121 }
122 }
123
124 private:
125 std::size_t rows_;
126 std::size_t cols_;
127 Dune::MatrixIndexSet pattern_;
128 };
129
130} // end namespace AMDiS
A general sparsity pattern implementation using the full pattern of the basis by adding all local ind...
Definition: SparsityPattern.hpp:15
std::size_t nnz() const
Total number of non-zeros.
Definition: SparsityPattern.hpp:68
std::size_t rows() const
Number of rows in the matrix.
Definition: SparsityPattern.hpp:40
std::size_t rowSize(std::size_t r) const
Number of non-zeros in row r.
Definition: SparsityPattern.hpp:62
std::size_t avgRowSize() const
Average number of non-zeros per row.
Definition: SparsityPattern.hpp:52
void applyTo(Matrix &matrix) const
Initialize a matrix with the non-zero structure.
Definition: SparsityPattern.hpp:75
std::size_t cols() const
Number of columns in the matrix.
Definition: SparsityPattern.hpp:46