AMDiS  0.1
The Adaptive Multi-Dimensional Simulation Toolbox
BiLinearForm.hpp
1 #pragma once
2 
3 #include <cmath>
4 
5 #include <amdis/LinearAlgebra.hpp>
6 #include <amdis/Observer.hpp>
7 #include <amdis/OperatorList.hpp>
8 #include <amdis/common/Concepts.hpp>
9 #include <amdis/common/ConceptsBase.hpp>
10 #include <amdis/common/FlatMatrix.hpp>
11 #include <amdis/common/TypeTraits.hpp>
12 #include <amdis/linearalgebra/SymmetryStructure.hpp>
13 #include <amdis/typetree/TreePath.hpp>
14 
15 namespace AMDiS
16 {
26  template <class RB, class CB, class T = double, class Traits = BackendTraits<RB>>
28  : public MatrixFacade<T, Traits::template MatrixImpl>
29  , private ObserverSequence<event::adapt,2>
30  {
31  using Self = BiLinearForm;
33 
34  public:
36  using RowBasis = RB;
37 
39  using ColBasis = CB;
40 
42  using CoefficientType = T;
43 
46 
48  using SparsityPattern = typename Traits::SparsityPattern;
49 
50  public:
52  BiLinearForm(std::shared_ptr<RB> const& rowBasis, std::shared_ptr<CB> const& colBasis)
53  : Super(*rowBasis, *colBasis)
54  , ObserverSequence<event::adapt,2>(*rowBasis, *colBasis)
55  , symmetry_(SymmetryStructure::unknown)
56  , rowBasis_(rowBasis)
57  , colBasis_(colBasis)
58  , updatePattern_(true)
59  {
60  auto const rowSize = rowBasis_->localView().maxSize();
61  auto const colSize = colBasis_->localView().maxSize();
62  elementMatrix_.resize(rowSize, colSize);
63  }
64 
66  template <class RB_, class CB_,
67  REQUIRES(Concepts::Similar<RB_,RB>),
68  REQUIRES(Concepts::Similar<CB_,CB>)>
69  BiLinearForm(RB_&& rowBasis, CB_&& colBasis)
70  : BiLinearForm(Dune::wrap_or_move(FWD(rowBasis)), Dune::wrap_or_move(FWD(colBasis)))
71  {}
72 
74  template <class RB_ = RB, class CB_ = CB,
75  REQUIRES(std::is_same_v<RB_,CB_>)>
76  explicit BiLinearForm(std::shared_ptr<RB> const& rowBasis)
77  : BiLinearForm(rowBasis, rowBasis)
78  {}
79 
81  template <class RB_,
82  REQUIRES(Concepts::Similar<RB_,RB>)>
83  explicit BiLinearForm(RB_&& rowBasis)
84  : BiLinearForm(Dune::wrap_or_move(FWD(rowBasis)))
85  {}
86 
87  std::shared_ptr<RowBasis const> const& rowBasis() const { return rowBasis_; }
88  std::shared_ptr<ColBasis const> const& colBasis() const { return colBasis_; }
89 
91 
110  template <class ContextTag, class Expr, class Row, class Col>
111  void addOperator(ContextTag contextTag, Expr const& expr, Row row, Col col);
112 
113  // Add an operator to be assembled on the elements of the grid
114  template <class Expr, class Row = RootTreePath, class Col = RootTreePath>
115  void addOperator(Expr const& expr, Row row = {}, Col col = {})
116  {
117  using E = typename RowBasis::LocalView::Element;
118  addOperator(tag::element_operator<E>{}, expr, row, col);
119  }
120 
121  // Add an operator to be assembled on the intersections of the grid
122  template <class Expr, class Row = RootTreePath, class Col = RootTreePath>
123  void addIntersectionOperator(Expr const& expr, Row row = {}, Col col = {})
124  {
125  using I = typename RowBasis::LocalView::GridView::Intersection;
126  addOperator(tag::intersection_operator<I>{}, expr, row, col);
127  }
129 
131 
135  void setSymmetryStructure(SymmetryStructure symm)
136  {
137  updatePattern_ = (symmetry_ != symm);
138  symmetry_ = symm;
139  }
140 
142  void setSymmetryStructure(std::string str)
143  {
144  setSymmetryStructure(symmetryStructure(str));
145  }
146 
148  template <class RowLocalView, class ColLocalView,
149  REQUIRES(Concepts::LocalView<RowLocalView>),
150  REQUIRES(Concepts::LocalView<ColLocalView>)>
151  void assemble(RowLocalView const& rowLocalView,
152  ColLocalView const& colLocalView);
153 
155  void assemble();
156 
159 
167  void init(bool forcePatternRebuild = false)
168  {
169  if (updatePattern_ || forcePatternRebuild)
170  {
171  Super::init(SparsityPattern{*rowBasis_, *colBasis_, symmetry_});
172  updatePattern_ = false;
173  }
174  else
175  {
176  Super::init();
177  }
178  }
179 
182  void updateImpl(event::adapt e, index_t<0> i) override { updatePattern_ = true; }
183  void updateImpl(event::adapt e, index_t<1> i) override { updatePattern_ = true; }
184 
185  protected:
187  SymmetryStructure symmetry_;
188 
191 
193  MatrixOperators<RowBasis,ColBasis,ElementMatrix> operators_;
194 
195  std::shared_ptr<RowBasis const> rowBasis_;
196  std::shared_ptr<ColBasis const> colBasis_;
197 
198  private:
199  // The pattern is rebuilt on calling init if this flag is set
200  bool updatePattern_;
201  };
202 
203 
204  // deduction guides
205  template <class RB, class CB>
206  BiLinearForm(RB&&, CB&&)
207  -> BiLinearForm<Underlying_t<RB>, Underlying_t<CB>>;
208 
209  template <class RB>
210  BiLinearForm(RB&&)
211  -> BiLinearForm<Underlying_t<RB>, Underlying_t<RB>>;
212 
213 
214  template <class T = double, class RB, class CB>
215  auto makeBiLinearForm(RB&& rowBasis, CB&& colBasis)
216  {
217  using BLF = BiLinearForm<Underlying_t<RB>, Underlying_t<CB>, T>;
218  return BLF{FWD(rowBasis), FWD(colBasis)};
219  }
220 
221  template <class T = double, class RB>
222  auto makeBiLinearForm(RB&& rowBasis)
223  {
224  using BLF = BiLinearForm<Underlying_t<RB>, Underlying_t<RB>, T>;
225  return BLF{FWD(rowBasis)};
226  }
227 
228 } // end namespace AMDiS
229 
230 #include <amdis/BiLinearForm.inc.hpp>
void addOperator(ContextTag contextTag, Expr const &expr, Row row, Col col)
Associate a local operator with this BiLinearForm.
void setSymmetryStructure(std::string str)
Set the symmetry structure using a string.
Definition: BiLinearForm.hpp:142
Definition: AdaptiveGrid.hpp:373
MatrixOperators< RowBasis, ColBasis, ElementMatrix > operators_
List of operators associated to row/col node.
Definition: BiLinearForm.hpp:193
BiLinearForm(std::shared_ptr< RB > const &rowBasis)
Constructor for rowBasis == colBasis.
Definition: BiLinearForm.hpp:76
Contains all classes needed for solving linear and non linear equation systems.
Definition: AdaptBase.hpp:6
void resize(size_type r, size_type c, value_type v={})
Resizes the container to contain r x c elements.
Definition: FlatMatrix.hpp:70
RB RowBasis
The type of the finite element space / basis of the row.
Definition: BiLinearForm.hpp:36
BiLinearForm(std::shared_ptr< RB > const &rowBasis, std::shared_ptr< CB > const &colBasis)
Constructor. Stores the row and column basis in a local shared_ptr to const.
Definition: BiLinearForm.hpp:52
BiLinearForm(RB_ &&rowBasis)
Wraps the passed row-basis into a (non-destroying) shared_ptr.
Definition: BiLinearForm.hpp:83
void setSymmetryStructure(SymmetryStructure symm)
Provide some additional information about the structure of the matrix pattern.
Definition: BiLinearForm.hpp:135
void assemble()
Assemble all matrix operators, TODO: incorporate boundary conditions.
Definition: BiLinearForm.inc.hpp:69
T CoefficientType
The type of the elements of the DOFVector.
Definition: BiLinearForm.hpp:42
Definition: BiLinearForm.hpp:27
void init(bool forcePatternRebuild=false)
Definition: BiLinearForm.hpp:167
Definition: Observer.hpp:25
BiLinearForm(RB_ &&rowBasis, CB_ &&colBasis)
Wraps the passed global bases into (non-destroying) shared_ptr.
Definition: BiLinearForm.hpp:69
std::integral_constant< std::size_t, I > index_t
A wrapper for std::size_t type.
Definition: Index.hpp:31
CB ColBasis
The type of the finite element space / basis of the column.
Definition: BiLinearForm.hpp:39
void updateImpl(event::adapt e, index_t< 0 > i) override
Definition: BiLinearForm.hpp:182
Definition: OperatorList.hpp:14
Definition: OperatorList.hpp:15
void init()
Initialize the matrix for insertion while keeping the pattern unchanged.
Definition: MatrixFacade.hpp:49
Definition: LinearSolverInterface.hpp:20
ElementMatrix elementMatrix_
Dense matrix to store coefficients during assemble()
Definition: BiLinearForm.hpp:190
SymmetryStructure symmetry_
The symmetry property if the bilinear form.
Definition: BiLinearForm.hpp:187
typename Traits::SparsityPattern SparsityPattern
Type of the sparsity pattern of the backend.
Definition: BiLinearForm.hpp:48