AMDiS  0.3
The Adaptive Multi-Dimensional Simulation Toolbox
LinearForm.hpp
1 #pragma once
2 
3 #include <dune/common/typeutilities.hh>
4 
5 #include <amdis/LinearAlgebra.hpp>
6 #include <amdis/OperatorList.hpp>
7 #include <amdis/common/Concepts.hpp>
8 #include <amdis/common/ConceptsBase.hpp>
9 #include <amdis/common/FlatVector.hpp>
10 #include <amdis/common/TypeTraits.hpp>
11 #include <amdis/typetree/TreePath.hpp>
12 
13 namespace AMDiS
14 {
16 
24  template <class GB, class T = double, class Traits = BackendTraits<GB>>
25  class LinearForm
26  : public VectorFacade<T, Traits::template VectorImpl>
27  {
29  using Self = LinearForm;
30 
31  public:
33  using GlobalBasis = GB;
34 
36  using CoefficientType = T;
37 
40 
41  public:
43  explicit LinearForm(std::shared_ptr<GB> const& basis)
44  : Super(*basis)
45  , basis_(basis)
46  {
47  auto const localSize = basis->localView().maxSize();
48  elementVector_.resize(localSize);
49  }
50 
52  template <class GB_,
53  REQUIRES(Concepts::Similar<GB_,GB>)>
54  explicit LinearForm(GB_&& basis)
55  : LinearForm(Dune::wrap_or_move(FWD(basis)))
56  {}
57 
58  GlobalBasis const& basis() const { return *basis_; }
59 
61 
78  template <class ContextTag, class Operator, class TreePath>
79  void addOperator(ContextTag contextTag, Operator const& op, TreePath path);
80 
81  // Add an operator to be assembled on the elements of the grid
82  template <class Operator, class TreePath = RootTreePath>
83  void addOperator(Operator const& op, TreePath path = {})
84  {
85  using E = typename GlobalBasis::LocalView::Element;
87  }
88 
89  // Add an operator to be assembled on the intersections of the grid
90  template <class Operator, class TreePath = RootTreePath>
91  void addIntersectionOperator(Operator const& op, TreePath path = {})
92  {
93  using I = typename GlobalBasis::LocalView::GridView::Intersection;
95  }
97 
98 
100  template <class LocalView,
101  REQUIRES(Concepts::LocalView<LocalView>)>
102  void assemble(LocalView const& localView);
103 
105  void assemble();
106 
107  private:
109  ElementVector elementVector_;
110 
112  VectorOperators<GlobalBasis,ElementVector> operators_;
113 
114  std::shared_ptr<GlobalBasis const> basis_;
115  };
116 
117 
118  // deduction guide
119  template <class GB>
120  LinearForm(GB&& basis)
122 
123 
124  template <class T = double, class GB>
125  auto makeLinearForm(GB&& basis)
126  {
127  using LF = LinearForm<Underlying_t<GB>, T>;
128  return LF{FWD(basis)};
129  }
130 
131 } // end namespace AMDiS
132 
133 #include <amdis/LinearForm.inc.hpp>
std::size_t localSize() const
Return the number of entries in the local part of the vector.
Definition: VectorFacade.hpp:82
void addOperator(ContextTag contextTag, Operator const &op, TreePath path)
Associate a local operator with this LinearForm.
The basic container that stores a base vector and a corresponding basis.
Definition: LinearSolverInterface.hpp:23
LinearForm(std::shared_ptr< GB > const &basis)
Constructor. Stores the shared_ptr of the basis and creates a new DataTransfer.
Definition: LinearForm.hpp:43
Definition: AdaptiveGrid.hpp:373
The basic container that stores a base vector and a corresponding basis.
Definition: LinearForm.hpp:25
Contains all classes needed for solving linear and non linear equation systems.
Definition: AdaptBase.hpp:6
constexpr bool LocalView
A Dune::Functions::LocalView type.
Definition: Concepts.hpp:182
typename GridView::template Codim< 0 >::Entity Element
Type of the grid element we are bound to.
Definition: LocalView.hpp:35
void assemble()
Assemble all vector operators added by addOperator().
Definition: LinearForm.inc.hpp:61
GB GlobalBasis
The type of the functionspace basis associated to this linearform.
Definition: LinearForm.hpp:33
LinearForm(GB_ &&basis)
Wraps the passed global basis into a (non-destroying) shared_ptr.
Definition: LinearForm.hpp:54
Definition: OperatorList.hpp:14
T CoefficientType
The type of the elements of the DOFVector.
Definition: LinearForm.hpp:36
Definition: OperatorList.hpp:15