AMDiS 2.10
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
7#include <amdis/GridFunctions.hpp>
8#include <amdis/Output.hpp>
9#include <amdis/common/Order.hpp>
10
11namespace AMDiS
12{
15
24 template <class Expr, class GV>
25 auto integrate(Expr&& expr, GV const& gridView, int localFctOrder = -1)
26 {
27 auto gridFct = makeGridFunction(FWD(expr), gridView);
28 auto localFct = localFunction(gridFct);
29
30 using LF = TYPEOF(localFct);
31 auto quadOrder = [localFctOrder](auto const& lf) {
32 DUNE_UNUSED_PARAMETER(localFctOrder);
33 if constexpr(Concepts::Polynomial<LF>)
34 return order(lf);
35 else
36 return localFctOrder;
37 };
38
39 test_exit(Concepts::Polynomial<LF> || localFctOrder >= 0, "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()`.");
40
41 using Range = typename decltype(localFct)::Range;
42 Range result(0);
43
44 using Rules = Dune::QuadratureRules<typename GV::ctype, GV::dimension>;
45 for (auto const& element : elements(gridView, Dune::Partitions::interior)) {
46 auto geometry = element.geometry();
47
48 localFct.bind(element);
49 auto const& quad = Rules::rule(element.type(), quadOrder(localFct));
50
51 for (auto const qp : quad)
52 result += localFct(qp.position()) * geometry.integrationElement(qp.position()) * qp.weight();
53 localFct.unbind();
54 }
55
56 return gridView.comm().sum(result);
57 }
58
59} // end namespace AMDiS
decltype(auto) makeGridFunction(PreGridFct const &preGridFct, GridView const &gridView)
Generator for Gridfunctions from Expressions (PreGridfunctions)
Definition: GridFunction.hpp:168