Assembler.hpp 4.9 KB
Newer Older
1
2
3
#pragma once

#include <memory>
4
#include <tuple>
5

6
7
8
#include <dune/common/fmatrix.hh>
#include <dune/common/fvector.hh>

9
#include <dune/amdis/FiniteElementSpaces.hpp>
10
11
#include <dune/amdis/LinearAlgebra.hpp>
#include <dune/amdis/common/Mpl.hpp>
12
#include <dune/amdis/common/TypeDefs.hpp>
13
14
15

namespace AMDiS
{
16
  template <class FeSpaces>
17
18
  class Assembler
  {
19
20
21
    template <std::size_t I>
    using FeSpace = std::tuple_element_t<I, FeSpaces>;

22
23
24
    /// Number of problem components
    static constexpr int nComponents = std::tuple_size<FeSpaces>::value;

25
26
27
    /// The grid view the global FE basis lives on
    using GridView = typename FeSpace<0>::GridView;

28
29
30
31
32
    using SystemMatrixType = SystemMatrix<FeSpaces>;
    using SystemVectorType = SystemVector<FeSpaces>;

    template <class T>
    using MatrixEntries = Dune::FieldMatrix<T, nComponents, nComponents>;
33

34
35
    template <class T>
    using VectorEntries = Dune::FieldVector<T, nComponents>;
36

37
38
39
    using ElementVector = Impl::ElementVector;
    using ElementMatrix = Impl::ElementMatrix;

40
  public:
41
    /// Constructor, stores a shared-pointer to the feSpaces
42
43
44
45
46
    Assembler(std::shared_ptr<FeSpaces> const& feSpaces, bool asmMatrix, bool asmVector)
      : globalBases_(feSpaces)
      , asmMatrix_(asmMatrix)
      , asmVector_(asmVector)
    {}
47

48
49
    /// Assemble the linear system
    template <class Operators>
50
    void assemble(
51
        GridView const& gv,
52
53
54
55
        SystemMatrixType& matrix,
        SystemVectorType& solution,
        SystemVectorType& rhs,
        MatrixEntries<Operators>& matrix_operators,
56
        VectorEntries<Operators>& rhs_operators);
57
58

  private:
59
    /// Sets the system to zero and initializes all operators and boundary conditions
60
    template <class Operators>
61
62
63
64
65
    void initMatrixVector(
        SystemMatrixType& matrix,
        SystemVectorType& solution,
        SystemVectorType& rhs,
        MatrixEntries<Operators>& matrix_operators,
66
        VectorEntries<Operators>& rhs_operators) const;
67
68


69
70
71
    /// Assemble a block-matrix of element-matrices and return a matrix of flags, whether
    /// a given block has received some entries.
    template <class Operators>
72
73
    MatrixEntries<bool> assembleElementMatrices(
        MatrixEntries<Operators>& operators,
74
        MatrixEntries<ElementMatrix>& elementMatrix) const;
75
76
77
78

    /// Assemble a block-vector of element-vectors and return a vector of flags, whether
    /// a given block has received some entries.
    template <class Operators>
79
80
    VectorEntries<bool> assembleElementVectors(
        VectorEntries<Operators>& operators,
81
        VectorEntries<ElementVector>& elementVector) const;
82
83


84
85
86
    /// Assemble one block of the block-element-matrix
    // The MatrixData argument stores all matrix-operators
    template <class Operators, class RowView, class ColView>
87
88
89
90
    bool assembleElementMatrix(
        Operators& operators,
        ElementMatrix& elementMatrix,
        RowView const& rowLocalView,
91
        ColView const& colLocalView) const;
92
93
94
95

    /// Assemble one block of the block-element-vector
    // The VectorData argument stores all vector-operators
    template <class Operators, class RowView>
96
97
98
    bool assembleElementVector(
        Operators& operators,
        ElementVector& elementVector,
99
        RowView const& rowLocalView) const;
100
101
102


    /// Add the block-element-matrix to the system-matrix
103
104
105
    void addElementMatrices(
        SystemMatrixType& dofmatrix,
        MatrixEntries<bool> const& addMat,
106
        MatrixEntries<ElementMatrix> const& elementMatrix) const;
107
108

    /// Add the block-element-vector to the system-vector
109
110
111
    void addElementVectors(
        SystemVectorType& dofvector,
        VectorEntries<bool> const& addVec,
112
        VectorEntries<ElementVector> const& elementVector) const;
113
114
115
116

    /// Finish insertion into the matrix and assembles boundary conditions
    /// Return the number of nonzeros assembled into the matrix
    template <class Operators>
117
118
119
120
121
    std::size_t finishMatrixVector(
        SystemMatrixType& matrix,
        SystemVectorType& solution,
        SystemVectorType& rhs,
        MatrixEntries<Operators>& matrix_operators,
122
        VectorEntries<Operators>& rhs_operators) const;
123
124
125
126
127


    /// Return whether the matrix-block needs to be assembled
    template <std::size_t R, std::size_t C, class Operators>
    bool assembleMatrix(const index_t<R>, const index_t<C>,
128
                        MatrixEntries<Operators> const& matrix_operators) const
129
130
131
132
133
134
    {
      return asmMatrix_ && (!matrix_operators[R][C].assembled || matrix_operators[R][C].changing);
    }

    /// Return whether the vector-block needs to be assembled
    template <std::size_t R, class Operators>
135
    bool assembleVector(const index_t<R>, VectorEntries<Operators> const& rhs_operators) const
136
137
138
    {
      return asmVector_ && (!rhs_operators[R].assembled || rhs_operators[R].changing);
    }
139
140

  private:
141
142
143
    FiniteElementSpaces<FeSpaces> globalBases_;
    bool asmMatrix_;
    bool asmVector_;
144
145
146
  };

} // end namespace AMDiS
147
148

#include "Assembler.inc.hpp"