Merge branch 'feature/gmsh-boundary-flag-for-alugrid' into 'master'

Feature/gmsh boundary flag for alugrid

See merge request !1672

dumux/common/boundaryflag.hh

@@ -25,24 +25,24 @@
 #define DUMUX_BOUNDARY_FLAG_HH
 
 #include <cstddef>
+#include <limits>
 
 namespace Dumux {
 
 /*!
- * \file
- * \ingroup Common
- * \brief Boundary flag to store e.g. in sub control volume faces
- * \note Can be specialized for each grid manager (in the gridmanager headers)
- * \tparam Grid the type of the grid
+ * \ingroup InputOutput
+ * \brief Class for accessing boundary flags
+ * \note this works for all grid managers with gmsh meshes.
  */
-template<class Grid>
-class BoundaryFlag
+class BoundarySegmentIndexFlag
 {
 public:
-    BoundaryFlag() : flag_(-1) {}
+    BoundarySegmentIndexFlag()
+    : flag_(std::numeric_limits<std::size_t>::max()) {}
 
     template<class Intersection>
-    BoundaryFlag(const Intersection& i) : flag_(-1)
+    BoundarySegmentIndexFlag(const Intersection& i)
+    : flag_(std::numeric_limits<std::size_t>::max())
     {
         if (i.boundary())
             flag_ = i.boundarySegmentIndex();
@@ -56,6 +56,17 @@ private:
     value_type flag_;
 };
 
+/*!
+ * \file
+ * \ingroup Common
+ * \brief Boundary flag to store e.g. in sub control volume faces
+ * \note Can be specialized for each grid manager (in the gridmanager headers)
+ * \tparam Grid the type of the grid
+ */
+template<class Grid>
+class BoundaryFlag : public BoundarySegmentIndexFlag
+{ using BoundarySegmentIndexFlag::BoundarySegmentIndexFlag; };
+
 }  // end namespace Dumux
 
 #endif

test/io/gridmanager/CMakeLists.txt

@@ -234,3 +234,18 @@ dune_add_test(NAME test_gridmanager_vtk_ug
                        --command ${CMAKE_CURRENT_BINARY_DIR}/test_gridmanager_vtk_ug
                        --files ${CMAKE_CURRENT_BINARY_DIR}/grids/lens.vtu
                                ${CMAKE_CURRENT_BINARY_DIR}/test-gridmanager-vtk-uggrid-0.vtu)
+
+
+dune_add_test(NAME test_gmshboundaryflag
+              SOURCES test_gmshboundaryflag.cc
+              COMPILE_DEFINITIONS ENABLE_CACHING=false
+              LABELS unit io
+              CMAKE_GUARD "( dune-alugrid_FOUND )"
+              )
+
+dune_add_test(NAME test_gmshboundaryflag_caching
+              SOURCES test_gmshboundaryflag.cc
+              COMPILE_DEFINITIONS ENABLE_CACHING=true
+              LABELS unit io
+              CMAKE_GUARD "( dune-alugrid_FOUND )"
+              )

test/io/gridmanager/gmshboundaryflagtest.hh

+/*****************************************************************************
+ *   See the file COPYING for full copying permissions.                      *
+ *                                                                           *
+ *   This program is free software: you can redistribute it and/or modify    *
+ *   it under the terms of the GNU General Public License as published by    *
+ *   the Free Software Foundation, either version 3 of the License, or       *
+ *   (at your option) any later version.                                     *
+ *                                                                           *
+ *   This program is distributed in the hope that it will be useful,         *
+ *   but WITHOUT ANY WARRANTY; without even the implied warranty of          *
+ *   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the            *
+ *   GNU General Public License for more details.                            *
+ *                                                                           *
+ *   You should have received a copy of the GNU General Public License       *
+ *   along with this program.  If not, see <http://www.gnu.org/licenses/>.   *
+ *****************************************************************************/
+/*!
+ * \file
+ *
+ * \brief Test whether GmshBoundaryFlag works as expected
+ * \note GmshBoundaryFlag only exists for ALUGrid
+ *
+ * This tests whether the boundary marker IDs accessed using boundary flags are
+ * correct when using ALUGrid, a Gmsh mesh file, and the class GmshBoundaryFlag.
+ */
+#ifndef DUMUX_TEST_IO_GMSHBOUNDARYFLAG_TEST_HH
+#define DUMUX_TEST_IO_GMSHBOUNDARYFLAG_TEST_HH
+
+#include <dumux/io/grid/gridmanager.hh>
+#include <dumux/discretization/method.hh>
+
+#include <iostream>
+
+namespace Dumux {
+
+template<class Grid>
+class GmshBoundaryFlagTest
+{
+    using GridView = typename Grid::LeafGridView;
+    using Scalar = double;
+    static const int dim = Grid::dimension;
+    using GridManager = typename Dumux::GridManager<Grid>;
+    using ReferenceElements = typename Dune::ReferenceElements<Scalar, dim>;
+
+public:
+
+    template<class FVGridGeometry>
+    static void testGmshBoundaryFlag(const GridView& leafGridView,
+                                     std::shared_ptr<const FVGridGeometry> fvGridGeometry,
+                                     std::shared_ptr<const GridData<Grid>> gridData)
+    {
+
+        for(const auto& element : elements(leafGridView))
+        {
+            auto fvGeometry = localView(*fvGridGeometry);
+            fvGeometry.bind(element);
+
+            for (auto&& scvf : scvfs(fvGeometry))
+            {
+                if (scvf.boundary())
+                {
+                    const auto boundaryMarkerId = gridData->getBoundaryDomainMarker(scvf.boundaryFlag());
+                    const auto& pos = scvf.center();
+                    std::cout << "z-coordinate: " << pos[dim-1] << ", actual ID = " << boundaryMarkerId << ", ";
+
+                    /* According to unitcube.geo:
+                     * top = 1, bottom = 2, side = 3
+                     *
+                     * According to unitcube.msh:
+                     * top = 1, bottom = 0, side = 2
+                     */
+                    const int topId = 1;
+                    const int bottomId = 0;
+                    const int sideId = 2;
+
+                    const bool isTop = pos[dim-1] > 1.0 - eps_;
+                    const bool isBottom = pos[dim-1] < eps_;
+                    const bool isSide = !isTop && !isBottom;
+
+                    if (isTop)
+                    {
+                        std::cout << "correct ID = " << topId << " (is top surface)" << std::endl;
+                        if (boundaryMarkerId != topId)
+                            DUNE_THROW(Dune::Exception, "BoundaryMarkerId for top is wrong!");
+                    }
+                    else if (isBottom)
+                    {
+                        std::cout << "correct ID = " << bottomId << " (is bottom surface)" << std::endl;
+                        if (boundaryMarkerId != bottomId)
+                            DUNE_THROW(Dune::Exception, "BoundaryMarkerId for bottom is wrong!");
+                    }
+                    else if (isSide)
+                    {
+                        std::cout << "correct ID = " << sideId << " (is side surface)" << std::endl;
+                        if (boundaryMarkerId != sideId)
+                            DUNE_THROW(Dune::Exception, "BoundaryMarkerId for side is wrong!");
+                    }
+                } // end if boundary
+            } // end scvf loop
+        } // end element loop
+    } // end testGmshBoundaryFlag
+
+private:
+    static constexpr Scalar eps_ = 1e-4;
+
+};
+
+
+} // end namespace Dumux
+
+
+#endif /* DUMUX_TEST_IO_GMSHBOUNDARYFLAG_TEST_HH */

test/io/gridmanager/test_gmshboundaryflag.cc

+/*****************************************************************************
+ *   See the file COPYING for full copying permissions.                      *
+ *                                                                           *
+ *   This program is free software: you can redistribute it and/or modify    *
+ *   it under the terms of the GNU General Public License as published by    *
+ *   the Free Software Foundation, either version 3 of the License, or       *
+ *   (at your option) any later version.                                     *
+ *                                                                           *
+ *   This program is distributed in the hope that it will be useful,         *
+ *   but WITHOUT ANY WARRANTY; without even the implied warranty of          *
+ *   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the            *
+ *   GNU General Public License for more details.                            *
+ *                                                                           *
+ *   You should have received a copy of the GNU General Public License       *
+ *   along with this program.  If not, see <http://www.gnu.org/licenses/>.   *
+ *****************************************************************************/
+/*!
+ * \file
+ *
+ * \brief Test whether BoundarySegmentIndexFlag works as expected with Box and CCTpfa and ALUGrid
+ * \note Alu currently defaults to a boundary flag that works for DGF files only
+ */
+#include <config.h>
+#include <iostream>
+
+#if HAVE_DUNE_ALUGRID
+
+#include <dune/alugrid/grid.hh>
+#include <dune/common/parallel/mpihelper.hh>
+#include <dumux/common/parameters.hh>
+#include <dumux/common/boundaryflag.hh>
+#include <dumux/io/grid/gridmanager.hh>
+
+#include <dumux/discretization/box.hh>
+#include <dumux/discretization/cctpfa.hh>
+
+#include "gmshboundaryflagtest.hh"
+
+namespace Dumux {
+
+// In order to use an alternative BoundaryFlag class, we have to adapt the GridGeometryTraits
+template<class GridView>
+struct MyBoxGridGeometryTraits : public BoxDefaultGridGeometryTraits<GridView>
+{
+    struct MyScvfTraits : public BoxDefaultScvfGeometryTraits<GridView>
+    { using BoundaryFlag = BoundarySegmentIndexFlag; };
+
+    using SubControlVolumeFace = BoxSubControlVolumeFace<GridView, MyScvfTraits>;
+};
+
+template<class GridView>
+struct MyCCTpfaGridGeometryTraits : public CCTpfaDefaultGridGeometryTraits<GridView>
+{
+    struct MyScvfTraits : public CCTpfaDefaultScvfGeometryTraits<GridView>
+    { using BoundaryFlag = BoundarySegmentIndexFlag; };
+
+    using SubControlVolumeFace = CCTpfaSubControlVolumeFace<GridView, MyScvfTraits>;
+};
+
+} // end namespace Dumux
+
+
+int main(int argc, char** argv) try
+{
+    using namespace Dumux;
+
+    Dune::MPIHelper::instance(argc, argv);
+
+    Parameters::init(argc, argv, "test_gmshboundaryflag.input");
+
+    // generate the grid manager and initialize the grid
+    using Grid = Dune::ALUGrid<3, 3, Dune::simplex, Dune::conforming>;
+    GridManager<Grid> gridManager;
+    gridManager.init();
+    const auto& leafGridView = gridManager.grid().leafGridView();
+    auto gridData = gridManager.getGridData();
+
+    ////////////////////////////////////////////////////////////////////////////////////
+    // Test using Box discretization
+    ////////////////////////////////////////////////////////////////////////////////////
+
+    // generate FV grid geometry
+    using BoxFVGridGeometry = BoxFVGridGeometry<double, typename Grid::LeafGridView,
+                                             ENABLE_CACHING,
+                                             MyBoxGridGeometryTraits<
+                                                 typename Grid::LeafGridView
+                                                 >
+                                             >;
+    auto boxFvGridGeometry = std::make_shared<BoxFVGridGeometry>(leafGridView);
+    boxFvGridGeometry->update();
+
+    // run the test
+    GmshBoundaryFlagTest<Grid>::testGmshBoundaryFlag<BoxFVGridGeometry>(leafGridView, boxFvGridGeometry, gridData);
+
+    ////////////////////////////////////////////////////////////////////////////////////
+    // Test using CCTpfa discretization
+    ////////////////////////////////////////////////////////////////////////////////////
+
+    // generate FV grid geometry
+    using CCTpfaFVGridGeometry = CCTpfaFVGridGeometry<typename Grid::LeafGridView,
+                                             ENABLE_CACHING,
+                                             MyCCTpfaGridGeometryTraits<
+                                                 typename Grid::LeafGridView
+                                                 >
+                                             >;
+    auto ccTpfaFvGridGeometry = std::make_shared<CCTpfaFVGridGeometry>(leafGridView);
+    ccTpfaFvGridGeometry->update();
+
+    // run the test
+    GmshBoundaryFlagTest<Grid>::testGmshBoundaryFlag<CCTpfaFVGridGeometry>(leafGridView, ccTpfaFvGridGeometry, gridData);
+
+    return 0;
+}
+catch (Dumux::ParameterException &e) {
+    std::cerr << e << ". Abort!\n";
+    return 1;
+}
+catch (Dune::Exception &e) {
+    std::cerr << "Dune reported error: " << e << std::endl;
+    return 3;
+}
+
+#endif /* HAVE_DUNE_ALUGRID */

test/io/gridmanager/test_gmshboundaryflag.input

+[Grid]
+File = ./grids/unitcube.msh
+Verbosity = 1 # verbose grid file parsing
+DomainMarkers = true # read domain markers (gmsh physical entities)
+BoundarySegments = true # read boundary markers (gmsh physical entities)