mmiteratorbase.hh

// -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
// vi: set et ts=4 sw=2 sts=2:
#ifndef DUNE_MULTIMESH_ITERATOR_BASE_HH
#define DUNE_MULTIMESH_ITERATOR_BASE_HH

#include <numeric>
#include <stack>
#include <type_traits>
#if DUNE_HAVE_CXX_VARIANT
#include <variant>
#endif

#include <dune/common/std/type_traits.hh>
#include <dune/grid/common/gridenums.hh>

#include "mmentity.hh"
#include "mmhierarchiciterator.hh"
#include "mmiteratorinterface.hh"

namespace Dune
{
  template <PartitionIteratorType pitype, class HostGrid, class Derived>
  class MultiMeshIteratorBase
      : public MultiMeshIteratorInterface<HostGrid,Derived>
  {
    using Super = MultiMeshIteratorInterface<HostGrid,Derived>;

  protected:
    // LevelIterator to the equivalent entity in the host grid
    using HostGridLevelIterator =
      typename HostGrid::template Codim<0>::template Partition<pitype>::LevelIterator;

    template <class GV>
    using IsGridView = Std::disjunction<
      std::is_same<std::decay_t<GV>,typename HostGrid::LeafGridView>,
      std::is_same<std::decay_t<GV>,typename HostGrid::LevelGridView> >;

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

    using EntityTest = typename MultiMeshHierarchicIterator<HostGrid>::EntityTest;
    using EntityStack = typename MultiMeshHierarchicIterator<HostGrid>::EntityStack;

  public:

    /// Constructor which creates the begin iterator
    /**
     * \param  multiMesh  Pointer to grid instance
     */
    template <class GridImp>
    MultiMeshIteratorBase (tag::begin_iterator, const GridImp* multiMesh)
      : maxLevel_(multiMesh->size())
    {
      for (std::size_t i = 0; i < multiMesh->size(); ++i) {
        macroIterators_.push_back(multiMesh->levelGridView(i,0).template begin<0,pitype>());
        macroEndIterators_.push_back(multiMesh->levelGridView(i,0).template end<0,pitype>());
      }

      // go to first leaf entity on all grids
      entityStacks_.reserve(multiMesh->size());
      for (std::size_t i = 0; i < multiMesh->size(); ++i) {
        maxLevel_[i] = multiMesh->maxLevel(i);
        entityStacks_.emplace_back(multiMesh->maxLevel(i));
      }
    }

    /// Constructor which creates the end iterator
    /**
     * \param  multiMesh  Pointer to grid instance
     */
    template <class GridImp>
    MultiMeshIteratorBase (tag::end_iterator, const GridImp* multiMesh)
    {
      for (std::size_t i = 0; i < multiMesh->size(); ++i)
        macroIterators_.push_back(multiMesh->levelGridView(i,0).template end<0,pitype>());
    }

    /// Construct an iterator from n gridViews of possibly different type
    template <class... GridViews,
      std::enable_if_t<Std::conjunction<IsGridView<GridViews>...>::value, int> = 0>
    MultiMeshIteratorBase (tag::begin_iterator, GridViews const&... gridViews)
      : maxLevel_{gridViews.grid().maxLevel()...}
      , macroIterators_{gridViews.grid().levelGridView(0).template begin<0,pitype>()...}
      , macroEndIterators_{gridViews.grid().levelGridView(0).template end<0,pitype>()...}
      , entityStacks_{EntityStack{gridViews.grid().maxLevel()}...}
    {}

    template <class... GridViews,
      std::enable_if_t<Std::conjunction<IsGridView<GridViews>...>::value, int> = 0>
    MultiMeshIteratorBase (tag::end_iterator, GridViews const&... gridViews)
      : macroIterators_{gridViews.grid().levelGridView(0).template end<0,pitype>()...}
    {}

#if DUNE_HAVE_CXX_VARIANT
    /// Construct an iterator from a vector of GridViews, where each GridView can be
    /// either a leaf- or a level-gridView, thus we need to use a variant.
    template <class... GV>
    MultiMeshIteratorBase (tag::begin_iterator, const std::vector<std::variant<GV...>>& gridViews)
    {
      entityStacks_.reserve(gridViews.size());
      maxLevel_.reserve(gridViews.size());

      for (auto const& gvs : gridViews) {
        macroIterators_.push_back(std::visit([](auto const& gv) {
          return gv.grid().levelGridView(0).template begin<0,pitype>(); }, gvs));
        macroEndIterators_.push_back(std::visit([](auto const& gv) {
          return gv.grid().levelGridView(0).template end<0,pitype>(); }, gvs));
        maxLevel_.push_back(std::visit([](auto const& gv) {
          return gv.grid().maxLevel(); }, gvs));
        entityStacks_.emplace_back(std::visit([](auto const& gv) {
          return gv.grid().maxLevel(); }, gvs));
      }
    }

    template <class... GV>
    MultiMeshIteratorBase (tag::end_iterator, const std::vector<std::variant<GV...>>& gridViews)
    {
      for (auto const& gvs : gridViews) {
        macroIterators_.push_back(std::visit([](auto const& gv) {
          return gv.grid().levelGridView(0).template end<0,pitype>(); }, gvs));
      }
    }
#endif

    /// equality
    bool equals_impl (const MultiMeshIteratorBase& that) const
    {
      return macroIterators_ == that.macroIterators_;
    }

    // got to next entity in grid i
    int 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 0;
      }

      // 3. go down in tree until leaf entity
      auto child = dereference(i);
      for (; !this->levelReached(i, child); child = dereference(i)) {
        entityStack.emplace(child);
        assert(int(entityStack.size()) <= maxLevel_[i]);
      }

      // 4. go up in tree again to the first regular entity, since
      // irregular element can not be traversed in a multi-mesh sense
      while (!child.isRegular() && !entityStack.empty()) {
        entityStack.pop();
        child = dereference(i);
      }

      return entityStack.size();
    }
    using Super::increment;

    /// 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,
        [size](bool allowed, auto const& entityStack) {
          return allowed && (entityStack.size() <= size || entityStack.finished(size));
        });
    }

    // go to first entity on grid i on the desired level
    int 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 0;

      // 1. go down in tree until desired level is reached
      auto child = dereference(i);
      for (; !this->levelReached(i, child); child = dereference(i)) {
        entityStack.emplace(child);
        assert(int(entityStack.size()) <= maxLevel_[i]);
      }

      // 2. go up in tree again to the first regular entity, since
      // irregular element can not be traversed in a multi-mesh sense
      while (!child.isRegular() && !entityStack.empty()) {
        entityStack.pop();
        child = dereference(i);
      }

      return entityStack.size();
    }
    using Super::initialIncrement;

    // Return the current entity the hierarchic iterator or macro iterator pointing to
    HostEntity dereference (std::size_t i) const
    {
      if (entityStacks_[i].empty()) {
        assert(macroIterators_[i] != macroEndIterators_[i]);
        return *macroIterators_[i];
      } else {
        assert(entityStacks_[i].top().it != entityStacks_[i].top().end);
        return *entityStacks_[i].top().it;
      }
    }
    using Super::dereference;

    std::size_t size () const
    {
      return entityStacks_.size();
    }

  private:
    std::vector<int> maxLevel_;

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

}  // end namespace Dune

#endif // DUNE_MULTIMESH_ITERATOR_BASE_HH