From fcd473a583f12b3be42780195c288264bc2d81d6 Mon Sep 17 00:00:00 2001
From: "Praetorius, Simon" <simon.praetorius@tu-dresden.de>
Date: Wed, 28 Aug 2019 11:11:12 +0200
Subject: [PATCH] changed grid construction to allow parallel grids to be
 created directly

---
 src/amdis/MeshCreator.hpp | 78 +++++++++++++++++++++++++++++++--------
 1 file changed, 63 insertions(+), 15 deletions(-)

diff --git a/src/amdis/MeshCreator.hpp b/src/amdis/MeshCreator.hpp
index 2dc00df3..41d1d0b9 100644
--- a/src/amdis/MeshCreator.hpp
+++ b/src/amdis/MeshCreator.hpp
@@ -22,6 +22,16 @@
 #include <amdis/common/TypeTraits.hpp>
 #include <amdis/utility/MacroGridFactory.hpp>
 
+namespace Dune
+{
+  // forward declarations
+  template <int dim, class Coordinates> class YaspGrid;
+  template <class ct, int dim> class EquidistantCoordinates;
+  template <class ct, int dim> class EquidistantOffsetCoordinates;
+  template <class ct, int dim> class TensorProductCoordinates;
+}
+
+
 namespace AMDiS
 {
   /// A creator class for dune grids.
@@ -49,6 +59,7 @@ namespace AMDiS
      * Otherwise reads the parameter `[meshName]->structured` and if set, creates:
      * - cube: a structured cube grid
      * - simplex: a structured simplex grid
+     * using a StructuredGridFactory
      *
      * Otherwise tries to create a grid depending on the grid type.
      **/
@@ -89,7 +100,7 @@ namespace AMDiS
     {
       return create_structured_grid([](auto&& lower, auto&& upper, auto&& numCells)
       {
-        return MacroGridFactory<Grid>::createSimplexGrid(lower, upper, numCells);
+        return Dune::StructuredGridFactory<Grid>::createSimplexGrid(lower, upper, numCells);
       });
     }
 
@@ -107,7 +118,7 @@ namespace AMDiS
 
   private:
     // use the structured grid factory to create the grid
-    template <class Factory>
+    template <class Size = unsigned int, class Factory>
     std::unique_ptr<Grid> create_structured_grid(Factory factory) const
     {
       // Lower left corner of the domain
@@ -115,7 +126,7 @@ namespace AMDiS
       // Upper right corner of the domain
       Dune::FieldVector<ctype,int(dimworld)> upper(1);
       // number of blocks in each coordinate direction
-      auto numCells = Dune::filledArray<std::size_t(dimension),unsigned int>(1);
+      auto numCells = Dune::filledArray<std::size_t(dimension),Size>(1);
 
       Parameters::get(name_ + "->min corner", lower);
       Parameters::get(name_ + "->max corner", upper);
@@ -136,7 +147,7 @@ namespace AMDiS
         return read_alberta_file<Grid>(filename, Dune::PriorityTag<42>{});
       }
       else {
-        error_exit("Can not read grid file. Unknown file extension.");
+        error_exit("Cannot read grid file. Unknown file extension.");
         return {};
       }
     }
@@ -149,7 +160,7 @@ namespace AMDiS
     // use GmshReader if GridFactory supports insertBoundarySegments
     template <class GridType = Grid,
       REQUIRES(Dune::Std::is_detected<SupportsGmshReader, GridType>::value)>
-    std::unique_ptr<GridType> read_gmsh_file(std::string const& filename, Dune::PriorityTag<2>) const
+    std::unique_ptr<GridType> read_gmsh_file(std::string const& filename, Dune::PriorityTag<1>) const
     {
       Dune::GmshReader<GridType> reader;
       return std::unique_ptr<GridType>{reader.read(filename, boundaryIds_, elementIds_)};
@@ -184,8 +195,10 @@ namespace AMDiS
     std::unique_ptr<GridType> read_alberta_file(std::string const& filename, Dune::PriorityTag<2>) const
     {
       Dune::GridFactory<GridType> factory;
-      Dune::AlbertaReader<GridType> reader;
-      reader.readGrid(filename, factory);
+      if (Environment::mpiRank() == 0) {
+        Dune::AlbertaReader<GridType> reader;
+        reader.readGrid(filename, factory);
+      }
       return std::unique_ptr<GridType>{factory.createGrid()};
     }
 
@@ -194,7 +207,7 @@ namespace AMDiS
       REQUIRES(GridType::dimensionworld != DIM_OF_WORLD)>
     std::unique_ptr<GridType> read_alberta_file(std::string const& filename, Dune::PriorityTag<1>) const
     {
-      error_exit("Alberta (and AlbertaReader) require WORLDDIM == Grid::dimensionworld. Change the cmake parameters!");
+      error_exit("AlbertaGrid (and AlbertaReader) require WORLDDIM == Grid::dimensionworld. Change the cmake parameters!");
       return {};
     }
 #else
@@ -202,7 +215,7 @@ namespace AMDiS
     template <class GridType>
     std::unique_ptr<GridType> read_alberta_file(std::string const& filename, Dune::PriorityTag<0>) const
     {
-      error_exit("Alberta (and AlbertaReader) not available. Set AlbertaFlags to your executable in cmake!");
+      error_exit("AlbertaGrid (and AlbertaReader) not available. Set AlbertaFlags to your executable in cmake!");
       return {};
     }
 #endif
@@ -214,7 +227,10 @@ namespace AMDiS
       REQUIRES(Dune::Std::is_detected<IsAlbertaGrid, GridType>::value)>
     std::unique_ptr<GridType> create_by_gridtype(Dune::PriorityTag<3>) const
     {
-      return create_simplex_grid();
+      return create_structured_grid([](auto&& lower, auto&& upper, auto&& numCells)
+      {
+        return MacroGridFactory<Grid>::createSimplexGrid(lower, upper, numCells);
+      });
     }
 #endif
 
@@ -223,9 +239,43 @@ namespace AMDiS
       class = typename GridType::YGrid>
     std::unique_ptr<GridType> create_by_gridtype(Dune::PriorityTag<2>) const
     {
-      return create_cube_grid();
+      int overlap = 0;
+      Parameters::get(name_ + "->overlap", overlap);
+      std::string periodic(dimension, '0');
+      Parameters::get(name_ + "->periodic", periodic); // e.g. 000 01 111
+
+      return create_yaspgrid(Types<GridType>{}, overlap, std::bitset<dimension>(periodic));
+    }
+
+    template <int dim, class ct>
+    std::unique_ptr<Grid> create_yaspgrid(Types<Dune::YaspGrid<dim,Dune::EquidistantCoordinates<ct,dim>>>,
+                                          int overlap, std::bitset<dimension> const& per) const
+    {
+      return create_structured_grid<int>([&](auto&& /*lower*/, auto&& upper, std::array<int,dimension> const& numCells)
+      {
+        return std::make_unique<Grid>(upper, numCells, per, overlap);
+      });
+    }
+
+    template <int dim, class ct>
+    std::unique_ptr<Grid> create_yaspgrid(Types<Dune::YaspGrid<dim,Dune::EquidistantOffsetCoordinates<ct,dim>>>,
+                                          int overlap, std::bitset<dimension> const& per) const
+    {
+      return create_structured_grid<int>([&](auto&& lower, auto&& upper, std::array<int,dimension> const& numCells)
+      {
+        return std::make_unique<Grid>(lower, upper, numCells, per, overlap);
+      });
     }
 
+    template <int dim, class ct>
+    std::unique_ptr<Grid> create_yaspgrid(Types<Dune::YaspGrid<dim,Dune::TensorProductCoordinates<ct,dim>>>,
+                                          int overlap, std::bitset<dimension> const& per) const
+    {
+      error_exit("MeshCreator cannot create YaspGrid with TensorProductCoordinates.");
+      return {};
+    }
+
+
 #if HAVE_DUNE_SPGRID
     // spgrid -> cube
     template <class GridType = Grid,
@@ -233,11 +283,9 @@ namespace AMDiS
       class = typename GridType::MultiIndex>
     std::unique_ptr<GridType> create_by_gridtype(Dune::PriorityTag<1>) const
     {
-      return create_structured_grid([](auto&& lower, auto&& upper, auto&& numCells)
+      return create_structured_grid<int>([](auto&& lower, auto&& upper, std::array<int,dimension> const& numCells)
       {
-        std::array<int,dimworld> cells;
-        std::copy(std::begin(numCells), std::end(numCells), std::begin(cells));
-        typename GridType::MultiIndex multiIndex(cells);
+        typename GridType::MultiIndex multiIndex(numCells);
         return std::make_unique<GridType>(lower, upper, multiIndex);
       });
     }
-- 
GitLab