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> 26 template <
class RB,
class CB,
class T =
double,
class Traits = BackendTraits<RB>>
29 ,
private ObserverSequence<event::adapt,2>
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)
58 , updatePattern_(true)
60 auto const rowSize = rowBasis_->localView().maxSize();
61 auto const colSize = colBasis_->localView().maxSize();
66 template <
class RB_,
class CB_,
67 REQUIRES(Concepts::Similar<RB_,RB>),
68 REQUIRES(Concepts::Similar<CB_,CB>)>
74 template <
class RB_ = RB,
class CB_ = CB,
75 REQUIRES(std::is_same_v<RB_,CB_>)>
82 REQUIRES(Concepts::Similar<RB_,RB>)>
87 std::shared_ptr<RowBasis const>
const& rowBasis()
const {
return rowBasis_; }
88 std::shared_ptr<ColBasis const>
const& colBasis()
const {
return colBasis_; }
110 template <
class ContextTag,
class Expr,
class Row,
class Col>
111 void addOperator(ContextTag contextTag, Expr
const& expr, Row row, Col col);
114 template <
class Expr,
class Row = RootTreePath,
class Col = RootTreePath>
115 void addOperator(Expr
const& expr, Row row = {}, Col col = {})
117 using E =
typename RowBasis::LocalView::Element;
122 template <
class Expr,
class Row = RootTreePath,
class Col = RootTreePath>
123 void addIntersectionOperator(Expr
const& expr, Row row = {}, Col col = {})
125 using I =
typename RowBasis::LocalView::GridView::Intersection;
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);
167 void init(
bool forcePatternRebuild =
false)
169 if (updatePattern_ || forcePatternRebuild)
172 updatePattern_ =
false;
195 std::shared_ptr<RowBasis const> rowBasis_;
196 std::shared_ptr<ColBasis const> colBasis_;
205 template <
class RB,
class CB>
214 template <
class T =
double,
class RB,
class CB>
215 auto makeBiLinearForm(RB&& rowBasis, CB&& colBasis)
218 return BLF{FWD(rowBasis), FWD(colBasis)};
221 template <
class T =
double,
class RB>
222 auto makeBiLinearForm(RB&& rowBasis)
225 return BLF{FWD(rowBasis)};
230 #include <amdis/BiLinearForm.inc.hpp>
Definition: AdaptiveGrid.hpp:373
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
Definition: Observer.hpp:25
std::integral_constant< std::size_t, I > index_t
A wrapper for std::size_t type.
Definition: Index.hpp:31
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