diff --git a/dune/fufem/boundarypatchprolongator.hh b/dune/fufem/boundarypatchprolongator.hh
index cb3a9c362b1b8ffb7e74274e76aec934363bf7c8..6be766c53f892059e327d1011999f70df26ca316 100644
--- a/dune/fufem/boundarypatchprolongator.hh
+++ b/dune/fufem/boundarypatchprolongator.hh
@@ -57,8 +57,7 @@ public:
     template<class BoundaryPatchPointerType>
     static void prolong(std::vector<BoundaryPatchPointerType>& patches)
     {
-        typedef typename GridType::template Codim<0>::Entity::EntityPointer ElementPointer;
-        typedef typename BoundaryPatch<typename GridType::LevelGridView>::iterator PatchIterator;
+        typedef typename GridType::template Codim<0>::Entity::Entity Element;
 
         if (not(patches[0]->isInitialized()))
             DUNE_THROW(Dune::Exception, "Coarse boundary patch has not been set up correctly!");
@@ -70,15 +69,13 @@ public:
         for (int i=1; i<=maxLevel; i++)
             patches[i]->setup(grid.levelGridView(i));
 
-        PatchIterator pIt = patches[0]->begin();
-        PatchIterator pEnd = patches[0]->end();
-        for (; pIt!=pEnd; ++pIt)
+        for (const auto& pIt : patches[0])
         {
-            ElementPointer inside = pIt->inside();
+            const Element& inside = pIt.inside();
 
-            if (not(inside->isLeaf()) and (inside->level()<maxLevel))
+            if (not(inside.isLeaf()) and (inside.level()<maxLevel))
             {
-                Face<GridType> face(grid, *inside, pIt->indexInInside());
+                Face<GridType> face(grid, inside, pIt.indexInInside());
 
                 typename Face<GridType>::HierarchicIterator it = face.hbegin(maxLevel);
                 typename Face<GridType>::HierarchicIterator end = face.hend(maxLevel);
@@ -102,7 +99,6 @@ public:
     static void prolong(const BoundaryPatch<typename GridType::LevelGridView>& coarsePatch,
                         BoundaryPatch<typename GridType::LeafGridView>& finePatch)
     {
-        typedef typename BoundaryPatch<typename GridType::LevelGridView>::iterator PatchIterator;
 
         if (not(coarsePatch.isInitialized()))
             DUNE_THROW(Dune::Exception, "Coarse boundary patch has not been set up correctly!");
@@ -112,16 +108,15 @@ public:
         // Set array sizes correctly
         finePatch.setup(grid.leafGridView());
 
-        PatchIterator pIt = coarsePatch.begin();
-        PatchIterator pEnd = coarsePatch.end();
-        for (; pIt!=pEnd; ++pIt)
+        for (const auto& pIt : coarsePatch)
         {
-            const auto inside = pIt->inside();
+            const auto& inside = pIt.inside();
+
             if (inside.isLeaf())
-                finePatch.addFace(inside, pIt->indexInInside());
+                finePatch.addFace(inside, pIt.indexInInside());
             else
             {
-                Face<GridType> face(grid, inside, pIt->indexInInside());
+                Face<GridType> face(grid, inside, pIt.indexInInside());
 
                 typename Face<GridType>::HierarchicIterator it = face.hbegin(grid.maxLevel());
                 typename Face<GridType>::HierarchicIterator end = face.hend(grid.maxLevel());
diff --git a/dune/fufem/differencenormsquared.hh b/dune/fufem/differencenormsquared.hh
index 9c39297589a4ff1239ce4bd6b20ab7553d15fa4b..9d62b0af67e87f9d6f132c4996582b92d03c748b 100644
--- a/dune/fufem/differencenormsquared.hh
+++ b/dune/fufem/differencenormsquared.hh
@@ -97,34 +97,35 @@ public:
         // handle each vertex only once
         Dune::BitSetVector<1> handled(targetGrid.size(dim), false);
 
-        typename GridType::template Codim<0>::LeafIterator eIt    = targetGrid.template leafbegin<0>();
-        typename GridType::template Codim<0>::LeafIterator eEndIt = targetGrid.template leafend<0>();
+        typename GridType::LeafGridView targetLeafView = targetGrid.leafGridView();
+        const auto& indexSet = targetLeafView.indexSet();
 
         // Loop over all vertices by looping over all elements
-        for (; eIt!=eEndIt; ++eIt) {
+        for (const auto& e : elements(targetLeafView)) {
 
-            for (int i=0; i<eIt->subEntities(dim); i++) {
+            for (int i=0; i<e.subEntities(dim); i++) {
 
-                typename GridType::template Codim<0>::EntityPointer element(eIt);
-                FieldVector<double,dim> pos = ReferenceElements<double, dim>::general(eIt->type()).position(i,dim);
+                typename GridType::template Codim<0>::Entity element(e);
+                FieldVector<double,dim> pos = ReferenceElements<double, dim>::general(e.type()).position(i,dim);
 
-                if (handled[targetGrid.leafIndexSet().subIndex(*eIt,i,dim)][0])
+                if (handled[indexSet.subIndex(e,i,dim)][0])
                     continue;
 
-                handled[targetGrid.leafIndexSet().subIndex(*eIt,i,dim)] = true;
+                handled[indexSet.subIndex(e,i,dim)] = true;
 
-                assert(checkInside(element->type(), pos, 1e-7));
+                assert(checkInside(element.type(), pos, 1e-7));
 
                 // ////////////////////////////////////////////////////////////////////
                 // Get an element on the coarsest grid which contains the vertex and
                 // its local coordinates there
                 // ////////////////////////////////////////////////////////////////////
-                while (element->level() != 0){
 
-                    pos = element->geometryInFather().global(pos);
-                    element = element->father();
+                while (element.level() != 0){
 
-                    assert(checkInside(element->type(), pos, 1e-7));
+                    pos = element.geometryInFather().global(pos);
+                    element = element.father();
+
+                    assert(checkInside(element.type(), pos, 1e-7));
                 }
 
                 // ////////////////////////////////////////////////////////////////////
@@ -133,30 +134,29 @@ public:
                 // ////////////////////////////////////////////////////////////////////
                 LevelMultipleCodimMultipleGeomTypeMapper<GridType,Dune::MCMGElementLayout> uniformP0Mapper (targetGrid,  0);
                 LevelMultipleCodimMultipleGeomTypeMapper<GridType,Dune::MCMGElementLayout> adaptiveP0Mapper(sourceGrid, 0);
-                int coarseIndex = uniformP0Mapper.map(*element);
+                int coarseIndex = uniformP0Mapper.map(element);
 
-                typename GridType::template Codim<0>::LevelIterator adaptEIt    = sourceGrid.template lbegin<0>(0);
-                typename GridType::template Codim<0>::LevelIterator adaptEEndIt = sourceGrid.template lend<0>(0);
+                typename GridType::LevelGridView sourceLevelView = sourceGrid.levelView(0);
 
-                for (; adaptEIt!=adaptEEndIt; ++adaptEIt)
-                    if (adaptiveP0Mapper.map(*adaptEIt) == coarseIndex)
+                for (auto&& adaptE : elements(sourceLevelView))
+                    if (adaptiveP0Mapper.map(adaptE) == coarseIndex) {
+                        element = std::move(adaptE);
                         break;
-
-                element = adaptEIt;
+                    }
 
                 // ////////////////////////////////////////////////////////////////////////
                 //   Find a corresponding point (not necessarily vertex) on the leaf level
                 //   of the adaptive grid.
                 // ////////////////////////////////////////////////////////////////////////
 
-                while (!element->isLeaf()) {
+                while (!element.isLeaf()) {
 
 
-                    typename GridType::template Codim<0>::Entity::HierarchicIterator hIt    = element->hbegin(element->level()+1);
-                    typename GridType::template Codim<0>::Entity::HierarchicIterator hEndIt = element->hend(element->level()+1);
+                    typename GridType::template Codim<0>::Entity::HierarchicIterator hIt    = element.hbegin(element.level()+1);
+                    typename GridType::template Codim<0>::Entity::HierarchicIterator hEndIt = element.hend(element.level()+1);
 
                     FieldVector<double,dim> childPos;
-                    assert(checkInside(element->type(), pos, 1e-7));
+                    assert(checkInside(element.type(), pos, 1e-7));
 
                     for (; hIt!=hEndIt; ++hIt) {
 
@@ -167,7 +167,7 @@ public:
                     }
 
                     assert(hIt!=hEndIt);
-                    element = hIt;
+                    element = *hIt;
                     pos     = childPos;
 
                 }
@@ -175,8 +175,8 @@ public:
                 // ////////////////////////////////////////////////////////////////////////
                 //   Sample adaptive function
                 // ////////////////////////////////////////////////////////////////////////
-                sourceFunction.evaluateLocal(*element, pos,
-                                            target[targetGrid.leafIndexSet().subIndex(*eIt,i,dim)]);
+                sourceFunction.evaluateLocal(element, pos,
+                                            target[indexSet.subIndex(*eIt,i,dim)]);
 
             }
 
diff --git a/dune/fufem/facehierarchy.hh b/dune/fufem/facehierarchy.hh
index f7f059e97853dbca376ecce981b3ed57a1379350..ccd04309fdc173f84daf1bcc600c94e18ca67324 100644
--- a/dune/fufem/facehierarchy.hh
+++ b/dune/fufem/facehierarchy.hh
@@ -60,7 +60,7 @@ struct FaceHierarchy
 
         static int findFatherFace(const GridType& grid, const Element& element, const int face)
         {
-            const Element father = element.father();
+            const Element& father = element.father();
             for(unsigned int i=0; i < father->subEntities(1); ++i)
                 if (isSubFace(father, i, element, face))
                     return i;
@@ -69,9 +69,10 @@ struct FaceHierarchy
 
         static LevelIntersectionIterator findFatherFaceIntersection(const GridType& grid, const Element& element, const int face)
         {
-            const Element father = element.father();
-            LevelIntersectionIterator nIt = grid.levelGridView(father->level()).ibegin(father);
-            LevelIntersectionIterator nEnd = grid.levelGridView(father->level()).iend(father);
+            const Element& father = element.father();
+
+            LevelIntersectionIterator nIt = grid.levelGridView(father.level()).ibegin(father);
+            LevelIntersectionIterator nEnd = grid.levelGridView(father.level()).iend(father);
             for(; nIt != nEnd; ++nIt)
             {
                 if (isSubFace(father, nIt->indexInInside(), element, face))
@@ -80,13 +81,7 @@ struct FaceHierarchy
             DUNE_THROW(Dune::Exception, "No father face found!");
         }
 
-        static
-#if DUNE_VERSION_NEWER(DUNE_GRID, 2, 4) || DOXYGEN
-        Element
-#else
-        ElementPointer
-#endif
-        finestCoarseNeighbor(const GridType& grid, const Element& element, int face)
+        static Element finestCoarseNeighbor(const GridType& grid, const Element& element, int face)
         {
             LevelIntersectionIterator nIt = grid.levelGridView(element.level()).ibegin(element);
             LevelIntersectionIterator nEnd = grid.levelGridView(element.level()).iend(element);
@@ -96,16 +91,16 @@ struct FaceHierarchy
                     return nIt->outside();
             }
 
-            ElementPointer checkElement(element);
-            while(checkElement->level()>0)
+            Element checkElement(element);
+            while(checkElement.level()>0)
             {
                 LevelIntersectionIterator fatherIntersection = findFatherFaceIntersection(grid, checkElement, face);
                 if (fatherIntersection->neighbor())
                     return fatherIntersection->outside();
-                checkElement = checkElement->father();
+                checkElement = checkElement.father();
                 face = fatherIntersection->indexInInside();
             }
-            return ElementPointer(element);
+            return element;
         }
 
 };
diff --git a/dune/fufem/improvegrid.hh b/dune/fufem/improvegrid.hh
index 8d9ae55fb9f4025cbde2a87b413455ae9f2114e0..3110fc6249c26e333a82a729efc08865f8bf2e92 100644
--- a/dune/fufem/improvegrid.hh
+++ b/dune/fufem/improvegrid.hh
@@ -93,7 +93,6 @@ void improveGrid(GridType& grid, double threshold)
     typedef typename GridType::template Codim<dim>::LevelIterator VertexIterator;
     typedef typename GridType::template Codim<0>::LevelIterator   ElementIterator;
     typedef typename GridType::template Codim<0>::Entity EntityType;
-    typedef typename EntityType::EntityPointer EntityPointer;
 
     VertexIterator vIt    = grid.template lbegin<dim>(maxLevel);
     VertexIterator vEndIt = grid.template lend<dim>(maxLevel);
@@ -127,7 +126,7 @@ void improveGrid(GridType& grid, double threshold)
 
     for (; eIt!=eEndIt; ++eIt) {
 
-        EntityPointer father = eIt->father();
+        const Entity& father = eIt->father();
 
         const typename GridType::template Codim<0>::LocalGeometry geometryInFather = eIt->geometryInFather();
         for (size_t i=0; i<eIt->subEntities(dim); i++) {
@@ -141,7 +140,7 @@ void improveGrid(GridType& grid, double threshold)
                      && std::abs(positionInFather[2]-edgeMidPoints[j][2]) < 0.1) {
 
                     isNewOnThisLevel[indexSet.subIndex(*eIt,i,dim)][0] = true;
-                    bisectionPosition[indexSet.subIndex(*eIt,i,dim)] = father->geometry().global(positionInFather);
+                    bisectionPosition[indexSet.subIndex(*eIt,i,dim)] = father.geometry().global(positionInFather);
                     break;
                 }
             }
diff --git a/dune/fufem/prolongboundarypatch.hh b/dune/fufem/prolongboundarypatch.hh
index 9dbbf05d3481b60f739f787818349d5f8a901b1c..1406e621baf6dab88c7e8abd864b4f968b0ec785 100644
--- a/dune/fufem/prolongboundarypatch.hh
+++ b/dune/fufem/prolongboundarypatch.hh
@@ -10,6 +10,7 @@
 #include <vector>
 
 #include <dune/grid/common/grid.hh>
+#include <dune/grid/common/rangegenerators.hh>
 #include <dune/fufem/boundarypatch.hh>
 #include <dune/fufem/facehierarchy.hh>
 
@@ -26,10 +27,7 @@ public:
     {
         using namespace Dune;
 
-        typedef typename GridType::template Codim<0>::LevelIterator ElementIterator;
         typedef typename GridType::template Codim<0>::Entity EntityType;
-        typedef typename EntityType::LevelIntersectionIterator NeighborIterator;
-        typedef typename EntityType::HierarchicIterator HierarchicIterator;
 
         const GridType& grid = patches[0].gridView().grid();
         if (!&grid)
@@ -41,21 +39,15 @@ public:
         for (int i=1; i<=grid.maxLevel(); i++)
             patches[i].setup(grid.levelGridView(i));
 
-
-        ElementIterator eIt    = grid.template lbegin<0>(0);
-        ElementIterator eEndIt = grid.template lend<0>(0);
-
         // Loop over the boundary on level 0
-        for (; eIt!=eEndIt; ++eIt) {
-
-            NeighborIterator nIt    = eIt->ilevelbegin();
-            NeighborIterator nEndIt = eIt->ilevelend();
+        typename GridType::LevelGridView levelView = grid.levelGridView(0);
+        for (const auto& e : elements(levelView)) {
 
             // Loop over all neighbors
-            for (; nIt!=nEndIt; ++nIt) {
+            for (const auto& i : intersections(levelView,e)) {
 
                 // if the element is a boundary element
-                if (patches[0].contains(nIt)) {
+                if (patches[0].contains(i)) {
 
                     // ///////////////////////////////////////////////////////////////
                     //   This implementation doesn't use the UGGrid-specific method
@@ -66,30 +58,25 @@ public:
                     const int dimworld = GridType::dimensionworld;
                     typedef typename GridType::ctype ctype;
                     FieldVector<ctype,dimworld-1> dummy;
-                    FieldVector<ctype,dimworld> level0Normal = nIt->unitOuterNormal(dummy);
+                    FieldVector<ctype,dimworld> level0Normal = i.unitOuterNormal(dummy);
 
                     // Find all boundary segments of the descendants of this element
-                    HierarchicIterator hIt    = eIt->hbegin(grid.maxLevel());
-                    HierarchicIterator hEndIt = eIt->hend(grid.maxLevel());
-
-                    for (; hIt!=hEndIt; ++hIt) {
+                    for (const auto& h : descendantElements(e,grid.maxLevel())) {
 
-                        if (hIt->level() == eIt->level())
+                        if (h.level() == e.level())
                             continue;
 
-                        NeighborIterator hnIt = hIt->ilevelbegin();
-                        NeighborIterator hnEndIt = hIt->ilevelend();
-
-                        for (; hnIt!=hnEndIt; ++hnIt) {
+                        typename GridType::LevelGridView hlevelView = grid.levelGridView(h.level());
+                        for (const auto& hi : intersections(hlevelView,h)) {
 
-                            if (hnIt->boundary()) {
+                            if (hi.boundary()) {
 
                                 // Test whether this boundary segment is contained in the
                                 // level 0 boundary segment by comparing the normals.
                                 /** \todo Isn't there a topological way to check this? */
                                 FieldVector<ctype,dim-1> dummy;
-                                if (0.95 < level0Normal*hnIt->unitOuterNormal(dummy))
-                                    patches[hIt->level()].addFace(hnIt);
+                                if (0.95 < level0Normal*hi.unitOuterNormal(dummy))
+                                    patches[h.level()].addFace(hi);
 
 
                             }
@@ -153,10 +140,7 @@ public:
     {
         using namespace Dune;
 
-        typedef typename GridType::template Codim<0>::LevelIterator ElementIterator;
         typedef typename GridType::template Codim<0>::Entity EntityType;
-        typedef typename GridType::template Codim<0>::Entity::EntityPointer EntityPointer;
-        typedef typename EntityType::LevelIntersectionIterator NeighborIterator;
 
         const GridType& grid = patches[0].gridView().grid();
         if (!&grid)
@@ -168,18 +152,12 @@ public:
         for (int i=1; i<=grid.maxLevel(); i++)
             patches[i].setup(grid.levelGridView(i));
 
-
-        ElementIterator eIt    = grid.template lbegin<0>(0);
-        ElementIterator eEndIt = grid.template lend<0>(0);
-
         // Loop over the boundary on level 0
-        for (; eIt!=eEndIt; ++eIt) {
-
-            NeighborIterator nIt    = grid.levelGridView(0).ibegin(*eIt);
-            NeighborIterator nEndIt = grid.levelGridView(0).iend(*eIt);
+        typename GridType::LevelGridView levelView = grid.levelGridView(0);
+        for (const auto& e : elements(levelView)) {
 
             // Loop over all neighbors
-            for (; nIt!=nEndIt; ++nIt) {
+            for (const auto& i : intersections(levelView,e)) {
 
                 // if the element is a boundary element
                 if (patches[0].contains(nIt)) {
@@ -189,14 +167,14 @@ public:
                     //   purely topologically.  For the other grids we need to compare
                     //   surface normals.
                     // /////////////////////////////////////////////////////////////////
-                    std::vector<EntityPointer> childElements(0,EntityPointer(grid.template lbegin<0>(0)));
+                    std::vector<EntityType> childElements(0,*levelView.template begin<0>());
                     std::vector<unsigned char> childElementSides;
-                    grid.getChildrenOfSubface(EntityPointer(eIt), nIt->indexInInside(), grid.maxLevel(),
+                    grid.getChildrenOfSubface(e, i.indexInInside(), grid.maxLevel(),
                                               childElements, childElementSides);
 
                     for (size_t i=0; i<childElements.size(); i++) {
-                        const int& level = childElements[i]->level();
-                        patches[level].addFace(*childElements[i], childElementSides[i]);
+                        const int& level = childElements[i].level();
+                        patches[level].addFace(childElements[i], childElementSides[i]);
                     }
 
                 }
@@ -252,10 +230,7 @@ public:
                         BoundaryPatch<typename GridType::LeafGridView>& finePatch)
 #endif
     {
-        typedef typename GridType::template Codim<0>::LevelIterator LevelElementIterator;
         typedef typename GridType::template Codim<0>::Entity EntityType;
-        typedef typename GridType::template Codim<0>::Entity::EntityPointer EntityPointer;
-        typedef typename EntityType::LevelIntersectionIterator NeighborIterator;
 
         const GridType& grid = coarsePatch.gridView().grid();
 
@@ -265,17 +240,12 @@ public:
         // Set array sizes correctly
         finePatch.setup(grid.leafGridView());
 
-        LevelElementIterator eIt    = grid.template lbegin<0>(0);
-        LevelElementIterator eEndIt = grid.template lend<0>(0);
-
         // Loop over the boundary on level 0
-        for (; eIt!=eEndIt; ++eIt) {
-
-            NeighborIterator nIt    = eIt->ilevelbegin();
-            NeighborIterator nEndIt = eIt->ilevelend();
+        typename GridType::LevelGridView levelView = grid.levelGridView(0);
+        for (const auto& e : elements(levelView)) {
 
             // Loop over all neighbors
-            for (; nIt!=nEndIt; ++nIt) {
+            for (const auto& i : intersections(levelView,e)) {
 
                 // if the element is a boundary element
                 if (coarsePatch.contains(nIt)) {
@@ -285,18 +255,18 @@ public:
                     //   purely topologically.  For the other grids we need to compare
                     //   surface normals.
                     // /////////////////////////////////////////////////////////////////
-                    if (eIt->isLeaf())
-                        finePatch.addFace(*eIt, nIt->indexInInside());
+                    if (e.isLeaf())
+                        finePatch.addFace(e, i.indexInInside());
 
-                    std::vector<EntityPointer> childElements(0, EntityPointer(grid.template lbegin<0>(0)));
+                    std::vector<EntityType> childElements(0, *levelView.template begin<0>());
                     std::vector<unsigned char> childElementSides;
-                    grid.getChildrenOfSubface(EntityPointer(eIt), nIt->indexInInside(), grid.maxLevel(),
+                    grid.getChildrenOfSubface(e, i.indexInInside(), grid.maxLevel(),
                                               childElements, childElementSides);
 
                     for (size_t i=0; i<childElements.size(); i++) {
 
-                        if (childElements[i]->isLeaf())
-                            finePatch.addFace(*childElements[i], childElementSides[i]);
+                        if (childElements[i].isLeaf())
+                            finePatch.addFace(childElements[i], childElementSides[i]);
 
                     }
 
diff --git a/dune/fufem/sampleonbitfield.hh b/dune/fufem/sampleonbitfield.hh
index 4d18d443af2802c3a905558804654f0661a55fc7..c2601d4cc4374df2b5737c1189f339938f5035f6 100644
--- a/dune/fufem/sampleonbitfield.hh
+++ b/dune/fufem/sampleonbitfield.hh
@@ -4,6 +4,7 @@
 #define DUNE_FUFEM_SAMPLE_ON_BITFIELD_HH
 
 #include <dune/grid/common/grid.hh>
+#include <dune/grid/common/rangegenerators.hh>
 #include <dune/common/bitsetvector.hh>
 #include <dune/geometry/referenceelements.hh>
 #include <dune/common/version.hh>
@@ -26,9 +27,7 @@ void sampleOnBitField(const GridType& grid,
     // We assume that VectorType is a block vector
     const int blocksize = VectorType::block_type::dimension;
 
-    typedef typename GridType::template Codim<0>::LevelIterator ElementIterator;
-    typedef typename GridType::template Codim<0>::Entity EntityType;
-    typedef typename EntityType::HierarchicIterator HierarchicIterator;
+    typedef typename GridType::template Codim<0>::Entity Element;
 
     const int dim = GridType::dimension;
 
@@ -47,40 +46,34 @@ void sampleOnBitField(const GridType& grid,
         dirichletValues[i] = 0;
     }
 
-    ElementIterator eIt    = grid.template lbegin<0>(0);
-    ElementIterator eEndIt = grid.template lend<0>(0);
-
     // Loop over the elements on level 0
-    for (; eIt!=eEndIt; ++eIt) {
+    for (const auto& e : elements(levelView)) {
 
         // Loop over all descendants of this element
-        HierarchicIterator hIt    = eIt->hbegin(grid.maxLevel());
-        HierarchicIterator hEndIt = eIt->hend(grid.maxLevel());
-
-        for (; hIt!=hEndIt; ++hIt) {
+        for (const auto& h : descendantElements(e,grid.maxLevel())) {
 
-            for (size_t i=0; i<hIt->subEntities(dim); i++) {
+            for (size_t i=0; i<h.subEntities(dim); i++) {
 
-                int fIdx = grid.levelIndexSet(hIt->level()).subIndex(*hIt, i, dim);
+                int fIdx = grid.levelIndexSet(h.level()).subIndex(h, i, dim);
 
                 // Logical position of the fine grid lagrange point
                 // in local coordinates of the coarse grid element
-                typename GridType::template Codim<0>::EntityPointer element(hIt);
-                Dune::FieldVector<double,dim> pos = Dune::ReferenceElements<double, dim>::general(hIt->type()).position(i,dim);
+                Element element(h);
+                Dune::FieldVector<double,dim> pos = Dune::ReferenceElements<double, dim>::general(h.type()).position(i,dim);
 
-                while (element->level() != 0){
+                while (element.level() != 0){
 
-                    pos = element->geometryInFather().global(pos);
-                    element = element->father();
+                    pos = element.geometryInFather().global(pos);
+                    element = element.father();
 
                 }
 
                 typename VectorType::value_type tmp;
-                coarseFunction.evaluateLocal(*eIt, pos, tmp);
+                coarseFunction.evaluateLocal(e, pos, tmp);
 
                 for (int j=0; j<blocksize; j++)
-                    if (dirichletNodes[hIt->level()][fIdx][j])
-                        dirichletValues[hIt->level()][fIdx][j] = tmp[j];
+                    if (dirichletNodes[h.level()][fIdx][j])
+                        dirichletValues[h.level()][fIdx][j] = tmp[j];
 
             }
 
@@ -107,8 +100,7 @@ void sampleOnBitField(const GridType& grid,
     // We assume that VectorType is a block vector
     const int blocksize = VectorType::block_type::dimension;
 
-    typedef typename GridType::template Codim<0>::LeafIterator ElementIterator;
-    typedef typename GridType::template Codim<0>::EntityPointer EntityPointer;
+    typedef typename GridType::template Codim<0>::Entity Element;
 
     const int dim = GridType::dimension;
 
@@ -124,30 +116,30 @@ void sampleOnBitField(const GridType& grid,
     fineDirichletValues.resize(grid.size(dim));
     fineDirichletValues = 0;
 
-    ElementIterator eIt    = grid.template leafbegin<0>();
-    ElementIterator eEndIt = grid.template leafend<0>();
+    typename GridType::LeafGridView leafView = grid.leafGridView();
+    const auto& indexSet = leafView.indexSet();
 
     // Loop over the leaf elements
-    for (; eIt!=eEndIt; ++eIt) {
+    for (const auto& e : elements(leafView)) {
 
-        for (size_t i=0; i<eIt->subEntities(dim); i++) {
+        for (size_t i=0; i<e.subEntities(dim); i++) {
 
-            int fIdx = grid.leafIndexSet().subIndex(*eIt, i, dim);
+            int fIdx = indexSet.subIndex(e, i, dim);
 
             // Logical position of the fine grid lagrange point
             // in local coordinates of the coarse grid element
-            typename GridType::template Codim<0>::EntityPointer element = EntityPointer(eIt);
-            Dune::FieldVector<double,dim> pos = Dune::ReferenceElements<double, dim>::general(eIt->type()).position(i,dim);
+            Element element(e);
+            Dune::FieldVector<double,dim> pos = Dune::ReferenceElements<double, dim>::general(e.type()).position(i,dim);
 
-            while (element->level() != 0){
+            while (element.level() != 0){
 
-                pos = element->geometryInFather().global(pos);
-                element = element->father();
+                pos = element.geometryInFather().global(pos);
+                element = element.father();
 
             }
 
             typename VectorType::value_type tmp;
-            coarseFunction.evaluateLocal(*element, pos, tmp);
+            coarseFunction.evaluateLocal(element, pos, tmp);
 
             for (int j=0; j<blocksize; j++)
                  if (dirichletNodes[eIt->level()][fIdx][j])
@@ -174,8 +166,7 @@ void sampleOnBitField(const GridType& grid,
 {
     // We assume that VectorType is a block vector
     const int blocksize = VectorType::block_type::dimension;
-
-    typedef typename GridType::template Codim<0>::LeafIterator ElementIterator;
+    typedef typename GridType::template Codim<0>::Entity Element;
 
     const int dim = GridType::dimension;
 
@@ -190,32 +181,30 @@ void sampleOnBitField(const GridType& grid,
     fineDirichletValues.resize(grid.size(dim));
     fineDirichletValues = 0;
 
-    ElementIterator eIt    = grid.template leafbegin<0>();
-    ElementIterator eEndIt = grid.template leafend<0>();
-
-    typedef typename GridType::template Codim<0>::EntityPointer EntityPointer;
+    typename GridType::LeafGridView leafView = grid.leafGridView();
+    const auto& indexSet = leafView.indexSet();
 
     // Loop over the leaf elements
-    for (; eIt!=eEndIt; ++eIt) {
+    for (const auto& e : elements(leafView)) {
 
-        for (size_t i=0; i<eIt->subEntities(dim); i++) {
+        for (size_t i=0; i<e.subEntities(dim); i++) {
 
-            int fIdx = grid.leafIndexSet().subIndex(*eIt, i, dim);
+            int fIdx = indexSet().subIndex(e, i, dim);
 
             // Logical position of the fine grid lagrange point
             // in local coordinates of the coarse grid element
-            typename GridType::template Codim<0>::EntityPointer element = EntityPointer(eIt);
-            Dune::FieldVector<double,dim> pos = Dune::ReferenceElements<double, dim>::general(eIt->type()).position(i,dim);
+            Element element(e);
+            Dune::FieldVector<double,dim> pos = Dune::ReferenceElements<double, dim>::general(e.type()).position(i,dim);
 
-            while (element->level() != 0){
+            while (element.level() != 0){
 
-                pos = element->geometryInFather().global(pos);
-                element = element->father();
+                pos = element.geometryInFather().global(pos);
+                element = element.father();
 
             }
 
             typename VectorType::value_type tmp;
-            coarseFunction.evaluateLocal(*element, pos, tmp);
+            coarseFunction.evaluateLocal(element, pos, tmp);
 
             for (int j=0; j<blocksize; j++)
                  if (fineDirichletNodes[fIdx][j])
diff --git a/dune/fufem/test/subgridxyfunctionalassemblertest.cc b/dune/fufem/test/subgridxyfunctionalassemblertest.cc
index 77bd907cd1834bdabce2786034bf8878252dfa29..501a1199289d6ee51d6117b183fb8e59de6d937d 100644
--- a/dune/fufem/test/subgridxyfunctionalassemblertest.cc
+++ b/dune/fufem/test/subgridxyfunctionalassemblertest.cc
@@ -43,7 +43,7 @@ void setupSubgridToHostgridTransfer(TransferOperatorType& matrix, const SubGridB
 {
     typedef typename SubGridBasis::GridView::Grid SubGridType;
     typedef typename SubGridType::HostGridType HostGridType;
-    typedef typename HostGridType::LeafGridView::template Codim<0>::EntityPointer HostElementPtr;
+    typedef typename HostGridType::LeafGridView::template Codim<0>::Entity HostElement;
     typedef typename SubGridType::LeafGridView::template Codim<0>::Iterator ElementIterator;
     typedef typename HostGridType::template Codim<0>::Entity::HierarchicIterator HierarchicIterator;
     typedef std::map<int, double> LinearCombination;
@@ -68,17 +68,17 @@ void setupSubgridToHostgridTransfer(TransferOperatorType& matrix, const SubGridB
     ElementIterator cEnd = subgrid.leafGridView().template end<0>();
     for(; cIt!=cEnd; ++cIt)
     {
-        HostElementPtr hostElement=subgrid.template getHostEntity<0>(*cIt);
+        HostElement hostElement = *subgrid.template getHostEntity<0>(*cIt);
 
         const CLFE& coarseFE = subgridbasis.getLocalFiniteElement(*cIt);
 
         // if element is leaf the transfer to the next level is locally the identity
-        if (hostElement->isLeaf())
+        if (hostElement.isLeaf())
         {
             for (size_t j=0; j<coarseFE.localBasis().size(); ++j)
             {
                 int coarseIndex = subgridbasis.index(*cIt, j);
-                int fineIndex = hostgridbasis.index(*hostElement, j);
+                int fineIndex = hostgridbasis.index(hostElement, j);
 
                 // visit each host node only once
                 if (processed[fineIndex][0])
@@ -100,8 +100,8 @@ void setupSubgridToHostgridTransfer(TransferOperatorType& matrix, const SubGridB
             std::vector<CFERange> valuesAtPosition(coarseFE.localBasis().size());
 
             // loop over all children on next level
-            HierarchicIterator fIt = hostElement->hbegin(level+1);
-            HierarchicIterator fEnd = hostElement->hend(level+1);
+            HierarchicIterator fIt = hostElement.hbegin(level+1);
+            HierarchicIterator fEnd = hostElement.hend(level+1);
             for (; fIt != fEnd; ++fIt)
             {
                 const FLFE& fineFE = hostgridbasis.getLocalFiniteElement(*fIt);
diff --git a/dune/fufem/utilities/hierarchicleafiterator.hh b/dune/fufem/utilities/hierarchicleafiterator.hh
index 05eff7df78f68522221f27940e6736c0b71e6ece..8304c9b02fcbe8c330f441a05c9e912ce732f2f0 100644
--- a/dune/fufem/utilities/hierarchicleafiterator.hh
+++ b/dune/fufem/utilities/hierarchicleafiterator.hh
@@ -21,7 +21,6 @@ class HierarchicLeafIterator :
     public Dune::ForwardIteratorFacade<HierarchicLeafIterator<GridImp>, const typename GridImp::template Codim<0>::Entity>
 {
     typedef typename GridImp::template Codim<0>::Entity Element;
-    typedef typename GridImp::template Codim<0>::EntityPointer EntityPointer;
     typedef typename GridImp::HierarchicIterator HierarchicIterator;
     enum {dim = GridImp::dimension};
     enum {dimworld = GridImp::dimensionworld};
@@ -37,19 +36,19 @@ public:
     {
 
         // if the element itself is leaf, then we don't have to iterate over the children
-        if (flag==begin && !element_->isLeaf()) {
-            hIt_ = element_->hbegin(maxlevel_);
+        if (flag==begin && !element_.isLeaf()) {
+            hIt_ = element_.begin(maxlevel_);
 
             //NOTE This class by now assumes that possible non-nestedness of the grid levels only arises
             // due to boundary parametrisation
             if (!nested_)  {
                 // Check if the element is a boundary element, and set the nested flag correspondingly
                 // If the element is not at the boundary, then we can neglect the nested flag
-                typename GridImp::LevelIntersectionIterator lIt = grid.levelGridView(element.level()).ibegin(element);
-                typename GridImp::LevelIntersectionIterator lEndIt = grid.levelGridView(element.level()).ibegin(element);
+                typename GridImp::LevelGridView levelView = grid.levelGridView(element_.level());
                 nested_ = true;
-                for (; lIt != lEndIt; ++lIt)
-                    if (lIt->boundary()) {
+
+                for (const auto lIt : intersections(levelView,element_)) {
+                    if (lIt.boundary()) {
                         nested_ = false;
                         break;
                     }
@@ -57,7 +56,7 @@ public:
 
             if (!hIt_->isLeaf())
                 increment();
-        }
+           }
     }
 
     //! Copy constructor
@@ -89,20 +88,22 @@ public:
 
         if (hIt_ == hEndIt)
             flag_ = end;
+        else
+            leafChild_ = *hIt_;
     }
 
     //! Dereferencing
     const Element& dereference() const {
         if (flag_ == end)
             DUNE_THROW(Dune::Exception,"HierarchicLeafIterator: Cannot dereference end iterator!");
-        return (element_->isLeaf()) ? *element_ : *hIt_;
+        return (element_.isLeaf()) ? element_ : leafChild_;
     }
 
     //! Compute the local geometry mapping from the leaf child to the starting element
     LocalGeometry geometryInAncestor() {
 
         //NOTE: We assume here that the type of the child and the ancestors are the same!
-        const Dune::ReferenceElement<typename GridImp::ctype, dim>& ref = Dune::ReferenceElements<typename GridImp::ctype,dim>::general(element_->type());
+        const Dune::ReferenceElement<typename GridImp::ctype, dim>& ref = Dune::ReferenceElements<typename GridImp::ctype,dim>::general(element_.type());
 
         // store local coordinates of the leaf child element within the coarse starting element
         std::vector<Dune::FieldVector<typename GridImp::ctype,dim> > corners(ref.size(dim));
@@ -110,22 +111,22 @@ public:
             corners[i] = ref.position(i,dim);
 
         // If the element is leaf, then return an Identity mapping
-        if (element_->isLeaf())
+        if (element_.isLeaf())
             return LocalGeometry(ref,corners);
 
         if (nested_) {
-            const typename Element::Geometry leafGeom = hIt_->geometry();
-            const typename Element::Geometry coarseGeom = element_->geometry();
+            const typename Element::Geometry leafGeom = leafChild_->geometry();
+            const typename Element::Geometry coarseGeom = element_.geometry();
 
             for (int i=0; i<corners.size();i++)
                 corners[i] = coarseGeom.local(leafGeom.corner(i));
             return LocalGeometry(ref,corners);
         }
 
-        EntityPointer father(*hIt_);
+        Element father(leafChild_);
         while (father != element_) {
 
-            const typename Element::LocalGeometry fatherGeom = hIt_->geometryInFather();
+            const typename Element::LocalGeometry fatherGeom = father.geometryInFather();
             father = father->father();
 
             for (int i=0; i<corners.size(); i++)
@@ -138,7 +139,7 @@ public:
 private:
 
     //! The starting element
-    const EntityPointer element_;
+    const Element& element_;
 
     //! The grids maxlevel
     int maxlevel_;
@@ -146,6 +147,9 @@ private:
     //! The actual hierarchic iterator
     HierarchicIterator hIt_;
 
+    //! The actual hierachic leaf element
+    Element leafChild_;
+
     //! Position flag
     PositionFlag flag_;