Add static create method to MeshCreator and implement boundaryIds extraction from gmsh grids

amdis/MeshCreator.hpp

@@ -66,6 +66,12 @@ namespace AMDiS
       : name_(name)
     {}
 
+    /// Static create mthod. See \ref create()
+    static std::shared_ptr<Grid> create(std::string name)
+    {
+      return MeshCreator{name}.create();
+    }
+
     /// Create a new grid by inspecting the intifile parameter group `[meshName]`
     /**
      * Reads first the parameter `[meshName]->macro file name` and if set
@@ -78,7 +84,7 @@ namespace AMDiS
      *
      * Otherwise tries to create a grid depending on the grid type.
      **/
-    std::shared_ptr<Grid> create() const
+    std::unique_ptr<Grid> create() const
     {
       auto filename = Parameters::get<std::string>(name_ + "->macro file name");
       auto structured = Parameters::get<std::string>(name_ + "->structured");
@@ -142,15 +148,15 @@ namespace AMDiS
     }
 
   private:
-    static std::shared_ptr<Grid> construct(std::unique_ptr<Grid> hostGrid)
+    static std::unique_ptr<Grid> construct(std::unique_ptr<Grid> hostGrid)
     {
       return std::move(hostGrid);
     }
 
     template <class HG>
-    static std::shared_ptr<Grid> construct(std::unique_ptr<HG> hostGrid)
+    static std::unique_ptr<Grid> construct(std::unique_ptr<HG> hostGrid)
     {
-      return std::make_shared<Grid>(std::move(hostGrid));
+      return std::make_unique<Grid>(std::move(hostGrid));
     }
 
     // use the structured grid factory to create the grid
@@ -204,13 +210,55 @@ namespace AMDiS
       std::declval<std::vector<unsigned int>>(),
       std::declval<std::shared_ptr<Dune::BoundarySegment<GridType::dimension, GridType::dimensionworld> >>()) );
 
+    template <class GridType, class LC>
+    using SupportsInsertionIndex = decltype(std::declval<Dune::GridFactory<GridType>>().insertionIndex(std::declval<LC>()));
+
     // use GmshReader if GridFactory supports insertBoundarySegments
     template <class GridType = HostGrid,
       REQUIRES(Dune::Std::is_detected<SupportsGmshReader, GridType>::value)>
     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_)};
+      Dune::GridFactory<GridType> factory;
+      std::vector<int> boundaryIds, elementIds;
+      reader.read(factory, filename, boundaryIds, elementIds);
+
+      using HasInsertionIndexEntity
+        = Dune::Std::is_detected<SupportsInsertionIndex, GridType, typename GridType::template Codim<0>::Entity>;
+      using HasInsertionIndexIntersection
+        = Dune::Std::is_detected<SupportsInsertionIndex, GridType, typename GridType::LeafIntersection>;
+
+      auto gridPtr = factory.createGrid();
+      if ((boundaryIds.empty() && elementIds.empty()) ||
+          (!HasInsertionIndexEntity::value && !HasInsertionIndexIntersection::value))
+        return std::unique_ptr<GridType>(std::move(gridPtr));
+
+      // map boundaryIds and elementIds read from file to grid indexing.
+      if (!boundaryIds.empty() && HasInsertionIndexIntersection::value)
+        boundaryIds_.resize(gridPtr->numBoundarySegments());
+      if (!elementIds.empty() && HasInsertionIndexEntity::value)
+        elementIds_.resize(gridPtr->size(0));
+
+      auto const& indexSet = gridPtr->leafIndexSet();
+      for (auto const& e : elements(gridPtr->leafGridView())) {
+        if constexpr(HasInsertionIndexEntity::value) {
+          if (!elementIds.empty())
+            elementIds_[indexSet.index(e)] = elementIds[factory.insertionIndex(e)];
+        }
+
+        if (boundaryIds.empty() || !e.hasBoundaryIntersections())
+          continue;
+
+        for (auto const& it : intersections(gridPtr->leafGridView(), e)) {
+          if constexpr(HasInsertionIndexIntersection::value) {
+            if (it.boundary())
+              boundaryIds_[it.boundarySegmentIndex()]
+                = boundaryIds[factory.insertionIndex(it)];
+          }
+        }
+      }
+
+      return std::unique_ptr<GridType>(std::move(gridPtr));
     }
 
     // fallback if GmshReader cannot be used