mmiterator.hh 9.43 KB
Newer Older
Praetorius, Simon's avatar
Praetorius, Simon committed
1 2 3 4 5 6 7 8 9
// -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
// vi: set et ts=4 sw=2 sts=2:
#ifndef DUNE_MULTIMESH_ITERATOR_HH
#define DUNE_MULTIMESH_ITERATOR_HH

#include <numeric>
#include <stack>

#include <dune/common/iteratorfacades.hh>
10
#include <dune/common/std/type_traits.hh>
Praetorius, Simon's avatar
Praetorius, Simon committed
11 12 13
#include <dune/grid/common/exceptions.hh>
#include <dune/grid/common/gridenums.hh>

Praetorius, Simon's avatar
Praetorius, Simon committed
14
#include "mmentity.hh"
Praetorius, Simon's avatar
Praetorius, Simon committed
15 16 17

namespace Dune
{
18 19 20 21 22
  namespace tag
  {
    struct with_gridview {};
  }

Praetorius, Simon's avatar
Praetorius, Simon committed
23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46
  /** \brief Iterator over all entities of a given codimension and level of a grid.
   *  \ingroup MultiMesh
   */
  template <int codim, PartitionIteratorType pitype, class GridImp>
  class MultiMeshIterator
  {
  public:
    template <class... Args>
    MultiMeshIterator(Args&&... args)
    {
      DUNE_THROW(NotImplemented, "Not implemented for entity of codim != 0!");
    }

    int operator*() const { return 0; }
    MultiMeshIterator& operator++() { return *this; }
    MultiMeshIterator operator++(int) { return *this; }
    bool operator!=(MultiMeshIterator const&) const { return false; }
  };


  // Implemented for codim 0 entities only
  template <PartitionIteratorType pitype, class GridImp>
  class MultiMeshIterator<0, pitype, GridImp>
      : public ForwardIteratorFacade<MultiMeshIterator<0,pitype,GridImp>,
Praetorius, Simon's avatar
Praetorius, Simon committed
47 48
            MultiEntity<GridImp>,
            MultiEntity<GridImp> const&>
Praetorius, Simon's avatar
Praetorius, Simon committed
49 50 51 52 53 54 55 56
  {
  private:
    // LevelIterator to the equivalent entity in the host grid
    using HostGridLevelIterator =
      typename GridImp::HostGridType::template Codim<0>::template Partition<pitype>::LevelIterator;

    using HostEntity = typename GridImp::HostGridType::template Codim<0>::Entity;

57 58
    using EntityTest = std::function<bool(HostEntity)>;

Praetorius, Simon's avatar
Praetorius, Simon committed
59 60
    struct EntityStackEntry
    {
61
      template <class Entity>
Praetorius, Simon's avatar
Praetorius, Simon committed
62
      explicit EntityStackEntry (Entity&& entity)
63 64
        : it(entity.hbegin(entity.level()+1))
        , end(entity.hend(entity.level()+1))
Praetorius, Simon's avatar
Praetorius, Simon committed
65 66 67 68 69 70 71 72 73 74 75 76
      {}

      typename HostEntity::HierarchicIterator it;
      typename HostEntity::HierarchicIterator end;
    };

    class EntityStack
        : public std::stack<EntityStackEntry, std::vector<EntityStackEntry>>
    {
      using Super = std::stack<EntityStackEntry, std::vector<EntityStackEntry>>;

    public:
77
      explicit EntityStack (int maxLevel)
Praetorius, Simon's avatar
Praetorius, Simon committed
78
      {
79
        c.reserve(maxLevel);
Praetorius, Simon's avatar
Praetorius, Simon committed
80 81
      }

82
      // return true if all levels >= l are finished, i.e. hierarchic iterators it+1 == end
Praetorius, Simon's avatar
Praetorius, Simon committed
83
      bool finished (std::size_t l = 0) const
Praetorius, Simon's avatar
Praetorius, Simon committed
84 85
      {
        bool f = true;
86 87
        for (auto tmp = c[l].it; f && l < c.size(); ++l) {
          tmp = c[l].it;
Praetorius, Simon's avatar
Praetorius, Simon committed
88
          f = f && (++tmp) == c[l].end;
89
        }
Praetorius, Simon's avatar
Praetorius, Simon committed
90 91 92 93 94 95 96 97 98
        return f;
      }

    protected:
      using Super::c;
    };

  public:
    /// Constructor. Stores a pointer to the grid
99 100 101 102
    MultiMeshIterator (const GridImp* multiMesh, int level)
      : incrementAllowed_(multiMesh->size(), true)
      , contains_(multiMesh->size(), EntityTest{[level](const HostEntity& entity) { return entity.level() == level; }})
      , maxLevel_(multiMesh->size(), level)
103
      , multiEntity_(multiMesh->size())
Praetorius, Simon's avatar
Praetorius, Simon committed
104
    {
105
      for (auto const& grid : multiMesh->grids_) {
Praetorius, Simon's avatar
Praetorius, Simon committed
106 107 108 109 110
        macroIterators_.push_back(grid->levelGridView(0).template begin<0,pitype>());
        macroEndIterators_.push_back(grid->levelGridView(0).template end<0,pitype>());
      }

      // go to first leaf entity on all grids
111 112 113 114
      entityStacks_.reserve(multiMesh->size());
      for (std::size_t i = 0; i < multiMesh->size(); ++i) {
        maxLevel_[i] = (*multiMesh)[i].maxLevel();
        entityStacks_.emplace_back((*multiMesh)[i].maxLevel());
Praetorius, Simon's avatar
Praetorius, Simon committed
115
        initialIncrement(i);
116
      }
Praetorius, Simon's avatar
Praetorius, Simon committed
117 118 119 120 121 122 123 124 125
    }

    /// Constructor which create the end iterator
    /**
     * \param  multiMesh  Pointer to grid instance
     * \param  endDummy   Here only to distinguish it from the other constructor
     */
    MultiMeshIterator (const GridImp* multiMesh, int level, bool endDummy)
    {
126
      for (auto const& grid : multiMesh->grids_)
Praetorius, Simon's avatar
Praetorius, Simon committed
127 128 129
        macroIterators_.push_back(grid->levelGridView(0).template end<0,pitype>());
    }

Praetorius, Simon's avatar
Praetorius, Simon committed
130
    /// Constructor which create the leaf-iterator
Praetorius, Simon's avatar
Praetorius, Simon committed
131 132 133 134
    MultiMeshIterator (const GridImp* multiMesh)
      : MultiMeshIterator{multiMesh, -1}
    {}

Praetorius, Simon's avatar
Praetorius, Simon committed
135
    /// Constructor which create the end leaf-iterator
Praetorius, Simon's avatar
Praetorius, Simon committed
136 137 138 139
    MultiMeshIterator (const GridImp* multiMesh, bool endDummy)
      : MultiMeshIterator{multiMesh, -1, endDummy}
    {}

140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163
    template <class... GridViews,
      std::enable_if_t<not Std::disjunction<std::is_same<GridViews,bool>...>::value, int> = 0>
    MultiMeshIterator (tag::with_gridview, GridViews const&... gridViews)
      : incrementAllowed_(sizeof...(GridViews), true)
      , maxLevel_{gridViews.grid().maxLevel()...}
      , multiEntity_(sizeof...(GridViews))
      , macroIterators_{gridViews.grid().levelGridView(0).template begin<0,pitype>()...}
      , macroEndIterators_{gridViews.grid().levelGridView(0).template end<0,pitype>()...}
      , entityStacks_{EntityStack{gridViews.grid().maxLevel()}...}
    {
      contains_.reserve(sizeof...(GridViews));
      Hybrid::forEach(std::forward_as_tuple(gridViews...), [this](auto const& gv) {
        contains_.emplace_back([gv](const HostEntity& entity) { return gv.contains(entity); });
      });

      for (std::size_t i = 0; i < sizeof...(GridViews); ++i)
        initialIncrement(i);
    }

    template <class... GridViews>
    MultiMeshIterator (tag::with_gridview, bool endDummy, GridViews const&... gridViews)
      : macroIterators_{gridViews.grid().levelGridView(0).template end<0,pitype>()...}
    {}

Praetorius, Simon's avatar
Praetorius, Simon committed
164 165 166 167 168 169 170 171 172 173 174 175 176 177

    /// prefix increment
    void increment ()
    {
      for (std::size_t i = 0; i < entityStacks_.size(); ++i)
        incrementAllowed_[i] = incrementAllowed(i);

      for (std::size_t i = 0; i < entityStacks_.size(); ++i) {
        if (incrementAllowed_[i])
          increment(i);
      }
    }

    /// dereferencing
Praetorius, Simon's avatar
Praetorius, Simon committed
178
    MultiEntity<GridImp> const& dereference () const
Praetorius, Simon's avatar
Praetorius, Simon committed
179
    {
180 181 182 183 184 185
      // update entries in multiEntity that have changed
      for (std::size_t i = 0; i < entityStacks_.size(); ++i) {
        if (incrementAllowed_[i])
          multiEntity_[i] = dereference(i);
      }
      return multiEntity_;
Praetorius, Simon's avatar
Praetorius, Simon committed
186 187 188 189 190 191 192 193 194 195 196 197 198 199 200 201 202 203 204 205 206 207 208 209 210 211 212 213 214 215 216 217 218 219 220
    }

    /// equality
    bool equals (const MultiMeshIterator& that) const
    {
      return macroIterators_ == that.macroIterators_;
    }

  protected:
    // got to next entity in grid i
    void increment (std::size_t i)
    {
      auto& entityStack = entityStacks_[i];
      auto& macroIt = macroIterators_[i];
      auto const& macroEnd = macroEndIterators_[i];

      // 1. go up in tree or to next entity on current level until we can go down again
      while (!entityStack.empty()) {
        auto& top = entityStack.top();
        ++top.it;
        if (top.it == top.end) {
          entityStack.pop();
        } else {
          break;
        }
      }

      // 2. if entityStack is empty, go to next macroElement
      if (entityStack.empty()) {
        ++macroIt;
        if (macroIt == macroEnd)
          return;
      }

      // 3. go down in tree until leaf entity
221
      for (auto child = dereference(i); !entityTest(i,child);
222
                child = dereference(i)) {
223
        entityStack.emplace(child);
224
        assert( entityStack.size() <= maxLevel_[i] );
225
      }
Praetorius, Simon's avatar
Praetorius, Simon committed
226 227 228 229 230 231 232 233 234 235 236 237 238 239 240 241 242 243 244 245 246 247 248 249
    }

    /// Return true, if all stacks with size > stack[i].size are finished
    bool incrementAllowed (std::size_t i) const
    {
      std::size_t size = entityStacks_[i].size();
      return std::accumulate(entityStacks_.begin(), entityStacks_.end(), true,
        [i,size](bool allowed, auto const& entityStack) {
          return allowed && (entityStack.size() <= size || entityStack.finished(size));
        });
    }

    // got to first leaf entity on grid i
    void initialIncrement (std::size_t i)
    {
      auto& entityStack = entityStacks_[i];
      auto& macroIt = macroIterators_[i];
      auto const& macroEnd = macroEndIterators_[i];

      assert( entityStack.empty() );
      if (macroIt == macroEnd)
        return;

      // 2. go down in tree until leaf entity
250
      for (auto child = dereference(i); !entityTest(i,child);
251
                child = dereference(i)) {
252
        entityStack.emplace(child);
253
        assert( entityStack.size() <= maxLevel_[i] );
254
      }
Praetorius, Simon's avatar
Praetorius, Simon committed
255 256 257 258
    }

    HostEntity dereference (std::size_t i) const
    {
259 260
      if (entityStacks_[i].empty()) {
        assert(macroIterators_[i] != macroEndIterators_[i]);
261
        return *macroIterators_[i];
262 263
      } else {
        assert(entityStacks_[i].top().it != entityStacks_[i].top().end);
264
        return *entityStacks_[i].top().it;
265
      }
Praetorius, Simon's avatar
Praetorius, Simon committed
266 267
    }

268
    bool entityTest (std::size_t i, HostEntity const& entity) const
Praetorius, Simon's avatar
Praetorius, Simon committed
269
    {
270
      return contains_[i](entity) || entity.isLeaf() || !entity.isRegular();
Praetorius, Simon's avatar
Praetorius, Simon committed
271 272 273
    }

  private:
274 275 276 277 278
    std::vector<bool> incrementAllowed_;
    std::vector<EntityTest> contains_;
    std::vector<int> maxLevel_;

    mutable MultiEntity<GridImp> multiEntity_;
Praetorius, Simon's avatar
Praetorius, Simon committed
279 280 281 282 283 284

    std::vector<HostGridLevelIterator> macroIterators_;
    std::vector<HostGridLevelIterator> macroEndIterators_;
    std::vector<EntityStack> entityStacks_;
  };

285 286 287 288 289 290 291 292 293 294 295 296 297
  template <class> class MultiMesh;

  template <class... GridViews>
  inline auto multi_elements(GridViews const&... gridViews)
  {
    using GridView0 = std::tuple_element_t<0,std::tuple<GridViews...>>;
    using GridImp = MultiMesh<typename GridView0::Grid>;
    using Iterator = MultiMeshIterator<0,All_Partition,GridImp>;
    using Range = IteratorRange<Iterator>;
    return Range{ Iterator{tag::with_gridview{}, gridViews...},
                  Iterator{tag::with_gridview{}, true, gridViews...} };
  }

Praetorius, Simon's avatar
Praetorius, Simon committed
298 299 300
}  // end namespace Dune

#endif // DUNE_MULTIMESH_ITERATOR_HH