AMDiS  0.3
The Adaptive Multi-Dimensional Simulation Toolbox
SecondOrderPartialTestPartialTrial.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 #include <amdis/utility/LocalToGlobalAdapter.hpp>
9 
10 namespace AMDiS
11 {
17  namespace tag
18  {
20  {
21  std::size_t comp_test; // i
22  std::size_t comp_trial; // j
23  };
24  }
25 
26 
28  template <class LC, class GridFct>
29  class GridFunctionOperator<tag::partialtest_partialtrial, LC, GridFct>
30  : public GridFunctionOperatorBase<GridFunctionOperator<tag::partialtest_partialtrial, LC, GridFct>, LC, GridFct>
31  {
32  using Self = GridFunctionOperator;
34 
35  static_assert( static_size_v<typename GridFct::Range> == 1, "Expression must be of scalar type." );
36 
37  public:
39  : Super(expr, 2)
40  , compTest_(tag.comp_test)
41  , compTrial_(tag.comp_trial)
42  {}
43 
44  template <class CG, class RN, class CN, class Mat>
45  void getElementMatrix(CG const& contextGeo, RN const& rowNode, CN const& colNode, Mat& elementMatrix)
46  {
47  static_assert(RN::isLeaf && CN::isLeaf, "Operator can be applied to Leaf-Nodes only.");
48 
49  auto const& quad = this->getQuadratureRule(contextGeo.type(), rowNode, colNode);
50 
51  LocalToGlobalBasisAdapter rowLocalBasis(rowNode, contextGeo.geometry());
52  LocalToGlobalBasisAdapter colLocalBasis(colNode, contextGeo.geometry());
53 
54  for (auto const& qp : quad) {
55  // Position of the current quadrature point in the reference element
56  decltype(auto) local = contextGeo.local(qp.position());
57 
58  // The multiplicative factor in the integral transformation formula
59  auto factor = contextGeo.integrationElement(qp.position()) * qp.weight();
60  factor *= Super::coefficient(local);
61 
62  auto const& rowPartial = rowLocalBasis.partialsAt(local, compTest_);
63  auto const& colPartial = colLocalBasis.partialsAt(local, compTrial_);
64 
65  for (std::size_t j = 0; j < colLocalBasis.size(); ++j) {
66  const auto local_j = colNode.localIndex(j);
67  const auto value = factor * colPartial[j];
68  for (std::size_t i = 0; i < rowLocalBasis.size(); ++i) {
69  const auto local_i = rowNode.localIndex(i);
70  elementMatrix[local_i][local_j] += value * rowPartial[i];
71  }
72  }
73  }
74 
75  }
76 
77  private:
78  std::size_t compTest_;
79  std::size_t compTrial_;
80  };
81 
82 
84  template <class Expr, class... QuadratureArgs>
85  auto sot_ij(Expr&& expr, std::size_t comp_test, std::size_t comp_trial, QuadratureArgs&&... args)
86  {
87  return makeOperator(tag::partialtest_partialtrial{comp_test, comp_trial}, FWD(expr), FWD(args)...);
88  }
89 
92 } // end namespace AMDiS
auto makeOperator(Tag tag, Expr &&expr, QuadratureArgs &&... args)
Store tag and expression into a PreGridFunctionOperator to create a GridFunctionOperator.
Definition: GridFunctionOperator.hpp:220
Convert a simple (scalar) local basis into a global basis.
Definition: LocalToGlobalAdapter.hpp:64
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
auto sot_ij(Expr &&expr, std::size_t comp_test, std::size_t comp_trial, QuadratureArgs &&... args)
Create a second-order term of partial derivatives.
Definition: SecondOrderPartialTestPartialTrial.hpp:85
Definition: SecondOrderPartialTestPartialTrial.hpp:19
The main implementation of an CRTP-base class for operators using a grid-function coefficient to be u...
Definition: GridFunctionOperator.hpp:39