AMDiS  2.10
The Adaptive Multi-Dimensional Simulation Toolbox
FirstOrderTestvecGradTrial.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 testvec_gradtrial {};
18  }
19 
20 
23  {
24  public:
26 
27  template <class CG, class RN, class CN, class Quad, class LocalFct, class Mat>
28  void assemble(CG const& contextGeo, RN const& rowNode, CN const& colNode,
29  Quad const& quad, LocalFct const& localFct, Mat& elementMatrix) const
30  {
31  static_assert(static_size_v<typename LocalFct::Range> == 1,
32  "Expression must be of scalar type.");
33  static_assert(RN::isPower && CN::isLeaf,
34  "col-node must be Leaf-Node and row-node must be a Power-Node.");
35 
36  assert( rowNode.degree() == CG::dow );
37 
38  std::size_t rowSize = rowNode.child(0).size();
39  std::size_t colSize = colNode.size();
40 
41  using RangeFieldType = typename CN::LocalBasis::Traits::RangeFieldType;
42  using WorldVector = FieldVector<RangeFieldType,CG::dow>;
43  std::vector<WorldVector> colGradients;
44 
45  for (auto const& qp : quad) {
46  // Position of the current quadrature point in the reference element
47  decltype(auto) local = contextGeo.coordinateInElement(qp.position());
48 
49  // The transposed inverse Jacobian of the map from the reference element to the element
50  const auto jacobian = contextGeo.elementGeometry().jacobianInverseTransposed(local);
51 
52  // The multiplicative factor in the integral transformation formula
53  const auto factor = localFct(local) * contextGeo.integrationElement(qp.position()) * qp.weight();
54 
55  // the values of the shape functions on the reference element at the quadrature point
56  auto const& shapeValues = rowNode.child(0).localBasisValuesAt(local);
57 
58  // The gradients of the shape functions on the reference element
59  auto const& shapeGradients = colNode.localBasisJacobiansAt(local);
60 
61  // Compute the shape function gradients on the real element
62  colGradients.resize(shapeGradients.size());
63 
64  for (std::size_t i = 0; i < colGradients.size(); ++i)
65  jacobian.mv(shapeGradients[i][0], colGradients[i]);
66 
67  for (std::size_t i = 0; i < rowSize; ++i) {
68  const auto value = factor * shapeValues[i];
69  for (std::size_t j = 0; j < colSize; ++j) {
70  const auto local_j = colNode.localIndex(j);
71  for (std::size_t k = 0; k < rowNode.degree(); ++k) {
72  const auto local_ki = rowNode.child(k).localIndex(i);
73  elementMatrix[local_ki][local_j] += value * colGradients[j][k];
74  }
75  }
76  }
77  }
78 
79  }
80  };
81 
82  template <class LC>
83  struct GridFunctionOperatorRegistry<tag::testvec_gradtrial, LC>
84  {
85  static constexpr int degree = 1;
87  };
88 
91 } // end namespace AMDiS
Definition: AdaptBase.hpp:6
Definition: FirstOrderTestvecGradTrial.hpp:17
first-order operator
Definition: FirstOrderTestvecGradTrial.hpp:22
Registry to specify a tag for each implementation type.
Definition: GridFunctionOperator.hpp:216