AMDiS 2.10
The Adaptive Multi-Dimensional Simulation Toolbox
BiLinearForm.hpp
1#pragma once
2
3#include <cmath>
4
5#include <amdis/Assembler.hpp>
6#include <amdis/LinearAlgebra.hpp>
7#include <amdis/Observer.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
15namespace AMDiS
16{
27 template <class RB, class CB, class T = double, class Traits = BackendTraits>
29 : public MatrixFacade<T, Traits::template Matrix<RB,CB>::template Impl>
30 , private ObserverSequence<event::adapt,2>
31 {
32 using Self = BiLinearForm;
34
35 public:
37 using RowBasis = RB;
38
40 using ColBasis = CB;
41
43 using CoefficientType = T;
44
47
49 using SparsityPattern = typename Traits::template Matrix<RB,CB>::SparsityPattern;
50
51 public:
53 BiLinearForm(std::shared_ptr<RB> const& rowBasis, std::shared_ptr<CB> const& colBasis)
54 : Super(*rowBasis, *colBasis)
55 , ObserverSequence<event::adapt,2>(*rowBasis, *colBasis)
56 , symmetry_(SymmetryStructure::unknown)
57 , rowBasis_(rowBasis)
58 , colBasis_(colBasis)
59 , updatePattern_(true)
60 {
61 auto const rowSize = rowBasis_->localView().maxSize();
62 auto const colSize = colBasis_->localView().maxSize();
63 elementMatrix_.resize(rowSize, colSize);
64 }
65
67 template <class RB_, class CB_,
68 REQUIRES(Concepts::Similar<RB_,RB>),
69 REQUIRES(Concepts::Similar<CB_,CB>)>
70 BiLinearForm(RB_&& rowBasis, CB_&& colBasis)
71 : BiLinearForm(Dune::wrap_or_move(FWD(rowBasis)), Dune::wrap_or_move(FWD(colBasis)))
72 {}
73
75 template <class RB_ = RB, class CB_ = CB,
76 REQUIRES(std::is_same_v<RB_,CB_>)>
77 explicit BiLinearForm(std::shared_ptr<RB> const& rowBasis)
78 : BiLinearForm(rowBasis, rowBasis)
79 {}
80
82 template <class RB_,
83 REQUIRES(Concepts::Similar<RB_,RB>)>
84 explicit BiLinearForm(RB_&& rowBasis)
85 : BiLinearForm(Dune::wrap_or_move(FWD(rowBasis)))
86 {}
87
88 RowBasis const& rowBasis() const { return *rowBasis_; }
89 ColBasis const& colBasis() const { return *colBasis_; }
90
92
111 template <class ContextTag, class Operator, class Row, class Col>
112 void addOperator(ContextTag contextTag, Operator const& op, Row row, Col col);
113
114 // Add an operator to be assembled on the elements of the grid
115 template <class Operator, class Row = RootTreePath, class Col = RootTreePath>
116 void addOperator(Operator const& op, Row row = {}, Col col = {})
117 {
118 using E = typename RowBasis::LocalView::Element;
119 addOperator(tag::element_operator<E>{}, op, row, col);
120 }
121
122 // Add an operator to be assembled on the intersections of the grid
123 template <class Operator, class Row = RootTreePath, class Col = RootTreePath>
124 void addIntersectionOperator(Operator const& op, Row row = {}, Col col = {})
125 {
126 using I = typename RowBasis::LocalView::GridView::Intersection;
127 addOperator(tag::intersection_operator<I>{}, op, row, col);
128 }
130
132
136 void setSymmetryStructure(SymmetryStructure symm)
137 {
138 updatePattern_ = updatePattern_ || (symmetry_ != symm);
139 symmetry_ = symm;
140 }
141
143 void setSymmetryStructure(std::string str)
144 {
145 setSymmetryStructure(symmetryStructure(str));
146 }
147
149 template <class RowLocalView, class ColLocalView,
150 REQUIRES(Concepts::LocalView<RowLocalView>),
151 REQUIRES(Concepts::LocalView<ColLocalView>)>
152 void assemble(RowLocalView const& rowLocalView,
153 ColLocalView const& colLocalView);
154
156 void assemble();
157
160
168 void init(bool forcePatternRebuild = false)
169 {
170 if (updatePattern_ || forcePatternRebuild)
171 {
172 Super::init(SparsityPattern{*rowBasis_, *colBasis_, symmetry_});
173 updatePattern_ = false;
174 }
175 else
176 {
177 Super::init();
178 }
179 }
180
183 void updateImpl(event::adapt, index_t<0>) override { updateImpl2(*rowBasis_); }
184 void updateImpl(event::adapt, index_t<1>) override { updateImpl2(*colBasis_); }
185
186 auto& assembler() { return assembler_; }
187
188 private:
189 template <class Basis>
190 void updateImpl2(Basis const& basis)
191 {
192 if (!updatePattern_)
193 assembler_.update(basis.gridView());
194 updatePattern_ = true;
195 }
196
197 protected:
199 SymmetryStructure symmetry_;
200
203
206
207 std::shared_ptr<RowBasis const> rowBasis_;
208 std::shared_ptr<ColBasis const> colBasis_;
209
210 private:
211 // The pattern is rebuilt on calling init if this flag is set
212 bool updatePattern_;
213 };
214
215
216 // deduction guides
217 template <class RB, class CB>
218 BiLinearForm(RB&&, CB&&)
219 -> BiLinearForm<Underlying_t<RB>, Underlying_t<CB>>;
220
221 template <class RB>
222 BiLinearForm(RB&&)
223 -> BiLinearForm<Underlying_t<RB>, Underlying_t<RB>>;
224
225
226 template <class T = double, class RB, class CB>
227 auto makeBiLinearForm(RB&& rowBasis, CB&& colBasis)
228 {
229 using BLF = BiLinearForm<Underlying_t<RB>, Underlying_t<CB>, T>;
230 return BLF{FWD(rowBasis), FWD(colBasis)};
231 }
232
233 template <class T = double, class RB>
234 auto makeBiLinearForm(RB&& rowBasis)
235 {
236 using BLF = BiLinearForm<Underlying_t<RB>, Underlying_t<RB>, T>;
237 return BLF{FWD(rowBasis)};
238 }
239
240} // end namespace AMDiS
241
242#include <amdis/BiLinearForm.inc.hpp>
An Assembler is the main driver for building the local element matrices and vectors by assembling ope...
Definition: Assembler.hpp:132
Definition: BiLinearForm.hpp:31
void setSymmetryStructure(std::string str)
Set the symmetry structure using a string.
Definition: BiLinearForm.hpp:143
SymmetryStructure symmetry_
The symmetry property if the bilinear form.
Definition: BiLinearForm.hpp:199
void addOperator(ContextTag contextTag, Operator const &op, Row row, Col col)
Associate a local operator with this BiLinearForm.
ElementMatrix elementMatrix_
Dense matrix to store coefficients during assemble()
Definition: BiLinearForm.hpp:202
MatrixAssembler< RowBasis, ColBasis, ElementMatrix > assembler_
List of operators associated to row/col node.
Definition: BiLinearForm.hpp:205
T CoefficientType
The type of the elements of the DOFVector.
Definition: BiLinearForm.hpp:43
typename Traits::template Matrix< RB, CB >::SparsityPattern SparsityPattern
Type of the sparsity pattern of the backend.
Definition: BiLinearForm.hpp:49
void setSymmetryStructure(SymmetryStructure symm)
Provide some additional information about the structure of the matrix pattern.
Definition: BiLinearForm.hpp:136
void updateImpl(event::adapt, index_t< 0 >) override
Definition: BiLinearForm.hpp:183
void init(bool forcePatternRebuild=false)
Definition: BiLinearForm.hpp:168
BiLinearForm(RB_ &&rowBasis, CB_ &&colBasis)
Wraps the passed global bases into (non-destroying) shared_ptr.
Definition: BiLinearForm.hpp:70
BiLinearForm(RB_ &&rowBasis)
Wraps the passed row-basis into a (non-destroying) shared_ptr.
Definition: BiLinearForm.hpp:84
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:53
void assemble()
Assemble all matrix operators, TODO: incorporate boundary conditions.
Definition: BiLinearForm.inc.hpp:63
CB ColBasis
The type of the finite element space / basis of the column.
Definition: BiLinearForm.hpp:40
RB RowBasis
The type of the finite element space / basis of the row.
Definition: BiLinearForm.hpp:37
BiLinearForm(std::shared_ptr< RB > const &rowBasis)
Constructor for rowBasis == colBasis.
Definition: BiLinearForm.hpp:77
void resize(size_type r, size_type c, value_type v={})
Resizes the container to contain r x c elements.
Definition: FlatMatrix.hpp:70
Definition: MatrixFacade.hpp:25
void init()
Initialize the matrix for insertion while keeping the pattern unchanged.
Definition: MatrixFacade.hpp:49
The base class for an operator to be used in an Assembler.
Definition: Operator.hpp:80
Definition: Observer.hpp:25