dune-fem 2.8.0
Loading...
Searching...
No Matches
maximum.hh
Go to the documentation of this file.
1#ifndef DUNE_FEM_MARKING_MAXIMUM_HH
2#define DUNE_FEM_MARKING_MAXIMUM_HH
3
4#include <cstddef>
5
6#include <algorithm>
7
8#include <dune/geometry/referenceelements.hh>
9
10#include <dune/grid/common/rangegenerators.hh>
13
15
16namespace Dune
17{
18
19 namespace Fem
20 {
21
22 // maximalMarking
23 // --------------
24
25 template <class Grid, class Indicator>
26 static inline
27 std::pair<int, int>
28 maximumMarking( Grid& grid, const Indicator& indicator,
29 double theta, int maxLevel = -1)
30 {
31 using Marker = GridAdaptation::Marker;
32
33 int refMarked = 0;
34 int crsMarked = 0;
35 maxLevel = ( maxLevel < 0 ) ? std::numeric_limits< int >::max() : maxLevel;
36 using std::max;
37
38 double maxError( 0 );
39 for( const auto &element : elements( grid.leafGridView() ) )
40 maxError = max( maxError, localError( element ) );
41 maxError = grid.comm().max( maxError );
42
43 ConstLocalFunction< Indicator > localIndicator( indicator );
44 typename Indicator :: RangeType value;
45
46 for (const auto &e : indicator.space())
47 {
48 if (!e.isLeaf()) continue;
49 const auto &gridEntity = Dune::Fem::gridEntity(e);
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 )
55 {
56 if (e.level()<maxLevel)
57 refMarked += grid.mark(Marker::refine, gridEntity);
58 else
59 grid.mark(Marker::keep, gridEntity);
60 }
61 }
62 return std::make_pair( grid.comm().sum(refMarked), grid.comm().sum(crsMarked) );
63 }
64
65 } // namespace Fem
66
67} // namespace Dune
68
69#endif // #ifndef DUNE_FEM_MARKING_MAXIMUM_HH
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