AMDiS  0.3
The Adaptive Multi-Dimensional Simulation Toolbox
ZeroOrderTest.hpp
1 #pragma once
2 
3 #include <type_traits>
4 
5 #include <amdis/GridFunctionOperator.hpp>
6 #include <amdis/common/StaticSize.hpp>
7 
8 namespace AMDiS
9 {
15  namespace tag
16  {
17  struct test {};
18  }
19 
20 
22  template <class LC, class GridFct>
23  class GridFunctionOperator<tag::test, LC, GridFct>
24  : public GridFunctionOperatorBase<GridFunctionOperator<tag::test, LC, GridFct>, LC, GridFct>
25  {
26  using Self = GridFunctionOperator;
28 
29  static_assert( static_size_v<typename GridFct::Range> == 1, "Expression must be of scalar type." );
30 
31  public:
33  : Super(expr, 0)
34  {}
35 
36  template <class CG, class Node, class Vec>
37  void getElementVector(CG const& contextGeo, Node const& node, Vec& elementVector)
38  {
39  static_assert(Node::isLeaf, "Operator can be applied to Leaf-Nodes only");
40 
41  auto const& quad = this->getQuadratureRule(contextGeo.type(), node);
42  std::size_t size = node.size();
43 
44  for (auto const& qp : quad) {
45  // Position of the current quadrature point in the reference element
46  auto&& local = contextGeo.local(qp.position());
47 
48  // The multiplicative factor in the integral transformation formula
49  auto factor = contextGeo.integrationElement(qp.position()) * qp.weight();
50  factor *= Super::coefficient(local);
51 
52  auto const& shapeValues = node.localBasisValuesAt(local);
53  for (std::size_t i = 0; i < size; ++i) {
54  const auto local_i = node.localIndex(i);
55  elementVector[local_i] += factor * shapeValues[i];
56  }
57  }
58  }
59  };
60 
63 } // end namespace AMDiS
The base-template for GridFunctionOperators.
Definition: GridFunctionOperator.hpp:242
Contains all classes needed for solving linear and non linear equation systems.
Definition: AdaptBase.hpp:6
Definition: ZeroOrderTest.hpp:17
The main implementation of an CRTP-base class for operators using a grid-function coefficient to be u...
Definition: GridFunctionOperator.hpp:39