BiLinearForm.hpp 5.14 KB
Newer Older
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
#pragma once

#include <cmath>

#include <amdis/LinearAlgebra.hpp>
#include <amdis/OperatorList.hpp>
#include <amdis/common/FlatMatrix.hpp>
#include <amdis/common/TypeTraits.hpp>
#include <amdis/linearalgebra/SymmetryStructure.hpp>
#include <amdis/typetree/TreePath.hpp>

namespace AMDiS
{
  /**
   * 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 T   Coefficient type to store in the matrix
   **/
Praetorius, Simon's avatar
Praetorius, Simon committed
23
  template <class RB, class CB, class T = double, class Traits = BackendTraits<RB,T>>
24
  class BiLinearForm
Praetorius, Simon's avatar
Praetorius, Simon committed
25
      : public MatrixFacade<T, typename Traits::SparsityPattern, Traits::template MatrixImpl>
26
27
  {
    using Self = BiLinearForm;
Praetorius, Simon's avatar
Praetorius, Simon committed
28
    using Super = MatrixFacade<T, typename Traits::SparsityPattern, Traits::template MatrixImpl>;
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45

  public:
    /// The type of the finite element space / basis of the row
    using RowBasis = RB;
    using RowLocalView = typename RowBasis::LocalView;

    /// The type of the finite element space / basis of the column
    using ColBasis = CB;
    using ColLocalView = typename ColBasis::LocalView;

    /// The type of the elements of the DOFVector
    using CoefficientType = typename Traits::CoefficientType;

    /// The type of the matrix filled on an element with local contributions
    using ElementMatrix = FlatMatrix<CoefficientType>;

  public:
46
47
    /// Constructor. Wraps the reference into a non-destroying shared_ptr or moves the basis into
    /// a new shared_ptr.
Praetorius, Simon's avatar
Praetorius, Simon committed
48
49
50
51
    BiLinearForm(std::shared_ptr<RB> const& rowBasis, std::shared_ptr<CB> const& colBasis)
      : Super(*rowBasis, *colBasis)
      , rowBasis_(rowBasis)
      , colBasis_(colBasis)
52
    {
Praetorius, Simon's avatar
Praetorius, Simon committed
53
      operators_.init(*rowBasis, *colBasis);
54

Praetorius, Simon's avatar
Praetorius, Simon committed
55
56
      auto const rowSize = rowBasis->localView().maxSize();
      auto const colSize = colBasis->localView().maxSize();
57
58
59
      elementMatrix_.resize(rowSize, colSize);
    }

Praetorius, Simon's avatar
Praetorius, Simon committed
60
61
62
    std::shared_ptr<RowBasis const> const& rowBasis() const { return rowBasis_; }
    std::shared_ptr<ColBasis const> const& colBasis() const { return colBasis_; }

63
    /// \brief Associate a local operator with this BiLinearForm
64
65
66
67
68
69
70
71
    /**
     * 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
72
73
     *                     or \ref tag::boundary_operator indicating where to assemble
     *                     this operator.
74
     * \tparam Expr        An pre-operator that can be bound to a gridView, or a valid
75
76
77
     *                     GridOperator.
     * \tparam row         A tree-path for the RowBasis
     * \tparam col         A tree-path for the ColBasis
78
79
80
     *
     * [[expects: row is valid tree-path in RowBasis]]
     * [[expects: col is valid tree-path in ColBasis]]
81
     * @{
82
     **/
83
84
    template <class ContextTag, class Expr, class Row, class Col>
    void addOperator(ContextTag contextTag, Expr const& expr, Row row, Col col);
85
86

    // Add an operator to be assembled on the elements of the grid
87
88
    template <class Expr, class Row = RootTreePath, class Col = RootTreePath>
    void addOperator(Expr const& expr, Row row = {}, Col col = {})
89
90
91
92
93
94
    {
      using E = typename RowLocalView::Element;
      addOperator(tag::element_operator<E>{}, expr, row, col);
    }

    // Add an operator to be assembled on the intersections of the grid
95
96
    template <class Expr, class Row = RootTreePath, class Col = RootTreePath>
    void addIntersectionOperator(Expr const& expr, Row row = {}, Col col = {})
97
98
99
100
101
    {
      using I = typename RowLocalView::GridView::Intersection;
      addOperator(tag::intersection_operator<I>{}, expr, row, col);
    }
    /// @}
102
103
104
105
106

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

Praetorius, Simon's avatar
Praetorius, Simon committed
107
    /// Assemble all matrix operators, TODO: incorporate boundary conditions
108
109
110
111
112
113
114
115
116
    void assemble(SymmetryStructure symmetry = SymmetryStructure::unknown);

  protected:

    /// Dense matrix to store coefficients during \ref assemble()
    ElementMatrix elementMatrix_;

    /// List of operators associated to row/col node
    MatrixOperators<RowBasis,ColBasis,ElementMatrix> operators_;
Praetorius, Simon's avatar
Praetorius, Simon committed
117
118
119

    std::shared_ptr<RowBasis const> rowBasis_;
    std::shared_ptr<ColBasis const> colBasis_;
120
121
  };

122
123

#if DUNE_HAVE_CXX_CLASS_TEMPLATE_ARGUMENT_DEDUCTION
124
125
  template <class RB, class CB>
  BiLinearForm(RB&&, CB&&)
126
127
128
129
    -> BiLinearForm<Underlying_t<RB>, Underlying_t<CB>>;
#endif

  template <class RB, class CB, class... Args>
130
  auto makeBiLinearForm(RB&& rowBasis, CB&& colBasis)
131
132
  {
    using BLF = BiLinearForm<Underlying_t<RB>, Underlying_t<CB>>;
133
    return BLF{FWD(rowBasis), FWD(colBasis)};
134
135
  }

136
137
138
} // end namespace AMDiS

#include <amdis/BiLinearForm.inc.hpp>