4 #ifndef DUNE_PDELAB_COMMON_FUNCTIONUTILITIES_HH 5 #define DUNE_PDELAB_COMMON_FUNCTIONUTILITIES_HH 11 #include <dune/common/debugstream.hh> 13 #include <dune/geometry/quadraturerules.hh> 14 #include <dune/geometry/type.hh> 16 #include <dune/grid/common/gridenums.hh> 17 #include <dune/grid/utility/hierarchicsearch.hh> 51 typename GF::Traits::RangeType& sum,
52 unsigned qorder = 1) {
53 typedef typename GF::Traits::GridViewType GV;
54 typedef typename GV::template Codim<0>::Geometry Geometry;
55 typedef typename GF::Traits::RangeType Range;
56 typedef typename GF::Traits::DomainFieldType DF;
57 static const int dimD = GF::Traits::dimDomain;
58 typedef Dune::QuadratureRule<DF,dimD> QR;
59 typedef Dune::QuadratureRules<DF,dimD> QRs;
63 for(
const auto& cell : elements(gf.getGridView(),Dune::Partitions::interior)) {
64 const Geometry& geo = cell.geometry();
65 Dune::GeometryType gt = geo.type();
66 const QR& rule = QRs::rule(gt,qorder);
67 for (
const auto& qip : rule) {
69 gf.evaluate(cell,qip.position(),val);
72 val *= qip.weight() * geo.integrationElement(qip.position());
90 typedef typename GF::Traits::GridViewType GV;
91 typedef typename GV::template Codim<0>::Entity Entity;
92 typedef typename GF::Traits::DomainType Domain;
93 typedef typename GF::Traits::RangeType Range;
106 template<
class GFHandle>
111 evalRank = gfp->getGridView().comm().size();
112 int myRank = gfp->getGridView().comm().rank();
115 (HierarchicSearch<typename GV::Grid, typename GV::IndexSet>
116 (gfp->getGridView().grid(), gfp->getGridView().indexSet()).
117 template findEntity<Interior_Partition>(xg)));
119 if(e->partitionType() == InteriorEntity)
122 catch(
const Dune::GridError&) { }
123 evalRank = gfp->getGridView().comm().min(evalRank);
124 if(myRank == evalRank)
125 xl = e->geometry().local(xg);
128 if(myRank == 0 && evalRank == gfp->getGridView().comm().size())
129 dwarn <<
"Warning: GridFunctionProbe at (" << xg <<
") is outside " 130 <<
"the grid" << std::endl;
178 typedef typename GF::Traits::RangeFieldType RF;
179 if(evalRank == gfp->getGridView().comm().size())
180 val = std::numeric_limits<RF>::quiet_NaN();
182 if(gfp->getGridView().comm().rank() == evalRank)
183 gfp->evaluate(*e, xl, val);
184 gfp->getGridView().comm().broadcast(&val,1,evalRank);
201 void eval(Range& val,
int rank = 0)
const {
206 std::shared_ptr<const GF> gfsp;
208 std::shared_ptr<Entity> e;
218 #endif // DUNE_PDELAB_COMMON_FUNCTIONUTILITIES_HH void setGridFunction(const GF &gf)
Set a new GridFunction.
Definition: functionutilities.hh:139
For backward compatibility – Do not use this!
Definition: adaptivity.hh:27
void setGridFunction(const GF *gf)
Set a new GridFunction.
Definition: functionutilities.hh:151
GridFunctionProbe(const GFHandle &gf, const Domain &xg)
Constructor.
Definition: functionutilities.hh:107
void setGridFunction(const std::shared_ptr< const GF > &gf)
Set a new GridFunction.
Definition: functionutilities.hh:165
void integrateGridFunction(const GF &gf, typename GF::Traits::RangeType &sum, unsigned qorder=1)
Integrate a GridFunction.
Definition: functionutilities.hh:50
Evaluate a GridFunction at a certain global coordinate.
Definition: functionutilities.hh:89
XG & xg
Definition: interpolate.hh:67
void eval(Range &val, int rank=0) const
evaluate the GridFunction and communicate result to the given rank
Definition: functionutilities.hh:201
void eval_all(Range &val) const
evaluate the GridFunction and broadcast result to all ranks
Definition: functionutilities.hh:177