From c740863ce7109ee0289836723fdfc18944166aa7 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?St=C3=A9phane=20Del=20Pino?= <stephane.delpino44@gmail.com>
Date: Tue, 5 Apr 2022 16:25:33 +0200
Subject: [PATCH] Change referenced item list policy

Now in 1D, referenced face and edge lists are defined according to
referenced node lists. And in 2D, referenced edge lists are defined
according to referenced face lists.

Note that lists are duplicated. This could be improved (sharing lists)
but it would require (significant?) changes. Since the extra memory
cost only occurs in 1D and 2D, it does not seem necessary.
---
 src/mesh/Connectivity.cpp       |  31 ++++
 src/mesh/MeshFaceBoundary.cpp   |  23 ++-
 tests/test_MeshFaceBoundary.cpp | 247 +++++++++++++++++++++++---------
 3 files changed, 234 insertions(+), 67 deletions(-)

diff --git a/src/mesh/Connectivity.cpp b/src/mesh/Connectivity.cpp
index 7b81eeb03..005dc1f65 100644
--- a/src/mesh/Connectivity.cpp
+++ b/src/mesh/Connectivity.cpp
@@ -95,6 +95,24 @@ Connectivity<Dimension>::_buildFrom(const ConnectivityDescriptor& descriptor)
     m_edge_number   = WeakEdgeValue<int>(*this, node_number_array);
     m_edge_owner    = WeakEdgeValue<int>(*this, node_owner_array);
     m_edge_is_owned = WeakEdgeValue<bool>(*this, node_is_owned_array);
+
+    // edge and face references are set equal to node references
+    m_ref_edge_list_vector.reserve(descriptor.template refItemListVector<ItemType::node>().size());
+    m_ref_face_list_vector.reserve(descriptor.template refItemListVector<ItemType::node>().size());
+    for (auto ref_node_list : descriptor.template refItemListVector<ItemType::node>()) {
+      const RefId ref_id            = ref_node_list.refId();
+      Array<const NodeId> node_list = ref_node_list.list();
+      Array<EdgeId> edge_list(node_list.size());
+      Array<FaceId> face_list(node_list.size());
+      for (size_t i = 0; i < node_list.size(); ++i) {
+        edge_list[i] = EdgeId::base_type{node_list[i]};
+        face_list[i] = FaceId::base_type{node_list[i]};
+      }
+
+      m_ref_edge_list_vector.emplace_back(RefItemList<ItemType::edge>(ref_id, edge_list));
+      m_ref_face_list_vector.emplace_back(RefItemList<ItemType::face>(ref_id, face_list));
+    }
+
   } else {
     m_item_to_item_matrix[itemTId(ItemType::face)][itemTId(ItemType::node)] = descriptor.face_to_node_vector;
 
@@ -134,6 +152,19 @@ Connectivity<Dimension>::_buildFrom(const ConnectivityDescriptor& descriptor)
       m_edge_owner    = WeakEdgeValue<int>(*this, face_owner_array);
       m_edge_is_owned = WeakEdgeValue<bool>(*this, face_is_owned_array);
 
+      // edge references are set equal to face references
+      m_ref_edge_list_vector.reserve(descriptor.template refItemListVector<ItemType::face>().size());
+      for (auto ref_face_list : descriptor.template refItemListVector<ItemType::face>()) {
+        const RefId ref_id            = ref_face_list.refId();
+        Array<const FaceId> face_list = ref_face_list.list();
+        Array<EdgeId> edge_list(face_list.size());
+        for (size_t i = 0; i < face_list.size(); ++i) {
+          edge_list[i] = EdgeId::base_type{face_list[i]};
+        }
+
+        m_ref_edge_list_vector.emplace_back(RefItemList<ItemType::edge>(ref_id, edge_list));
+      }
+
     } else {
       m_item_to_item_matrix[itemTId(ItemType::edge)][itemTId(ItemType::node)] = descriptor.edge_to_node_vector;
 
diff --git a/src/mesh/MeshFaceBoundary.cpp b/src/mesh/MeshFaceBoundary.cpp
index 751bcf962..4e718f251 100644
--- a/src/mesh/MeshFaceBoundary.cpp
+++ b/src/mesh/MeshFaceBoundary.cpp
@@ -9,6 +9,7 @@ MeshFaceBoundary<Dimension>::MeshFaceBoundary(const Mesh<Connectivity<Dimension>
   : m_face_list(ref_face_list.list()), m_boundary_name(ref_face_list.refId().tagName())
 {}
 
+template MeshFaceBoundary<1>::MeshFaceBoundary(const Mesh<Connectivity<1>>&, const RefFaceList&);
 template MeshFaceBoundary<2>::MeshFaceBoundary(const Mesh<Connectivity<2>>&, const RefFaceList&);
 template MeshFaceBoundary<3>::MeshFaceBoundary(const Mesh<Connectivity<3>>&, const RefFaceList&);
 
@@ -21,14 +22,32 @@ getMeshFaceBoundary(const Mesh<Connectivity<Dimension>>& mesh, const IBoundaryDe
     const auto& ref_face_list = mesh.connectivity().template refItemList<ItemType::face>(i_ref_face_list);
     const RefId& ref          = ref_face_list.refId();
     if (ref == boundary_descriptor) {
+      auto face_list           = ref_face_list.list();
+      auto face_to_cell_matrix = mesh.connectivity().faceToCellMatrix();
+
+      bool is_bad = false;
+      parallel_for(face_list.size(), [=, &is_bad](int l) {
+        const auto& face_cells = face_to_cell_matrix[face_list[l]];
+        if (face_cells.size() > 1) {
+          is_bad = true;
+        }
+      });
+
+      if (parallel::allReduceOr(is_bad)) {
+        std::ostringstream ost;
+        ost << "invalid boundary " << rang::fgB::yellow << boundary_descriptor << rang::style::reset
+            << ": inner faces cannot be used to define mesh boundaries";
+        throw NormalError(ost.str());
+      }
+
       return MeshFaceBoundary<Dimension>{mesh, ref_face_list};
     }
   }
 
   std::ostringstream ost;
-  ost << "cannot find surface with name " << rang::fgB::red << boundary_descriptor << rang::style::reset << '\n';
+  ost << "cannot find face list with name " << rang::fgB::red << boundary_descriptor << rang::style::reset << '\n';
   ost << "The mesh contains " << mesh.connectivity().template numberOfRefItemList<ItemType::face>()
-      << " face boundaries: ";
+      << " referenced face lists: ";
   for (size_t i_ref_face_list = 0; i_ref_face_list < mesh.connectivity().template numberOfRefItemList<ItemType::face>();
        ++i_ref_face_list) {
     const auto& ref_face_list = mesh.connectivity().template refItemList<ItemType::face>(i_ref_face_list);
diff --git a/tests/test_MeshFaceBoundary.cpp b/tests/test_MeshFaceBoundary.cpp
index 5b390e17f..3ab99dcd3 100644
--- a/tests/test_MeshFaceBoundary.cpp
+++ b/tests/test_MeshFaceBoundary.cpp
@@ -15,12 +15,34 @@ TEST_CASE("MeshFaceBoundary", "[mesh]")
 {
   auto is_same = [](const auto& a, const auto& b) -> bool {
     if (a.size() > 0 and b.size() > 0) {
-      return (a[0] == b[0]);
+      return (a.size() == b.size()) and (&(a[0]) == &(b[0]));
     } else {
       return (a.size() == b.size());
     }
   };
 
+  auto get_face_list_from_tag = [](const size_t tag, const auto& connectivity) -> Array<const FaceId> {
+    for (size_t i = 0; i < connectivity.template numberOfRefItemList<ItemType::face>(); ++i) {
+      const auto& ref_face_list = connectivity.template refItemList<ItemType::face>(i);
+      const RefId ref_id        = ref_face_list.refId();
+      if (ref_id.tagNumber() == tag) {
+        return ref_face_list.list();
+      }
+    }
+    return {};
+  };
+
+  auto get_face_list_from_name = [](const std::string& name, const auto& connectivity) -> Array<const FaceId> {
+    for (size_t i = 0; i < connectivity.template numberOfRefItemList<ItemType::face>(); ++i) {
+      const auto& ref_face_list = connectivity.template refItemList<ItemType::face>(i);
+      const RefId ref_id        = ref_face_list.refId();
+      if (ref_id.tagName() == name) {
+        return ref_face_list.list();
+      }
+    }
+    return {};
+  };
+
   SECTION("1D")
   {
     static constexpr size_t Dimension = 1;
@@ -35,20 +57,28 @@ TEST_CASE("MeshFaceBoundary", "[mesh]")
 
       const ConnectivityType& connectivity = mesh.connectivity();
 
-      for (size_t i = 0; i < connectivity.numberOfRefItemList<ItemType::face>(); ++i) {
-        const auto& ref_face_list = connectivity.refItemList<ItemType::face>(i);
-        const RefId ref_id        = ref_face_list.refId();
+      {
+        const std::set<size_t> tag_set = {0, 1};
 
-        {
-          NumberedBoundaryDescriptor numbered_boundary_descriptor(ref_id.tagNumber());
+        for (auto tag : tag_set) {
+          NumberedBoundaryDescriptor numbered_boundary_descriptor(tag);
           const auto& face_boundary = getMeshFaceBoundary(mesh, numbered_boundary_descriptor);
-          REQUIRE(is_same(face_boundary.faceList(), ref_face_list.list()));
+
+          auto face_list = get_face_list_from_tag(tag, connectivity);
+          REQUIRE(is_same(face_boundary.faceList(), face_list));
         }
+      }
 
-        {
-          NamedBoundaryDescriptor named_boundary_descriptor(ref_id.tagName());
+      {
+        const std::set<std::string> name_set = {"XMIN", "XMAX"};
+
+        for (auto name : name_set) {
+          NamedBoundaryDescriptor named_boundary_descriptor(name);
           const auto& face_boundary = getMeshFaceBoundary(mesh, named_boundary_descriptor);
-          REQUIRE(is_same(face_boundary.faceList(), ref_face_list.list()));
+
+          auto face_list = get_face_list_from_name(name, connectivity);
+
+          REQUIRE(is_same(face_boundary.faceList(), face_list));
         }
       }
     }
@@ -60,20 +90,27 @@ TEST_CASE("MeshFaceBoundary", "[mesh]")
 
       const ConnectivityType& connectivity = mesh.connectivity();
 
-      for (size_t i = 0; i < connectivity.numberOfRefItemList<ItemType::face>(); ++i) {
-        const auto& ref_face_list = connectivity.refItemList<ItemType::face>(i);
-        const RefId ref_id        = ref_face_list.refId();
+      {
+        const std::set<size_t> tag_set = {1, 2};
 
-        {
-          NumberedBoundaryDescriptor numbered_boundary_descriptor(ref_id.tagNumber());
+        for (auto tag : tag_set) {
+          NumberedBoundaryDescriptor numbered_boundary_descriptor(tag);
           const auto& face_boundary = getMeshFaceBoundary(mesh, numbered_boundary_descriptor);
-          REQUIRE(is_same(face_boundary.faceList(), ref_face_list.list()));
+
+          auto face_list = get_face_list_from_tag(tag, connectivity);
+          REQUIRE(is_same(face_boundary.faceList(), face_list));
         }
+      }
+
+      {
+        const std::set<std::string> name_set = {"XMIN", "XMAX"};
 
-        {
-          NamedBoundaryDescriptor named_boundary_descriptor(ref_id.tagName());
+        for (auto name : name_set) {
+          NamedBoundaryDescriptor named_boundary_descriptor(name);
           const auto& face_boundary = getMeshFaceBoundary(mesh, named_boundary_descriptor);
-          REQUIRE(is_same(face_boundary.faceList(), ref_face_list.list()));
+
+          auto face_list = get_face_list_from_name(name, connectivity);
+          REQUIRE(is_same(face_boundary.faceList(), face_list));
         }
       }
     }
@@ -93,20 +130,27 @@ TEST_CASE("MeshFaceBoundary", "[mesh]")
 
       const ConnectivityType& connectivity = mesh.connectivity();
 
-      for (size_t i = 0; i < connectivity.numberOfRefItemList<ItemType::face>(); ++i) {
-        const auto& ref_face_list = connectivity.refItemList<ItemType::face>(i);
-        const RefId ref_id        = ref_face_list.refId();
+      {
+        const std::set<size_t> tag_set = {0, 1, 2, 3};
 
-        {
-          NumberedBoundaryDescriptor numbered_boundary_descriptor(ref_id.tagNumber());
+        for (auto tag : tag_set) {
+          NumberedBoundaryDescriptor numbered_boundary_descriptor(tag);
           const auto& face_boundary = getMeshFaceBoundary(mesh, numbered_boundary_descriptor);
-          REQUIRE(is_same(face_boundary.faceList(), ref_face_list.list()));
+
+          auto face_list = get_face_list_from_tag(tag, connectivity);
+          REQUIRE(is_same(face_boundary.faceList(), face_list));
         }
+      }
 
-        {
-          NamedBoundaryDescriptor named_boundary_descriptor(ref_id.tagName());
-          const auto& face_boundary = getMeshFaceBoundary(mesh, named_boundary_descriptor);
-          REQUIRE(is_same(face_boundary.faceList(), ref_face_list.list()));
+      {
+        const std::set<std::string> name_set = {"XMIN", "XMAX", "YMIN", "YMAX"};
+
+        for (auto name : name_set) {
+          NamedBoundaryDescriptor numbered_boundary_descriptor(name);
+          const auto& face_boundary = getMeshFaceBoundary(mesh, numbered_boundary_descriptor);
+
+          auto face_list = get_face_list_from_name(name, connectivity);
+          REQUIRE(is_same(face_boundary.faceList(), face_list));
         }
       }
     }
@@ -118,20 +162,27 @@ TEST_CASE("MeshFaceBoundary", "[mesh]")
 
       const ConnectivityType& connectivity = mesh.connectivity();
 
-      for (size_t i = 0; i < connectivity.numberOfRefItemList<ItemType::face>(); ++i) {
-        const auto& ref_face_list = connectivity.refItemList<ItemType::face>(i);
-        const RefId ref_id        = ref_face_list.refId();
+      {
+        const std::set<size_t> tag_set = {1, 2, 3, 4};
 
-        {
-          NumberedBoundaryDescriptor numbered_boundary_descriptor(ref_id.tagNumber());
+        for (auto tag : tag_set) {
+          NumberedBoundaryDescriptor numbered_boundary_descriptor(tag);
           const auto& face_boundary = getMeshFaceBoundary(mesh, numbered_boundary_descriptor);
-          REQUIRE(is_same(face_boundary.faceList(), ref_face_list.list()));
+
+          auto face_list = get_face_list_from_tag(tag, connectivity);
+          REQUIRE(is_same(face_boundary.faceList(), face_list));
         }
+      }
 
-        {
-          NamedBoundaryDescriptor named_boundary_descriptor(ref_id.tagName());
-          const auto& face_boundary = getMeshFaceBoundary(mesh, named_boundary_descriptor);
-          REQUIRE(is_same(face_boundary.faceList(), ref_face_list.list()));
+      {
+        const std::set<std::string> name_set = {"XMIN", "YMIN", "XMAX", "YMIN"};
+
+        for (auto name : name_set) {
+          NamedBoundaryDescriptor numbered_boundary_descriptor(name);
+          const auto& face_boundary = getMeshFaceBoundary(mesh, numbered_boundary_descriptor);
+
+          auto face_list = get_face_list_from_name(name, connectivity);
+          REQUIRE(is_same(face_boundary.faceList(), face_list));
         }
       }
     }
@@ -151,21 +202,27 @@ TEST_CASE("MeshFaceBoundary", "[mesh]")
 
       const ConnectivityType& connectivity = mesh.connectivity();
 
-      for (size_t i = 0; i < connectivity.numberOfRefItemList<ItemType::face>(); ++i) {
-        const auto& ref_face_list = connectivity.refItemList<ItemType::face>(i);
-        const RefId ref_id        = ref_face_list.refId();
+      {
+        const std::set<size_t> tag_set = {0, 1, 2, 3, 4, 5};
+
+        for (auto tag : tag_set) {
+          NumberedBoundaryDescriptor numbered_boundary_descriptor(tag);
+          const auto& face_boundary = getMeshFaceBoundary(mesh, numbered_boundary_descriptor);
 
-        {
-          NumberedBoundaryDescriptor numbered_boundary_descriptor(ref_id.tagNumber());
-          // force copy for tests
-          auto face_boundary = getMeshFaceBoundary(mesh, numbered_boundary_descriptor);
-          REQUIRE(is_same(face_boundary.faceList(), ref_face_list.list()));
+          auto face_list = get_face_list_from_tag(tag, connectivity);
+          REQUIRE(is_same(face_boundary.faceList(), face_list));
         }
+      }
 
-        {
-          NamedBoundaryDescriptor named_boundary_descriptor(ref_id.tagName());
-          const auto& face_boundary = getMeshFaceBoundary(mesh, named_boundary_descriptor);
-          REQUIRE(is_same(face_boundary.faceList(), ref_face_list.list()));
+      {
+        const std::set<std::string> name_set = {"XMIN", "XMAX", "YMIN", "YMAX", "ZMIN", "ZMAX"};
+
+        for (auto name : name_set) {
+          NamedBoundaryDescriptor numbered_boundary_descriptor(name);
+          const auto& face_boundary = getMeshFaceBoundary(mesh, numbered_boundary_descriptor);
+
+          auto face_list = get_face_list_from_name(name, connectivity);
+          REQUIRE(is_same(face_boundary.faceList(), face_list));
         }
       }
     }
@@ -177,21 +234,27 @@ TEST_CASE("MeshFaceBoundary", "[mesh]")
 
       const ConnectivityType& connectivity = mesh.connectivity();
 
-      for (size_t i = 0; i < connectivity.numberOfRefItemList<ItemType::face>(); ++i) {
-        const auto& ref_face_list = connectivity.refItemList<ItemType::face>(i);
-        const RefId ref_id        = ref_face_list.refId();
+      {
+        const std::set<size_t> tag_set = {22, 23, 24, 25, 26, 27};
 
-        {
-          NumberedBoundaryDescriptor numbered_boundary_descriptor(ref_id.tagNumber());
-          // force copy for tests
-          auto face_boundary = getMeshFaceBoundary(mesh, numbered_boundary_descriptor);
-          REQUIRE(is_same(face_boundary.faceList(), ref_face_list.list()));
+        for (auto tag : tag_set) {
+          NumberedBoundaryDescriptor numbered_boundary_descriptor(tag);
+          const auto& face_boundary = getMeshFaceBoundary(mesh, numbered_boundary_descriptor);
+
+          auto face_list = get_face_list_from_tag(tag, connectivity);
+          REQUIRE(is_same(face_boundary.faceList(), face_list));
         }
+      }
 
-        {
-          NamedBoundaryDescriptor named_boundary_descriptor(ref_id.tagName());
-          const auto& face_boundary = getMeshFaceBoundary(mesh, named_boundary_descriptor);
-          REQUIRE(is_same(face_boundary.faceList(), ref_face_list.list()));
+      {
+        const std::set<std::string> name_set = {"XMIN", "XMAX", "YMIN", "YMAX", "ZMIN", "ZMAX"};
+
+        for (auto name : name_set) {
+          NamedBoundaryDescriptor numbered_boundary_descriptor(name);
+          const auto& face_boundary = getMeshFaceBoundary(mesh, numbered_boundary_descriptor);
+
+          auto face_list = get_face_list_from_name(name, connectivity);
+          REQUIRE(is_same(face_boundary.faceList(), face_list));
         }
       }
     }
@@ -212,9 +275,63 @@ TEST_CASE("MeshFaceBoundary", "[mesh]")
       NamedBoundaryDescriptor named_boundary_descriptor("invalid_boundary");
 
       REQUIRE_THROWS_WITH(getMeshFaceBoundary(mesh, named_boundary_descriptor),
-                          "error: cannot find surface with name \"invalid_boundary\"\nThe mesh contains 8 face "
-                          "boundaries: XMIN(22), XMAX(23), ZMAX(24), ZMIN(25), YMAX(26), YMIN(27), INTERFACE1(55), "
+                          "error: cannot find face list with name \"invalid_boundary\"\nThe mesh contains 8 referenced "
+                          "face lists: XMIN(22), XMAX(23), ZMAX(24), ZMIN(25), YMAX(26), YMIN(27), INTERFACE1(55), "
                           "INTERFACE2(56)");
     }
+
+    SECTION("surface is inside")
+    {
+      SECTION("1D")
+      {
+        static constexpr size_t Dimension = 1;
+
+        using ConnectivityType = Connectivity<Dimension>;
+        using MeshType         = Mesh<ConnectivityType>;
+
+        std::shared_ptr p_mesh = MeshDataBaseForTests::get().unordered1DMesh();
+        const MeshType& mesh   = *p_mesh;
+
+        NamedBoundaryDescriptor named_boundary_descriptor("INTERFACE");
+
+        REQUIRE_THROWS_WITH(getMeshFaceBoundary(mesh, named_boundary_descriptor),
+                            "error: invalid boundary \"INTERFACE\": inner faces cannot be used to define mesh "
+                            "boundaries");
+      }
+
+      SECTION("2D")
+      {
+        static constexpr size_t Dimension = 2;
+
+        using ConnectivityType = Connectivity<Dimension>;
+        using MeshType         = Mesh<ConnectivityType>;
+
+        std::shared_ptr p_mesh = MeshDataBaseForTests::get().hybrid2DMesh();
+        const MeshType& mesh   = *p_mesh;
+
+        NamedBoundaryDescriptor named_boundary_descriptor("INTERFACE");
+
+        REQUIRE_THROWS_WITH(getMeshFaceBoundary(mesh, named_boundary_descriptor),
+                            "error: invalid boundary \"INTERFACE\": inner faces cannot be used to define mesh "
+                            "boundaries");
+      }
+
+      SECTION("3D")
+      {
+        static constexpr size_t Dimension = 3;
+
+        using ConnectivityType = Connectivity<Dimension>;
+        using MeshType         = Mesh<ConnectivityType>;
+
+        std::shared_ptr p_mesh = MeshDataBaseForTests::get().hybrid3DMesh();
+        const MeshType& mesh   = *p_mesh;
+
+        NamedBoundaryDescriptor named_boundary_descriptor("INTERFACE1");
+
+        REQUIRE_THROWS_WITH(getMeshFaceBoundary(mesh, named_boundary_descriptor),
+                            "error: invalid boundary \"INTERFACE1\": inner faces cannot be used to define mesh "
+                            "boundaries");
+      }
+    }
   }
 }
-- 
GitLab