Commit a1edbfb1 authored by Praetorius, Simon's avatar Praetorius, Simon

MultiGridView and MultiIndexSet added

parent 3cdbc33e
set(modules "DuneMultimeshMacros.cmake")
install(FILES ${modules} DESTINATION ${DUNE_INSTALL_MODULEDIR})
install(FILES
DuneMultimeshMacros.cmake
CXXFeatures.cmake
DESTINATION ${DUNE_INSTALL_MODULEDIR})
include(CheckCXXSourceCompiles)
# support for C++17's variant implementation
check_cxx_source_compiles("
#include <variant>
int main()
{
std::variant<int,double> a;
a = 1;
a = 2.0;
bool positive = std::visit([](auto&& x) { return x > 0; }, a);
}
" DUNE_HAVE_CXX_VARIANT )
# File for module specific CMake tests.
include(CXXFeatures)
......@@ -40,6 +40,9 @@
/* Define to the revision of dune-multimesh */
#define DUNE_MULTIMESH_VERSION_REVISION @DUNE_MULTIMESH_VERSION_REVISION@
/* some detected compiler features may be used in this module */
#cmakedefine DUNE_HAVE_CXX_VARIANT 1
/* end dune-multimesh
Everything below here will be overwritten
*/
......@@ -12,23 +12,23 @@
namespace Dune
{
template <class GridImp>
template <class HostGrid>
class MultiEntity
: public std::vector<typename GridImp::HostGridType::Traits::template Codim<0>::Entity>
: public std::vector<typename HostGrid::Traits::template Codim<0>::Entity>
{
private:
using ctype = typename GridImp::ctype;
using ctype = typename HostGrid::ctype;
public:
enum { dimension = GridImp::dimension };
enum { dimension = HostGrid::dimension };
/// The type of a local geometry
using LocalGeometry = Dune::Geometry<dimension, dimension, const GridImp, MultiMeshLocalGeometry>;
using LocalGeometry = MultiMeshLocalGeometry<dimension, dimension, HostGrid>;
using EntitySeed = typename GridImp::Traits::template Codim<0>::EntitySeed;
using EntitySeed = MultiMeshEntitySeed<0, HostGrid>;
/// Entity in the host grid
using HostGridEntity = typename GridImp::HostGridType::Traits::template Codim<0>::Entity;
using HostGridEntity = typename HostGrid::Traits::template Codim<0>::Entity;
/// Containertype
using Super = std::vector<HostGridEntity>;
......@@ -39,15 +39,13 @@ namespace Dune
/// Return a local geometry of source in target
static LocalGeometry localGeometry (HostGridEntity const& source, HostGridEntity const& target)
{
using LocalGeometryImp = MultiMeshLocalGeometry<dimension, dimension, const GridImp>;
return LocalGeometry{LocalGeometryImp{source, target}};
return {source, target};
}
/// Return a local geometry of source_i'th entity in target_t'th entity
LocalGeometry localGeometry (std::size_t source_i, std::size_t target_i) const
{
using LocalGeometryImp = MultiMeshLocalGeometry<dimension, dimension, const GridImp>;
return LocalGeometry{LocalGeometryImp{(*this)[source_i], (*this)[target_i]}};
return {(*this)[source_i], (*this)[target_i]};
}
/// \brief Return the entity seed which contains sufficient information
......@@ -58,8 +56,7 @@ namespace Dune
**/
EntitySeed seed (std::size_t entity_i) const
{
using EntitySeedImp = MultiMeshEntitySeed<0, const GridImp>;
return EntitySeedImp{(*this)[entity_i], entity_i};
return {(*this)[entity_i], entity_i};
}
/// Return the entity with maximal level
......
......@@ -9,23 +9,23 @@
namespace Dune
{
template <int mydim, int coorddim, class GridImp>
template <int mydim, int coorddim, class HostGrid>
class MultiMeshLocalGeometry
: public GeometryDefaultImplementation<mydim, coorddim, GridImp, MultiMeshLocalGeometry>
: public GeometryDefaultImplementation<mydim, coorddim, HostGrid, MultiMeshLocalGeometry>
{
private:
using ctype = typename GridImp::ctype;
using ctype = typename HostGrid::ctype;
public:
// The codimension of this entitypointer wrt the host grid
enum { codimension = GridImp::HostGridType::dimension - mydim };
enum { dimensionworld = GridImp::HostGridType::dimensionworld };
enum { codimension = HostGrid::dimension - mydim };
enum { dimensionworld = HostGrid::dimensionworld };
// select appropriate hostgrid geometry via typeswitch
using HostLocalGeometry = typename GridImp::HostGridType::Traits::template Codim<codimension>::LocalGeometry;
using HostLocalGeometry = typename HostGrid::Traits::template Codim<codimension>::LocalGeometry;
/// entity in the host grid
using HostGridEntity = typename GridImp::HostGridType::Traits::template Codim<codimension>::Entity;
using HostGridEntity = typename HostGrid::Traits::template Codim<codimension>::Entity;
//! type of jacobian transposed
using JacobianInverseTransposed = typename HostLocalGeometry::JacobianInverseTransposed;
......
......@@ -5,6 +5,7 @@
#include <numeric>
#include <stack>
#include <variant>
#include <dune/common/iteratorfacades.hh>
#include <dune/common/std/type_traits.hh>
......@@ -23,7 +24,7 @@ namespace Dune
/** \brief Iterator over all entities of a given codimension and level of a grid.
* \ingroup MultiMesh
*/
template <int codim, PartitionIteratorType pitype, class GridImp>
template <int codim, PartitionIteratorType pitype, class HostGrid>
class MultiMeshIterator
{
public:
......@@ -41,18 +42,18 @@ namespace Dune
// Implemented for codim 0 entities only
template <PartitionIteratorType pitype, class GridImp>
class MultiMeshIterator<0, pitype, GridImp>
: public ForwardIteratorFacade<MultiMeshIterator<0,pitype,GridImp>,
MultiEntity<GridImp>,
MultiEntity<GridImp> const&>
template <PartitionIteratorType pitype, class HostGrid>
class MultiMeshIterator<0, pitype, HostGrid>
: public ForwardIteratorFacade<MultiMeshIterator<0,pitype,HostGrid>,
MultiEntity<HostGrid>,
MultiEntity<HostGrid> const&>
{
private:
// LevelIterator to the equivalent entity in the host grid
using HostGridLevelIterator =
typename GridImp::HostGridType::template Codim<0>::template Partition<pitype>::LevelIterator;
typename HostGrid::template Codim<0>::template Partition<pitype>::LevelIterator;
using HostEntity = typename GridImp::HostGridType::template Codim<0>::Entity;
using HostEntity = typename HostGrid::template Codim<0>::Entity;
using EntityTest = std::function<bool(HostEntity)>;
......@@ -96,6 +97,7 @@ namespace Dune
public:
/// Constructor. Stores a pointer to the grid
template <class GridImp>
MultiMeshIterator (const GridImp* multiMesh, int level)
: incrementAllowed_(multiMesh->size(), true)
, contains_(multiMesh->size(), EntityTest{[level](const HostEntity& entity) { return entity.level() == level; }})
......@@ -121,6 +123,7 @@ namespace Dune
* \param multiMesh Pointer to grid instance
* \param endDummy Here only to distinguish it from the other constructor
*/
template <class GridImp>
MultiMeshIterator (const GridImp* multiMesh, int level, bool endDummy)
{
for (auto const& grid : multiMesh->grids_)
......@@ -128,11 +131,13 @@ namespace Dune
}
/// Constructor which create the leaf-iterator
template <class GridImp>
MultiMeshIterator (const GridImp* multiMesh)
: MultiMeshIterator{multiMesh, -1}
{}
/// Constructor which create the end leaf-iterator
template <class GridImp>
MultiMeshIterator (const GridImp* multiMesh, bool endDummy)
: MultiMeshIterator{multiMesh, -1, endDummy}
{}
......@@ -161,6 +166,41 @@ namespace Dune
: macroIterators_{gridViews.grid().levelGridView(0).template end<0,pitype>()...}
{}
#if DUNE_HAVE_CXX_VARIANT
template <class... GV>
MultiMeshIterator (const std::vector<std::variant<GV...>>& gridViews)
: incrementAllowed_(gridViews.size(), true)
, multiEntity_(gridViews.size())
{
contains_.reserve(gridViews.size());
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));
contains_.emplace_back(std::visit([](auto const& gv) {
return EntityTest{[gv](const HostEntity& entity) { return gv.contains(entity); }}; }, 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));
}
// go to first leaf entity on all grids
for (std::size_t i = 0; i < entityStacks_.size(); ++i)
initialIncrement(i);
}
template <class... GV>
MultiMeshIterator (const std::vector<std::variant<GV...>>& gridViews, bool endDummy)
{
for (auto const& gvs : gridViews) {
macroIterators_.push_back(std::visit([](auto const& gv) {
return gv.grid().levelGridView(0).template end<0,pitype>(); }, gvs));
}
}
#endif
/// prefix increment
void increment ()
......@@ -175,7 +215,7 @@ namespace Dune
}
/// dereferencing
MultiEntity<GridImp> const& dereference () const
MultiEntity<HostGrid> const& dereference () const
{
// update entries in multiEntity that have changed
for (std::size_t i = 0; i < entityStacks_.size(); ++i) {
......@@ -275,7 +315,7 @@ namespace Dune
std::vector<EntityTest> contains_;
std::vector<int> maxLevel_;
mutable MultiEntity<GridImp> multiEntity_;
mutable MultiEntity<HostGrid> multiEntity_;
std::vector<HostGridLevelIterator> macroIterators_;
std::vector<HostGridLevelIterator> macroEndIterators_;
......
// -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
// vi: set et ts=4 sw=2 sts=2:
#ifndef DUNE_MULTI_GRIDVIEW_HH
#define DUNE_MULTI_GRIDVIEW_HH
#if ! DUNE_HAVE_CXX_VARIANT
#error "Require C++17 variant!"
#endif
#include <variant>
#include <dune/common/typetraits.hh>
#include <dune/common/exceptions.hh>
#include <dune/common/std/type_traits.hh>
#include <dune/grid/common/capabilities.hh>
#include <dune/grid/common/gridview.hh>
#include "mmiterator.hh"
namespace Dune
{
template <class HostGrid>
class MultiIndexSet
{
public:
using LeafIndexSet = typename HostGrid::Traits::LeafIndexSet;
using LevelIndexSet = typename HostGrid::Traits::LevelIndexSet;
using IndexType = std::common_type_t<typename LeafIndexSet::IndexType, typename LevelIndexSet::IndexType>;
using IndexSetTypes = std::variant<LeafIndexSet const*, LevelIndexSet const*>;
template <class IndexSet>
MultiIndexSet(IndexSet const* indexSet)
: indexSets_(indexSet)
{}
MultiIndexSet(IndexSetTypes const& indexSet)
: indexSets_(indexSet)
{}
template <class Entity>
IndexType index (const Entity& e) const
{
return std::visit([&e](auto const* is) -> IndexType { return is->index(e); }, indexSets_);
}
template <class Entity>
IndexType subIndex (const Entity& e, int i, unsigned int codim) const
{
return std::visit([&e,i,codim](auto const* is) -> IndexType { return is->subIndex(e,i,codim); }, indexSets_);
}
private:
IndexSetTypes indexSets_;
};
// forward declaration
template <class HostGrid>
class MultiGridView;
template <class HostGrid>
struct MultiGridViewTraits
{
using Grid = typename std::remove_const<HostGrid>::type;
using IndexSet = MultiIndexSet<HostGrid>;
using Intersection = std::variant<
typename HostGrid::Traits::LeafIntersection, typename HostGrid::Traits::LevelIntersection>;
using IntersectionIterator = std::variant<
typename HostGrid::Traits::LeafIntersectionIterator, typename HostGrid::Traits::LevelIntersectionIterator>;
using CollectiveCommunication = typename HostGrid::Traits::CollectiveCommunication;
template <int cd>
struct Codim
{
using Entity = typename Grid::Traits::template Codim<cd>::Entity;
using Geometry = typename Grid::template Codim<cd>::Geometry;
using LocalGeometry = typename Grid::template Codim<cd>::LocalGeometry;
template <PartitionIteratorType pit>
struct Partition
{
using Iterator = MultiMeshIterator<cd,pit,HostGrid>;
};
};
enum { conforming = Capabilities::isLevelwiseConforming<HostGrid>::v };
};
template <class HostGrid>
class MultiGridView
{
public:
using Traits = MultiGridViewTraits<HostGrid>;
/// The MultiMesh GridType
using Grid = typename Traits::Grid;
using GridViewTypes = std::variant<typename HostGrid::LeafGridView, typename HostGrid::LevelGridView>;
using IndexSet = typename Traits::IndexSet;
using IndexSetTypes = std::variant<const typename HostGrid::LeafIndexSet*, const typename HostGrid::LevelIndexSet*>;
using IntersectionIterator = typename Traits::IntersectionIterator;
using CollectiveCommunication = typename Traits::CollectiveCommunication;
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> >;
template <int cd>
struct Codim : public Traits::template Codim<cd> {};
enum { dimension = HostGrid::dimension };
enum { dimensionworld = HostGrid::dimensionworld };
enum { conforming = Traits::conforming };
public:
/// Constructor.
template <class... GridViews,
std::enable_if_t<Std::conjunction<IsGridView<GridViews>...>::value, int> = 0>
MultiGridView (GridViews&&... gridViews)
: gridViews_{std::forward<GridViews>(gridViews)...}
{}
template <class GridView,
std::enable_if_t<IsGridView<GridView>::value, int> = 0>
MultiGridView (const std::vector<GridView>& gridViews)
: gridViews_(gridViews.begin(), gridViews.end())
{}
template <class Iter,
std::enable_if_t<IsGridView<typename std::iterator_traits<Iter>::value_type>::value, int> = 0>
MultiGridView (Iter first, Iter last)
: gridViews_(first, last)
{}
template <class Iter,
std::enable_if_t<IsGridView<typename std::iterator_traits<Iter>::value_type>::value, int> = 0>
MultiGridView (const IteratorRange<Iter>& gridViews)
: gridViews_(gridViews.begin(), gridViews.end())
{}
/// Obtain a const reference to the underlying hierarchic grid
const HostGrid& grid (std::size_t idx) const
{
return std::visit([](auto const& gv) { return gv.grid(); }, gridViews_[idx]);
}
/// Obtain the level-indexSet
// NOTE: IndexSet is one of {LeafIndexSet, LevelIndexSet}
// NOTE: do not return a reference, but a reference wrapper
IndexSet indexSet (std::size_t idx) const
{
return std::visit([](auto const& gv) -> IndexSetTypes { return &gv.indexSet(); }, gridViews_[idx]);
}
/// Obtain number of entities in a given codimension
int size (std::size_t idx, int codim) const
{
return std::visit([codim](auto const& gv) { return gv.size(codim); }, gridViews_[idx]);
}
/// Obtain number of entities with a given geometry type
int size (std::size_t idx, const GeometryType& type) const
{
return std::visit([&type](auto const& gv) { return gv.size(type); }, gridViews_[idx]);
}
/// Obtain begin iterator for this view
template <int cd, PartitionIteratorType pit = All_Partition>
typename Codim<cd>::template Partition<pit>::Iterator begin () const
{
static_assert(cd == 0, "Implemented for codim == 0 only");
return MultiMeshIterator<cd,pit,HostGrid>(gridViews_);
}
/// Obtain end iterator for this view
template <int cd, PartitionIteratorType pit = All_Partition>
typename Codim<cd>::template Partition<pit>::Iterator end () const
{
static_assert(cd == 0, "Implemented for codim == 0 only");
return MultiMeshIterator<cd,pit,HostGrid>(gridViews_, true);
}
/// Obtain begin intersection iterator with respect to this view
// NOTE: IntersectionIterator is one of {LeafIntersectionIterator, LevelIntersectionIterator}
IntersectionIterator ibegin (std::size_t idx, const typename Codim<0>::Entity& entity) const
{
return std::visit([&entity](auto const& gv) { return gv.ibegin(entity); }, gridViews_[idx]);
}
/// Obtain end intersection iterator with respect to this view
// NOTE: IntersectionIterator is one of {LeafIntersectionIterator, LevelIntersectionIterator}
IntersectionIterator iend (std::size_t idx, const typename Codim<0>::Entity& entity) const
{
return std::visit([&entity](auto const& gv) { return gv.iend(entity); }, gridViews_[idx]);
}
/// Obtain collective communication object
const CollectiveCommunication& comm (std::size_t idx) const
{
return std::visit([](auto const& gv) -> const auto& { return gv.comm(); }, gridViews_[idx]);
}
/// Return size of the overlap region for a given codim on the grid view.
int overlapSize (std::size_t idx, int codim) const
{
return std::visit([codim](auto const& gv) { return gv.overlapSize(codim); }, gridViews_[idx]);
}
/// Return size of the ghost region for a given codim on the grid view.
int ghostSize (std::size_t idx, int codim) const
{
return std::visit([codim](auto const& gv) { return gv.ghostSize(codim); }, gridViews_[idx]);
}
/// Communicate data on this view
template <class DataHandleImp, class DataType>
void communicate (std::size_t idx,
CommDataHandleIF<DataHandleImp, DataType>& data,
InterfaceType iftype,
CommunicationDirection dir) const
{
std::visit([&data,iftype,dir](auto const& gv) { gv.communicate(data,iftype,dir); }, gridViews_[idx]);
}
private:
std::vector<GridViewTypes> gridViews_;
};
} // end namespace Dune
#endif // DUNE_MULTI_GRIDVIEW_HH
......@@ -46,11 +46,11 @@ namespace Dune
{
/// The type of the iterator over the level entities of this codim on this partition.
// using LevelIterator = Dune::EntityIterator<cd, const Grid, MultiMeshIterator<cd, pitype, const Grid> >;
using LevelIterator = MultiMeshIterator<cd, pitype, const Grid>;
using LevelIterator = MultiMeshIterator<cd, pitype, HostGrid>;
/// The type of the iterator over the leaf entities of this codim on this partition.
// using LeafIterator = Dune::EntityIterator<cd, const Grid, MultiMeshIterator<cd, pitype, const Grid> >;
using LeafIterator = MultiMeshIterator<cd, pitype, const Grid>;
using LeafIterator = MultiMeshIterator<cd, pitype, HostGrid>;
};
/// The type of the iterator over all leaf entities of this codim.
......@@ -164,28 +164,28 @@ namespace Dune
template <int codim, PartitionIteratorType PiType = All_Partition>
typename Traits::template Codim<codim>::template Partition<PiType>::LevelIterator lbegin (int level) const
{
return MultiMeshIterator<codim,PiType, const MultiMesh<HostGrid> >(this, level);
return MultiMeshIterator<codim,PiType, HostGrid>(this, level);
}
/// one past the end on this level
template <int codim, PartitionIteratorType PiType = All_Partition>
typename Traits::template Codim<codim>::template Partition<PiType>::LevelIterator lend (int level) const
{
return MultiMeshIterator<codim,PiType, const MultiMesh<HostGrid> >(this, level, true);
return MultiMeshIterator<codim,PiType, HostGrid>(this, level, true);
}
/// Iterator to first leaf entity of given codim
template <int codim, PartitionIteratorType PiType = All_Partition>
typename Traits::template Codim<codim>::template Partition<PiType>::LeafIterator leafbegin () const
{
return MultiMeshIterator<codim,PiType, const MultiMesh<HostGrid> >(this);
return MultiMeshIterator<codim,PiType, HostGrid>(this);
}
/// one past the end of the sequence of leaf entities
template <int codim, PartitionIteratorType PiType = All_Partition>
typename Traits::template Codim<codim>::template Partition<PiType>::LeafIterator leafend () const
{
return MultiMeshIterator<codim,PiType, const MultiMesh<HostGrid> >(this, true);
return MultiMeshIterator<codim,PiType, HostGrid>(this, true);
}
......
......@@ -116,7 +116,7 @@ namespace Dune
* \param upperRight Upper right corner of the grid
* \param elements Number of elements in each coordinate direction
**/
static std::unique_ptr<GridType> createSimplexGrid (
static void createSimplexGrid (
GridFactory<GridType>& factory,
const FieldVector<ctype,dimworld>& lowerLeft,
const FieldVector<ctype,dimworld>& upperRight,
......@@ -166,9 +166,6 @@ namespace Dune
permutation.end()));
}
}
// Create the grid and hand it to the calling method
return std::unique_ptr<GridType>(factory.createGrid());
}
static std::unique_ptr<GridType> createSimplexGrid (
......@@ -177,7 +174,8 @@ namespace Dune
const std::array<unsigned int,dim>& elements)
{
GridFactory<GridType> factory;
return createCubeGrid(factory, lowerLeft, upperRight, elements);
createSimplexGrid(factory, lowerLeft, upperRight, elements);
return std::unique_ptr<GridType>(factory.createGrid());
}
private:
......
......@@ -9,3 +9,6 @@ target_link_dune_default_libraries("localrefinement")
add_executable("uggrid" uggrid.cc)
target_link_dune_default_libraries("uggrid")
add_executable("multigridview" multigridview.cc)
target_link_dune_default_libraries("multigridview")
......@@ -69,8 +69,8 @@ int main(int argc, char** argv)
using Factory = StructuredGridBuilder<MultiMesh<HostGrid> >;
GridFactory<MultiMesh<HostGrid> > gridFactory(3);
auto gridPtr = Factory::createSimplexGrid(gridFactory, lower_left, bbox, num_elements);
Factory::createSimplexGrid(gridFactory, lower_left, bbox, num_elements);
auto gridPtr = gridFactory.createGrid();
auto& grid = *gridPtr;
std::cout << "number of grids = " << grid.size() << "\n";
......
// -*- 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 <iostream>
#include <dune/common/parallel/mpihelper.hh> // An initializer of MPI
#include <dune/common/exceptions.hh> // We use exceptions
#include <dune/common/timer.hh>
#if HAVE_DUNE_ALUGRID
#include <dune/alugrid/grid.hh>
#endif
#include <dune/multimesh/multigridview.hh>
#include <dune/multimesh/utility/structuredgridbuilder.hh>
using namespace Dune;
template <class GridView>
void printGrid(GridView const& gv)
{
volatile std::size_t n = 0;
Dune::Timer t;
for (auto const& entity : elements(gv))
n += gv.indexSet().index(entity);
std::cout << n << "\n";
std::cout << "time: " << t.elapsed() << "\n";
}
template <class GridView>
void printGrid2(GridView const& gv)
{
volatile std::size_t n = 0;
Dune::Timer t;
for (auto const& entities : elements(gv)) {
n += gv.indexSet(0).index(entities[0]);
std::cout << "indices = [";
for (std::size_t i = 0; i < entities.size(); ++i) {
std::cout << (i > 0 ? ", " : "") << gv.indexSet(i).index(entities[i]);
}
std::cout << "]\n";
}
std::cout << "n = " << n << "\n";
std::cout << "time: " << t.elapsed() << "\n";
}
int main(int argc, char** argv)
{
MPIHelper::instance(argc, argv);
#if HAVE_DUNE_ALUGRID
FieldVector<double,2> lower_left = {-1.5, -1.5};
FieldVector<double,2> bbox = {1.5, 1.5};
std::array<unsigned int,2> num_elements = {2, 2};
using Grid = Dune::ALUGrid<2, 2, Dune::simplex, Dune::conforming>;
using Factory = StructuredGridBuilder<Grid>;
auto gridPtr = Factory::createSimplexGrid(lower_left, bbox, num_elements);
auto& grid = *gridPtr;
grid.globalRefine(4);
printGrid(grid.leafGridView());
printGrid(grid.levelGridView(2));
using GridView = MultiGridView<Grid>;
GridView gridView(grid.leafGridView(), grid.levelGridView(2));
printGrid2(gridView);
#endif
}
......@@ -70,8 +70,8 @@ int main(int argc, char** argv)
using Factory = StructuredGridBuilder<MultiMesh<HostGrid> >;
GridFactory<MultiMesh<HostGrid> > gridFactory(3);
auto gridPtr = Factory::createSimplexGrid(gridFactory, lower_left, bbox, num_elements);
Factory::createSimplexGrid(gridFactory, lower_left, bbox, num_elements);
auto gridPtr = gridFactory.createGrid();
auto& grid = *gridPtr;
#elif GRID == 3 && HAVE_ALBERTA
......
......@@ -30,6 +30,7 @@ int main(int argc, char** argv)
using Factory = StructuredGridBuilder<MultiMesh<HostGrid> >;
GridFactory<MultiMesh<HostGrid> > factory(1);
auto gridPtr = Factory::createSimplexGrid(factory, lower_left, bbox, num_elements);
Factory::createSimplexGrid(factory, lower_left, bbox, num_elements);
auto gridPtr = factory.createGrid();
#endif
}
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