1#ifndef DUNE_FEM_MARKING_MAXIMUM_HH
2#define DUNE_FEM_MARKING_MAXIMUM_HH
8#include <dune/geometry/referenceelements.hh>
10#include <dune/grid/common/rangegenerators.hh>
25 template <
class Gr
id,
class Indicator>
29 double theta,
int maxLevel = -1)
35 maxLevel = ( maxLevel < 0 ) ? std::numeric_limits< int >::max() : maxLevel;
39 for(
const auto &element : elements( grid.leafGridView() ) )
40 maxError =
max( maxError, localError( element ) );
41 maxError = grid.comm().max( maxError );
44 typename Indicator :: RangeType value;
46 for (
const auto &e : indicator.space())
48 if (!e.isLeaf())
continue;
50 localIndicator.bind(e);
51 const auto& center = Dune::referenceElement< typename Grid::ctype, Grid::dimension>( e.type() ).position(0,0);
52 localIndicator.evaluate(center,value);
53 double eta = value[0];
54 if( eta > theta*maxError )
56 if (e.level()<maxLevel)
57 refMarked += grid.mark(Marker::refine,
gridEntity);
62 return std::make_pair( grid.comm().sum(refMarked), grid.comm().sum(crsMarked) );
double max(const Dune::Fem::Double &v, const double p)
Definition: double.hh:965
Definition: bindguard.hh:11
const GridEntityAccess< Entity >::GridEntityType & gridEntity(const Entity &entity)
Definition: gridpart.hh:412
typename Impl::ConstLocalFunction< GridFunction >::Type ConstLocalFunction
Definition: const.hh:604
static double max(const Double &v, const double p)
Definition: double.hh:398
static std::pair< int, int > maximumMarking(Grid &grid, const Indicator &indicator, double theta, int maxLevel=-1)
Definition: maximum.hh:28
Definition: marking/default.hh:17