AMDiS  0.3
The Adaptive Multi-Dimensional Simulation Toolbox
FirstOrderGradTest.hpp
1 #pragma once
2 
3 #include <type_traits>
4 
5 #include <amdis/GridFunctionOperator.hpp>
6 #include <amdis/Output.hpp>
7 #include <amdis/common/StaticSize.hpp>
8 
9 namespace AMDiS
10 {
16  namespace tag
17  {
18  struct gradtest {};
19  }
20 
21 
23  template <class LC, class GridFct>
24  class GridFunctionOperator<tag::gradtest, LC, GridFct>
25  : public GridFunctionOperatorBase<GridFunctionOperator<tag::gradtest, LC, GridFct>, LC, GridFct>
26  {
27  static const int dow = ContextGeometry<LC>::dow;
28  using Self = GridFunctionOperator;
30 
31  static_assert( static_size_v<typename GridFct::Range> == dow, "Expression must be of vector type." );
32 
33  public:
35  : Super(expr, 1)
36  {}
37 
38  template <class CG, class Node, class Vec>
39  void getElementVector(CG const& contextGeo, Node const& node, Vec& elementVector)
40  {
41  static_assert(Node::isLeaf, "Node must be Leaf-Node.");
42 
43  auto const& quad = this->getQuadratureRule(contextGeo.type(), node);
44  std::size_t feSize = node.size();
45 
46  using RangeFieldType = typename Node::LocalBasis::Traits::RangeFieldType;
47  using WorldVector = FieldVector<RangeFieldType,CG::dow>;
48  std::vector<WorldVector> gradients;
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 transposed inverse Jacobian of the map from the reference element to the element
55  const auto jacobian = contextGeo.geometry().jacobianInverseTransposed(local);
56 
57  // The multiplicative factor in the integral transformation formula
58  const auto factor = Super::coefficient(local);
59  const auto dx = contextGeo.integrationElement(qp.position()) * qp.weight();
60 
61  // The gradients of the shape functions on the reference element
62  auto const& shapeGradients = node.localBasisJacobiansAt(local);
63 
64  // Compute the shape function gradients on the real element
65  gradients.resize(shapeGradients.size());
66 
67  for (std::size_t i = 0; i < gradients.size(); ++i)
68  jacobian.mv(shapeGradients[i][0], gradients[i]);
69 
70  for (std::size_t i = 0; i < feSize; ++i) {
71  const auto local_i = node.localIndex(i);
72  elementVector[local_i] += dx * (factor * gradients[i]);
73  }
74  }
75 
76  }
77  };
78 
81 } // end namespace AMDiS
The base-template for GridFunctionOperators.
Definition: GridFunctionOperator.hpp:242
Definition: FirstOrderGradTest.hpp:18
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
The main implementation of an CRTP-base class for operators using a grid-function coefficient to be u...
Definition: GridFunctionOperator.hpp:39