DOFMatrixBase.hpp 5.35 KB
Newer Older
1
2
#pragma once

3
#include <cmath>
4
5
6

#include <dune/common/timer.hh>

7
#include <amdis/OperatorList.hpp>
8
#include <amdis/common/Math.hpp>
9
10
#include <amdis/typetree/MultiIndex.hpp>
#include <amdis/typetree/TreePath.hpp>
11

12
13
namespace AMDiS
{
14
15
16
17
18
19
20
21
22
  /**
   * Basis implementation of DOFMatrix, i.e. a sparse matrix storing all the
   * assembled Operators indexed with DOF indices. The matrix data is associated
   * to a row and column global basis.
   *
   * \tparam RB  Basis of the matrix rows
   * \tparam CB  Basis of matrix columns
   * \tparam Backend  A linear-algebra backend for the matrix storage
   **/
23
  template <class RowBasisType, class ColBasisType, class Backend>
Praetorius, Simon's avatar
Praetorius, Simon committed
24
  class DOFMatrixBase
25
  {
Praetorius, Simon's avatar
Praetorius, Simon committed
26
27
28
  public:
    /// The type of the finite element space / basis of the row
    using RowBasis = RowBasisType;
29
    using RowLocalView = typename RowBasis::LocalView;
Praetorius, Simon's avatar
Praetorius, Simon committed
30
31
32

    /// The type of the finite element space / basis of the column
    using ColBasis = ColBasisType;
33
34
35
36
    using ColLocalView = typename ColBasis::LocalView;

    using Element = typename RowLocalView::Element;
    using Geometry = typename Element::Geometry;
Praetorius, Simon's avatar
Praetorius, Simon committed
37
38
39
40

    /// The index/size - type
    using size_type = typename RowBasis::size_type;

41
    using BaseMatrix = typename Backend::BaseMatrix;
42
    using ElementMatrix = Dune::DynamicMatrix<double>;
43

Praetorius, Simon's avatar
Praetorius, Simon committed
44
45
46
47
48
  public:
    /// Constructor. Constructs new BaseVector.
    DOFMatrixBase(RowBasis const& rowBasis, ColBasis const& colBasis)
      : rowBasis_(&rowBasis)
      , colBasis_(&colBasis)
49
    {
50
      operators_.init(rowBasis, colBasis);
51
    }
Praetorius, Simon's avatar
Praetorius, Simon committed
52

53
54
55
    DOFMatrixBase(DOFMatrixBase const&) = delete;
    DOFMatrixBase(DOFMatrixBase&&) = delete;

Praetorius, Simon's avatar
Praetorius, Simon committed
56
57
58
59
60
61
62
63
64
65
66
67
    /// Return the row-basis \ref rowBasis of the matrix
    RowBasis const& rowBasis() const
    {
      return *rowBasis_;
    }

    /// Return the col-basis \ref colBasis of the matrix
    ColBasis const& colBasis() const
    {
      return *colBasis_;
    }

68
69
70
71
72
73
74
75
76
77
78
79
    /// Return the data-matrix
    BaseMatrix const& matrix() const
    {
      return backend_.matrix();
    }

    /// Return the data-matrix
    BaseMatrix& matrix()
    {
      return backend_.matrix();
    }

Praetorius, Simon's avatar
Praetorius, Simon committed
80
81
82
83
84
85
86
87
88
89
90
91
    /// Return the size of the \ref rowBasis_
    size_type rows() const
    {
      return rowBasis_->dimension();
    }

    /// Return the size of the \ref colBasis_
    size_type cols() const
    {
      return colBasis_->dimension();
    }

92
93
    /// Initialize the matrix for insertion, e.g. allocate the non-zero pattern
    /// If \p setToZero is true, the matrix is set to 0
94
    void init(bool asmMatrix);
95
96

    /// Finish the matrix insertion, e.g. cleanup or final insertion
97
    void finish(bool asmMatrix);
98
99

    /// Insert a block of values into the matrix (add to existing values)
100
101
    void insert(RowLocalView const& rowLocalView,
                ColLocalView const& colLocalView,
102
                ElementMatrix const& elementMatrix);
103
104
105
106
107
108
109
110

    /// Insert a single value into the matrix (add to existing value)
    template <class RowIndex, class ColIndex>
    void insert(RowIndex row, ColIndex col, typename Backend::value_type const& value)
    {
      backend_.insert(flatMultiIndex(row), flatMultiIndex(col), value);
    }

111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
    /// \brief Associate a local operator with this DOFMatrix
    /**
     * Stores an operator in a list that gets assembled during a call to \ref assemble().
     * The operator may be assigned to a specific context, i.e. either an element
     * operator, an intersection operator, or a boundary operator.
     * The \p row and \p col tree paths specify the sub-basis for test and trial
     * functions the operator is applied to.
     *
     * \tparam ContextTag  One of \ref tag::element_operator, \ref tag::intersection_operator
     *                     or \ref tag::boundary_operator indicating where to assemble this operator.
     * \tparam PreOperator  An pre-operator that can be bound to a gridView, or a valid
     *                      GridOperator.
     * \tparam row  A tree-path for the RowBasis
     * \tparam col  A tree-path for the ColBasis
     *
     * [[expects: row is valid tree-path in RowBasis]]
     * [[expects: col is valid tree-path in ColBasis]]
     **/
129
130
131
    template <class ContextTag, class Operator,
              class RowTreePath = RootTreePath, class ColTreePath = RootTreePath>
    void addOperator(ContextTag contextTag, Operator const& preOp,
132
                     RowTreePath row = {}, ColTreePath col = {});
133

134
135
136
    /// Assemble the matrix operators on the bound element.
    void assemble(RowLocalView const& rowLocalView,
                  ColLocalView const& colLocalView);
137

138
139
    /// Assemble all matrix operators, TODO: incooperate boundary conditions
    void assemble();
140

141
142
143
144
145
146
147
148
149
150
151
152
    /// \brief Deletes all rows with \p dirichletNodes != 0 and if \p apply is true, adds
    /// a one on the diagonal of the matrix.
    auto applyDirichletBC(std::vector<char> const& dirichletNodes)
    {
      return backend_.applyDirichletBC(dirichletNodes);
    }

    size_type nnz() const
    {
      return backend_.nnz();
    }

Praetorius, Simon's avatar
Praetorius, Simon committed
153
  protected:
154
155
156
    /// The finite element space / basis associated with the rows/columns
    RowBasis const* rowBasis_;
    ColBasis const* colBasis_;
157
158

    /// Data backend
159
    Backend backend_;
160

161
162
163
164
165
    /// Dense matrix to store coefficients during \ref assemble()
    ElementMatrix elementMatrix_;

    /// List of operators associated to row/col node
    MatrixOperators<RowBasis,ColBasis> operators_;
166
167
168
  };

} // end namespace AMDiS
169
170

#include "DOFMatrixBase.inc.hpp"