dune-fem 2.8.0
Loading...
Searching...
No Matches
marking/default.hh
Go to the documentation of this file.
1#ifndef DUNE_FEMPY_PY_DEFAULT_MARKING_HH
2#define DUNE_FEMPY_PY_DEFAULT_MARKING_HH
3
4#include <array>
5#include <memory>
6#include <dune/geometry/referenceelements.hh>
9
10namespace Dune
11{
12 namespace Fem
13 {
14 namespace GridAdaptation
15 {
16 struct Marker
17 {
18 static const int refine = 1;
19 static const int keep = 0;
20 static const int coarsen = -1;
21 };
22
23 // Indicator is of type DiscreteFunction
24 template <class Grid, class Indicator>
25 static inline
26 std::pair< int, int >
27 mark( Grid& grid, Indicator& indicator,
28 const double refineTolerance, const double coarsenTolerance,
29 int minLevel = 0, int maxLevel = -1,
30 const bool markNeighbors = false )
31 {
32 Dune::Fem::ConstLocalFunction<Indicator> localIndicator(indicator);
33 typename Indicator::RangeType value;
34
35 std::array< int, 2 > sumMarks = {{ 0,0 }};
36
37 int& refMarked = sumMarks[ 0 ];
38 int& crsMarked = sumMarks[ 1 ];
39
40 if ( maxLevel < 0 )
41 maxLevel = std::numeric_limits<int>::max();
42
43 const auto& gridPart = indicator.space().gridPart();
44
45 for (const auto &e : indicator.space())
46 {
47 if (!e.isLeaf()) continue;
48
49 const auto &gridEntity = Dune::Fem::gridEntity(e);
50 int marked = grid.getMark(gridEntity);
51 if (marked==1) continue;
52
53 localIndicator.bind(e);
54 const auto &center = Dune::referenceElement< typename Grid::ctype, Grid::dimension>( e.type() ).position(0,0);
55 localIndicator.evaluate(center,value);
56 double eta = std::abs(value[0]);
57 const int level = e.level();
58 if (eta>refineTolerance && level<maxLevel)
59 {
60 refMarked += grid.mark( Marker::refine, gridEntity);
61 if( markNeighbors )
62 {
63 const auto iEnd = gridPart.iend( e );
64 for( auto it = gridPart.ibegin( e ); it != iEnd; ++it )
65 {
66 const auto& intersection = *it;
67 if( intersection.neighbor() )
68 {
69 const auto& outside = intersection.outside();
70 if (outside.level()<maxLevel)
71 refMarked += grid.mark( Marker::refine, Dune::Fem::gridEntity(outside) );
72 }
73 }
74
75 }
76 }
77 else if (eta<coarsenTolerance && level>minLevel)
78 crsMarked += grid.mark( Marker::coarsen, gridEntity);
79 else
80 grid.mark(Marker::keep, gridEntity);
81 }
82
83 grid.comm().sum( sumMarks.data(), 2 );
84 return std::make_pair( refMarked, crsMarked );
85 } // end mark
86
87 } // end namespace GridAdaptation
88
89 } // end namespace Fem
90} // end namespace Dune
91#endif
Dune::Fem::Double abs(const Dune::Fem::Double &a)
Definition: double.hh:942
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 std::pair< int, int > mark(Grid &grid, Indicator &indicator, const double refineTolerance, const double coarsenTolerance, int minLevel=0, int maxLevel=-1, const bool markNeighbors=false)
Definition: marking/default.hh:27
Definition: marking/default.hh:17
static const int coarsen
Definition: marking/default.hh:20
static const int keep
Definition: marking/default.hh:19
static const int refine
Definition: marking/default.hh:18