4 #include <dune/geometry/quadraturerules.hh> 5 #include <dune/grid/common/partitionset.hh> 6 #include <amdis/GridFunctions.hpp> 7 #include <amdis/common/Order.hpp> 13 template <
class GF,
class GV,
class QP>
14 auto integrateImpl(GF&& gf, GV
const& gridView, QP makeQuad)
16 auto localFct = localFunction(FWD(gf));
18 using GridFct = remove_cvref_t<GF>;
19 using Range =
typename GridFct::Range;
22 for (
auto const& element : elements(gridView, Dune::Partitions::interior)) {
23 auto geometry = element.geometry();
25 localFct.bind(element);
26 auto const& quad = makeQuad(element.type(), localFct);
28 for (
auto const qp : quad)
29 result += localFct(qp.position()) * geometry.integrationElement(qp.position()) * qp.weight();
33 return gridView.comm().sum(result);
48 template <
class Expr,
class Gr
idView>
49 auto integrate(Expr&& expr, GridView
const& gridView)
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()`.");
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)); });
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)
78 return Impl::integrateImpl(FWD(gridFct), gridView,
79 [&](
auto&&,
auto&&) {
return quad; });
90 template <
class Expr,
class Gr
idView>
91 auto integrate(Expr&& expr, GridView
const& gridView,
int degree,
92 Dune::QuadratureType::Enum qt = Dune::QuadratureType::GaussLegendre)
94 using Rules = Dune::QuadratureRules<typename GridView::ctype, GridView::dimension>;
97 return Impl::integrateImpl(FWD(gridFct), gridView,
98 [&](
auto const& type,
auto&&) {
return Rules::rule(type, degree, qt); });
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