Commit 3cdbc33e authored by Praetorius, Simon's avatar Praetorius, Simon

iterator for tuples of gridViews

parent f809a5b4
......@@ -7,6 +7,7 @@
#include <stack>
#include <dune/common/iteratorfacades.hh>
#include <dune/common/std/type_traits.hh>
#include <dune/grid/common/exceptions.hh>
#include <dune/grid/common/gridenums.hh>
......@@ -14,6 +15,11 @@
namespace Dune
{
namespace tag
{
struct with_gridview {};
}
/** \brief Iterator over all entities of a given codimension and level of a grid.
* \ingroup MultiMesh
*/
......@@ -48,6 +54,8 @@ namespace Dune
using HostEntity = typename GridImp::HostGridType::template Codim<0>::Entity;
using EntityTest = std::function<bool(HostEntity)>;
struct EntityStackEntry
{
template <class Entity>
......@@ -88,21 +96,22 @@ namespace Dune
public:
/// Constructor. Stores a pointer to the grid
explicit MultiMeshIterator (const GridImp* multiMesh, int level)
: multiMesh_(multiMesh)
, incrementAllowed_(multiMesh->size(), true)
, level_(level)
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)
, multiEntity_(multiMesh->size())
{
for (auto const& grid : multiMesh_->grids_) {
for (auto const& grid : multiMesh->grids_) {
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
entityStacks_.reserve(multiMesh_->size());
for (std::size_t i = 0; i < multiMesh_->size(); ++i) {
entityStacks_.emplace_back((*multiMesh_)[i].maxLevel());
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());
initialIncrement(i);
}
}
......@@ -113,13 +122,9 @@ namespace Dune
* \param endDummy Here only to distinguish it from the other constructor
*/
MultiMeshIterator (const GridImp* multiMesh, int level, bool endDummy)
: multiMesh_(multiMesh)
, level_(level)
{
for (auto const& grid : multiMesh_->grids_) {
for (auto const& grid : multiMesh->grids_)
macroIterators_.push_back(grid->levelGridView(0).template end<0,pitype>());
macroEndIterators_.push_back(grid->levelGridView(0).template end<0,pitype>());
}
}
/// Constructor which create the leaf-iterator
......@@ -132,6 +137,30 @@ namespace Dune
: MultiMeshIterator{multiMesh, -1, endDummy}
{}
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>()...}
{}
/// prefix increment
void increment ()
......@@ -189,10 +218,10 @@ namespace Dune
}
// 3. go down in tree until leaf entity
for (auto child = dereference(i); !entityTest(child);
for (auto child = dereference(i); !entityTest(i,child);
child = dereference(i)) {
entityStack.emplace(child);
assert( entityStack.size() <= (*multiMesh_)[i].maxLevel() );
assert( entityStack.size() <= maxLevel_[i] );
}
}
......@@ -218,10 +247,10 @@ namespace Dune
return;
// 2. go down in tree until leaf entity
for (auto child = dereference(i); !entityTest(child);
for (auto child = dereference(i); !entityTest(i,child);
child = dereference(i)) {
entityStack.emplace(child);
assert( entityStack.size() <= (*multiMesh_)[i].maxLevel() );
assert( entityStack.size() <= maxLevel_[i] );
}
}
......@@ -236,24 +265,36 @@ namespace Dune
}
}
bool entityTest (HostEntity const& entity)
bool entityTest (std::size_t i, HostEntity const& entity) const
{
return entity.level() == level_ || entity.isLeaf() || !entity.isRegular();
return contains_[i](entity) || entity.isLeaf() || !entity.isRegular();
}
private:
const GridImp* multiMesh_;
std::vector<bool> incrementAllowed_;
std::vector<EntityTest> contains_;
std::vector<int> maxLevel_;
mutable MultiEntity<GridImp> multiEntity_;
std::vector<HostGridLevelIterator> macroIterators_;
std::vector<HostGridLevelIterator> macroEndIterators_;
std::vector<EntityStack> entityStacks_;
std::vector<bool> incrementAllowed_;
int level_ = -1;
mutable MultiEntity<GridImp> multiEntity_;
};
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...} };
}
} // end namespace Dune
#endif // DUNE_MULTIMESH_ITERATOR_HH
dune_add_test(SOURCES testvolume.cc)
dune_add_test(SOURCES testnumentities.cc)
dune_add_test(SOURCES testtransform.cc)
\ No newline at end of file
dune_add_test(SOURCES testtransform.cc)
dune_add_test(SOURCES testgridviews.cc)
\ No newline at end of file
// -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
// vi: set et ts=4 sw=2 sts=2:
#ifdef HAVE_CONFIG_H
# include "config.h"
#endif
#include <functional>
#include <iostream>
#include <numeric>
#include <dune/common/filledarray.hh>
#include <dune/common/test/testsuite.hh>
#include <dune/common/parallel/mpihelper.hh>
#include <dune/grid/yaspgrid.hh>
#include <dune/multimesh/multimesh.hh>
using namespace Dune;
template <std::size_t dim, class Test>
void test_dim(Test& test)
{
FieldVector<double,dim> lower; lower = -1.5;
FieldVector<double,dim> upper; upper = 1.5;
auto num_elements = filledArray<dim>(2);
using HostGrid = YaspGrid<dim, EquidistantOffsetCoordinates<double,dim>>;
MultiMesh<HostGrid> grid(3, lower, upper, num_elements);
grid[0].globalRefine(1);
grid[1].globalRefine(2);
grid[2].globalRefine(3);
for (auto const& entities : multi_elements(grid[0].leafGridView(), grid[1].levelGridView(1))) {
test.check(entities.size() == 2);
test.check(entities[0].isLeaf());
test.check(entities[1].level() == 1);
}
}
int main(int argc, char** argv)
{
MPIHelper::instance(argc, argv);
Dune::TestSuite test;
test_dim<1>(test);
test_dim<2>(test);
test_dim<3>(test);
return test.exit();
}
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment