AMDiS  0.3
The Adaptive Multi-Dimensional Simulation Toolbox
ZeroOrderTestTrialvec.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_trialvec {};
18  }
19 
20 
22  template <class LC, class GridFct>
23  class GridFunctionOperator<tag::test_trialvec, LC, GridFct>
24  : public GridFunctionOperatorBase<GridFunctionOperator<tag::test_trialvec, LC, GridFct>, LC, GridFct>
25  {
26  static const int dow = ContextGeometry<LC>::dow;
27  using Self = GridFunctionOperator;
29 
30  static_assert( static_size_v<typename GridFct::Range> == dow, "Expression must be of vector type." );
31 
32  public:
34  : Super(expr, 0)
35  {}
36 
37  template <class CG, class RN, class CN, class Mat>
38  void getElementMatrix(CG const& contextGeo, RN const& rowNode, CN const& colNode, Mat& elementMatrix)
39  {
40  static_assert(RN::isLeaf && CN::isPower,
41  "RN must be Leaf-Node and CN must be a Power-Node.");
42 
43  static const std::size_t CHILDREN = CN::CHILDREN;
44  static_assert( static_size_v<typename GridFct::Range> == CHILDREN, "" );
45 
46  auto const& quad = this->getQuadratureRule(contextGeo.type(), rowNode, colNode);
47  std::size_t rowSize = rowNode.size();
48  std::size_t colSize = colNode.child(0).size();
49 
50  for (auto const& qp : quad) {
51  // Position of the current quadrature point in the reference element
52  auto&& local = contextGeo.local(qp.position());
53 
54  // The multiplicative factor in the integral transformation formula
55  const auto factor = contextGeo.integrationElement(qp.position()) * qp.weight();
56  const auto b = Super::coefficient(local);
57 
58  auto const& rowShapeValues = rowNode.localBasisValuesAt(local);
59  auto const& colShapeValues = colNode.child(0).localBasisValuesAt(local);
60 
61  for (std::size_t i = 0; i < rowSize; ++i) {
62  const auto local_i = rowNode.localIndex(i);
63 
64  for (std::size_t j = 0; j < colSize; ++j) {
65  const auto value = b * (factor * rowShapeValues[i] * colShapeValues[j]);
66 
67  for (std::size_t k = 0; k < CHILDREN; ++k) {
68  const auto local_kj = colNode.child(k).localIndex(j);
69  elementMatrix[local_i][local_kj] += at(value,k);
70  }
71  }
72  }
73  }
74  }
75  };
76 
79 } // 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
Wrapper class for element and geometry.
Definition: ContextGeometry.hpp:43
Definition: ZeroOrderTestTrialvec.hpp:17
The main implementation of an CRTP-base class for operators using a grid-function coefficient to be u...
Definition: GridFunctionOperator.hpp:39