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>
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>
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)
59 , updatePattern_(true)
61 auto const rowSize = rowBasis_->localView().maxSize();
62 auto const colSize = colBasis_->localView().maxSize();
67 template <
class RB_,
class CB_,
68 REQUIRES(Concepts::Similar<RB_,RB>),
69 REQUIRES(Concepts::Similar<CB_,CB>)>
71 :
BiLinearForm(Dune::wrap_or_move(FWD(rowBasis)), Dune::wrap_or_move(FWD(colBasis)))
75 template <
class RB_ = RB,
class CB_ = CB,
76 REQUIRES(std::is_same_v<RB_,CB_>)>
83 REQUIRES(Concepts::Similar<RB_,RB>)>
88 RowBasis const& rowBasis()
const {
return *rowBasis_; }
89 ColBasis const& colBasis()
const {
return *colBasis_; }
111 template <
class ContextTag,
class Operator,
class Row,
class Col>
115 template <
class Operator,
class Row = RootTreePath,
class Col = RootTreePath>
118 using E =
typename RowBasis::LocalView::Element;
119 addOperator(tag::element_operator<E>{}, op, row, col);
123 template <
class Operator,
class Row = RootTreePath,
class Col = RootTreePath>
124 void addIntersectionOperator(Operator
const& op, Row row = {}, Col col = {})
126 using I =
typename RowBasis::LocalView::GridView::Intersection;
127 addOperator(tag::intersection_operator<I>{}, op, row, col);
138 updatePattern_ = updatePattern_ || (
symmetry_ != symm);
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);
168 void init(
bool forcePatternRebuild =
false)
170 if (updatePattern_ || forcePatternRebuild)
173 updatePattern_ =
false;
189 template <
class Basis>
190 void updateImpl2(Basis
const& basis)
194 updatePattern_ =
true;
207 std::shared_ptr<RowBasis const> rowBasis_;
208 std::shared_ptr<ColBasis const> colBasis_;
217 template <
class RB,
class CB>
226 template <
class T =
double,
class RB,
class CB>
227 auto makeBiLinearForm(RB&& rowBasis, CB&& colBasis)
230 return BLF{FWD(rowBasis), FWD(colBasis)};
233 template <
class T =
double,
class RB>
234 auto makeBiLinearForm(RB&& rowBasis)
236 using BLF = BiLinearForm<Underlying_t<RB>, Underlying_t<RB>, T>;
237 return BLF{FWD(rowBasis)};
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
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