diff --git a/dune/fufem/estimators/geometricmarking.hh b/dune/fufem/estimators/geometricmarking.hh index 3c2fc9c68fd6e673725b5dceb1aadc23c552a6f7..d26b07b8bf6940942415368f19e46f66d54f48cc 100644 --- a/dune/fufem/estimators/geometricmarking.hh +++ b/dune/fufem/estimators/geometricmarking.hh @@ -2,7 +2,6 @@ #define DUNE_GEOMETRIC_ERROR_ESTIMATOR_HH #include <dune/common/fvector.hh> -#include <dune/istl/bvector.hh> /** \brief Error estimator which marks an element for refinement if it fulfills a geometric condition */ @@ -13,64 +12,49 @@ class GeometricEstimator { enum {dimworld = GridType::dimensionworld}; public: + static void estimate(GridType& grid, + bool (*refinementCondition)(const Dune::FieldVector<double, dim>& pos)) { - void estimate(GridType& grid, - bool (*refinementCondition)(const Dune::FieldVector<double, GridType::dimension>& pos) ) { - - typedef typename GridType::template Codim<0>::LeafIterator LeafIterator; - - LeafIterator lIt = grid.template leafbegin<0>(); - LeafIterator lEndIt = grid.template leafend<0>(); - - for (; lIt!=lEndIt; ++lIt) { + for (const auto& e : elements(grid.leafGridView())) { // Compute element's center of mass - Dune::FieldVector<double, GridType::dimension> center(0); - int numCorners = lIt->template count<GridType::dimension>(); + const auto& geometry = e.geometry(); + auto center = geometry.corner(0); - for (int i=0; i<numCorners; i++) - center += lIt->geometry().corner(i); + for (int i=1; i<geometry.corners(); i++) + center += geometry.corner(i); - center /= (double)numCorners; + center /= (double)geometry.corners(); // Mark if the condition is fulfilled - if ( (*refinementCondition)(center) ) - grid.mark(1, *lIt); - + if ((*refinementCondition)(center) ) + grid.mark(1, e); } - } - void estimate(GridType& grid, + template <class BlockVectorType> + static void estimate(GridType& grid, bool (*refinementCondition)(const Dune::FieldVector<double, dimworld>& pos), - const Dune::BlockVector<Dune::FieldVector<double, dimworld> >& deformation) { + const BlockVectorType& deformation) { - const typename GridType::Traits::LeafIndexSet& indexSet = grid.leafIndexSet(); - typedef typename GridType::template Codim<0>::LeafIterator LeafIterator; + const auto& indexSet = grid.leafIndexSet(); - LeafIterator lIt = grid.template leafbegin<0>(); - LeafIterator lEndIt = grid.template leafend<0>(); - - for (; lIt!=lEndIt; ++lIt) { + for (const auto& e : elements(grid.leafGridView())) { // Compute element's center of mass Dune::FieldVector<double, dimworld> center(0); - int numCorners = lIt->template count<GridType::dimension>(); + const auto& geometry = e.geometry(); - for (int i=0; i<numCorners; i++) { - center += lIt->geometry()[i] + deformation[indexSet.template subIndex<dim>(*lIt,i)]; - } + for (int i=0; i<geometry.corners(); i++) + center += geometry.corner(i) + deformation[indexSet.template subIndex<dim>(e,i)]; - center /= numCorners; + center /= (double) geometry.corners(); // Mark if the condition is fulfilled if ( (*refinementCondition)(center) > 0) - grid.mark(1, lIt); - + grid.mark(1, e); } - } - }; #endif