BiLinearForm.hpp 4.04 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
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
#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
   **/
  template <class RB, class CB, class T = double>
  class BiLinearForm
      : public MatrixBase<RB,CB,MatrixBackend<BackendTraits<RB,T>>>
  {
    using Self = BiLinearForm;
    using Traits = BackendTraits<RB,T>;
    using Backend = MatrixBackend<Traits>;
    using Super = MatrixBase<RB,CB,Backend>;

  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:
    /// Constructor. Wraps the reference into a non-destroying shared_ptr or moves the basis into a new shared_ptr.
    template <class RB_, class CB_, class Comm_>
    BiLinearForm(RB_&& rowBasis, CB_&& colBasis, Comm_&& comm)
      : Super(FWD(rowBasis), FWD(colBasis), FWD(comm))
    {
      operators_.init(*Super::rowBasis(), *Super::colBasis());
    }

    /// Prepare the matrix for insertion. Clears all the entries.
    void init(SymmetryStructure symmetry = SymmetryStructure::unknown)
    {
      Super::init(symmetry);

      auto const rowSize = this->rowBasis()->localView().maxSize();
      auto const colSize = this->colBasis()->localView().maxSize();
      elementMatrix_.resize(rowSize, colSize);
    }

    /// \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 Expr        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]]
     **/
    // TODO: add method without contextTag.
    template <class ContextTag, class Expr,
              class RowTreePath = RootTreePath, class ColTreePath = RootTreePath>
    void addOperator(ContextTag contextTag, Expr const& expr,
                     RowTreePath row = {}, ColTreePath col = {});

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

    /// Assemble all matrix operators, TODO: incooperate boundary conditions
    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_;
  };

} // end namespace AMDiS

#include <amdis/BiLinearForm.inc.hpp>