AMDiS  0.3
The Adaptive Multi-Dimensional Simulation Toolbox
Integrate.hpp
1 #pragma once
2 
3 #include <type_traits>
4 #include <dune/geometry/quadraturerules.hh>
5 #include <dune/grid/common/partitionset.hh>
6 #include <amdis/GridFunctions.hpp>
7 #include <amdis/common/Order.hpp>
8 
9 namespace AMDiS
10 {
11  namespace Impl
12  {
13  template <class GF, class GV, class QP>
14  auto integrateImpl(GF&& gf, GV const& gridView, QP makeQuad)
15  {
16  auto localFct = localFunction(FWD(gf));
17 
18  using GridFct = remove_cvref_t<GF>;
19  using Range = typename GridFct::Range;
20  Range result(0);
21 
22  for (auto const& element : elements(gridView, Dune::Partitions::interior)) {
23  auto geometry = element.geometry();
24 
25  localFct.bind(element);
26  auto const& quad = makeQuad(element.type(), localFct);
27 
28  for (auto const qp : quad)
29  result += localFct(qp.position()) * geometry.integrationElement(qp.position()) * qp.weight();
30  localFct.unbind();
31  }
32 
33  return gridView.comm().sum(result);
34  }
35 
36  } // end namespace Impl
37 
38 
40 
48  template <class Expr, class GridView>
49  auto integrate(Expr&& expr, GridView const& gridView)
50  {
51  auto&& gridFct = makeGridFunction(FWD(expr), gridView);
52 
53  using LocalFct = TYPEOF(localFunction(gridFct));
54  static_assert(Concepts::Polynomial<LocalFct>, "Polynomial degree of expression can not be deduced. You need to provide an explicit value for the quadrature degree or a quadrature rule in `integrate()`.");
55 
56  using Rules = Dune::QuadratureRules<typename GridView::ctype, GridView::dimension>;
57  if constexpr(Concepts::Polynomial<LocalFct>)
58  return Impl::integrateImpl(FWD(gridFct), gridView,
59  [](auto&& t, auto&& lf) { return Rules::rule(t, order(lf)); });
60  else
61  return 0.0;
62  }
63 
64 
66 
73  template <class Expr, class GridView,
74  class QuadratureRule = Dune::QuadratureRule<typename GridView::ctype, GridView::dimension>>
75  auto integrate(Expr&& expr, GridView const& gridView, QuadratureRule const& quad)
76  {
77  auto&& gridFct = makeGridFunction(FWD(expr), gridView);
78  return Impl::integrateImpl(FWD(gridFct), gridView,
79  [&](auto&&, auto&&) { return quad; });
80  }
81 
82 
84 
90  template <class Expr, class GridView>
91  auto integrate(Expr&& expr, GridView const& gridView, int degree,
92  Dune::QuadratureType::Enum qt = Dune::QuadratureType::GaussLegendre)
93  {
94  using Rules = Dune::QuadratureRules<typename GridView::ctype, GridView::dimension>;
95 
96  auto&& gridFct = makeGridFunction(FWD(expr), gridView);
97  return Impl::integrateImpl(FWD(gridFct), gridView,
98  [&](auto const& type, auto&&) { return Rules::rule(type, degree, qt); });
99  }
100 
101 } // end namespace AMDiS
decltype(auto) makeGridFunction(PreGridFct const &preGridFct, GridView const &gridView)
Generator for Gridfunctions from Expressions (PreGridfunctions)
Definition: GridFunction.hpp:154
Contains all classes needed for solving linear and non linear equation systems.
Definition: AdaptBase.hpp:6
auto order(F const &f) -> decltype(&F::operator(), f.order())
polynomial order of functions
Definition: Order.hpp:11
auto integrate(Expr &&expr, GridView const &gridView)
Integrate expression with quadrature rule determined by polynomial order of expression.
Definition: Integrate.hpp:49