AMDiS  0.3
The Adaptive Multi-Dimensional Simulation Toolbox
FirstOrderTestGradTrial.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 test_gradtrial {};
19  struct grad_trial {};
20  }
21 
22 
24  template <class LC, class GridFct>
25  class GridFunctionOperator<tag::test_gradtrial, LC, GridFct>
26  : public GridFunctionOperatorBase<GridFunctionOperator<tag::test_gradtrial, LC, GridFct>, LC, GridFct>
27  {
28  static const int dow = ContextGeometry<LC>::dow;
29  using Self = GridFunctionOperator;
31 
32  static_assert( static_size_v<typename GridFct::Range> == dow, "Expression must be of vector type." );
33 
34  public:
36  : Super(expr, 1)
37  {}
38 
39  template <class CG, class RN, class CN, class Mat>
40  void getElementMatrix(CG const& contextGeo, RN const& rowNode, CN const& colNode, Mat& elementMatrix)
41  {
42  static_assert(RN::isLeaf && CN::isLeaf,
43  "Operator can be applied to Leaf-Nodes only.");
44 
45  auto const& quad = this->getQuadratureRule(contextGeo.type(), rowNode, colNode);
46  std::size_t rowSize = rowNode.size();
47  std::size_t colSize = colNode.size();
48 
49  using RangeFieldType = typename CN::LocalBasis::Traits::RangeFieldType;
50  using WorldVector = FieldVector<RangeFieldType,CG::dow>;
51  std::vector<WorldVector> colGradients;
52 
53  for (auto const& qp : quad) {
54  // Position of the current quadrature point in the reference element
55  auto&& local = contextGeo.local(qp.position());
56 
57  // The transposed inverse Jacobian of the map from the reference element to the element
58  const auto jacobian = contextGeo.geometry().jacobianInverseTransposed(local);
59 
60  // The multiplicative factor in the integral transformation formula
61  const auto factor = contextGeo.integrationElement(qp.position()) * qp.weight();
62  const auto b = Super::coefficient(local);
63 
64  // the values of the shape functions on the reference element at the quadrature point
65  auto const& shapeValues = rowNode.localBasisValuesAt(local);
66 
67  // The gradients of the shape functions on the reference element
68  auto const& shapeGradients = colNode.localBasisJacobiansAt(local);
69 
70  // Compute the shape function gradients on the real element
71  colGradients.resize(shapeGradients.size());
72 
73  for (std::size_t i = 0; i < colGradients.size(); ++i)
74  jacobian.mv(shapeGradients[i][0], colGradients[i]);
75 
76  for (std::size_t j = 0; j < colSize; ++j) {
77  const auto local_j = colNode.localIndex(j);
78  const auto value = factor * (b * colGradients[j]);
79  for (std::size_t i = 0; i < rowSize; ++i) {
80  const auto local_i = rowNode.localIndex(i);
81  elementMatrix[local_i][local_j] += value * shapeValues[i];
82  }
83  }
84  }
85 
86  }
87  };
88 
89 
91  template <class Expr, class... QuadratureArgs>
92  auto fot(Expr&& expr, tag::grad_trial, QuadratureArgs&&... args)
93  {
94  return makeOperator(tag::test_gradtrial{}, FWD(expr), FWD(args)...);
95  }
96 
99 } // end namespace AMDiS
Definition: FirstOrderTestGradTrial.hpp:18
auto makeOperator(Tag tag, Expr &&expr, QuadratureArgs &&... args)
Store tag and expression into a PreGridFunctionOperator to create a GridFunctionOperator.
Definition: GridFunctionOperator.hpp:220
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
auto fot(Expr &&expr, tag::grad_test, QuadratureArgs &&... args)
Create a first-order term with derivative on trial-function.
Definition: FirstOrderGradTestTrial.hpp:40
The main implementation of an CRTP-base class for operators using a grid-function coefficient to be u...
Definition: GridFunctionOperator.hpp:39
Definition: FirstOrderTestGradTrial.hpp:19