diff --git a/src/mesh/CMakeLists.txt b/src/mesh/CMakeLists.txt
index afd3cb28a354332c35b0ca2757499c1e8da47d09..113e2527dbe4c3896979c3bd7f68a00369e62029 100644
--- a/src/mesh/CMakeLists.txt
+++ b/src/mesh/CMakeLists.txt
@@ -24,6 +24,7 @@ add_library(
MeshDataManager.cpp
MeshEdgeBoundary.cpp
MeshFaceBoundary.cpp
+ MeshFlatEdgeBoundary.cpp
MeshFlatFaceBoundary.cpp
MeshFlatNodeBoundary.cpp
MeshRelaxer.cpp
diff --git a/src/mesh/MeshFlatEdgeBoundary.cpp b/src/mesh/MeshFlatEdgeBoundary.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..d4c1ec155a0e8c4d863facbf6280d71265944933
--- /dev/null
+++ b/src/mesh/MeshFlatEdgeBoundary.cpp
@@ -0,0 +1,20 @@
+#include <mesh/MeshFlatEdgeBoundary.hpp>
+
+#include <mesh/Connectivity.hpp>
+#include <mesh/Mesh.hpp>
+#include <mesh/MeshFlatNodeBoundary.hpp>
+
+template <size_t Dimension>
+MeshFlatEdgeBoundary<Dimension>
+getMeshFlatEdgeBoundary(const Mesh<Connectivity<Dimension>>& mesh, const IBoundaryDescriptor& boundary_descriptor)
+{
+ MeshEdgeBoundary<Dimension> mesh_edge_boundary = getMeshEdgeBoundary(mesh, boundary_descriptor);
+ MeshFlatNodeBoundary<Dimension> mesh_flat_node_boundary = getMeshFlatNodeBoundary(mesh, boundary_descriptor);
+
+ return MeshFlatEdgeBoundary<Dimension>{mesh, mesh_edge_boundary.refEdgeList(),
+ mesh_flat_node_boundary.outgoingNormal()};
+}
+
+template MeshFlatEdgeBoundary<1> getMeshFlatEdgeBoundary(const Mesh<Connectivity<1>>&, const IBoundaryDescriptor&);
+template MeshFlatEdgeBoundary<2> getMeshFlatEdgeBoundary(const Mesh<Connectivity<2>>&, const IBoundaryDescriptor&);
+template MeshFlatEdgeBoundary<3> getMeshFlatEdgeBoundary(const Mesh<Connectivity<3>>&, const IBoundaryDescriptor&);
diff --git a/src/mesh/MeshFlatEdgeBoundary.hpp b/src/mesh/MeshFlatEdgeBoundary.hpp
new file mode 100644
index 0000000000000000000000000000000000000000..a104d1413a7319f08b74bf1ef103a5157b803f18
--- /dev/null
+++ b/src/mesh/MeshFlatEdgeBoundary.hpp
@@ -0,0 +1,46 @@
+#ifndef MESH_FLAT_EDGE_BOUNDARY_HPP
+#define MESH_FLAT_EDGE_BOUNDARY_HPP
+
+#include <mesh/MeshEdgeBoundary.hpp>
+
+template <size_t Dimension>
+class MeshFlatEdgeBoundary final : public MeshEdgeBoundary<Dimension> // clazy:exclude=copyable-polymorphic
+{
+ public:
+ using Rd = TinyVector<Dimension, double>;
+
+ private:
+ const Rd m_outgoing_normal;
+
+ public:
+ const Rd&
+ outgoingNormal() const
+ {
+ return m_outgoing_normal;
+ }
+
+ MeshFlatEdgeBoundary& operator=(const MeshFlatEdgeBoundary&) = default;
+ MeshFlatEdgeBoundary& operator=(MeshFlatEdgeBoundary&&) = default;
+
+ template <size_t MeshDimension>
+ friend MeshFlatEdgeBoundary<MeshDimension> getMeshFlatEdgeBoundary(const Mesh<Connectivity<MeshDimension>>& mesh,
+ const IBoundaryDescriptor& boundary_descriptor);
+
+ private:
+ template <typename MeshType>
+ MeshFlatEdgeBoundary(const MeshType& mesh, const RefEdgeList& ref_edge_list, const Rd& outgoing_normal)
+ : MeshEdgeBoundary<Dimension>(mesh, ref_edge_list), m_outgoing_normal(outgoing_normal)
+ {}
+
+ public:
+ MeshFlatEdgeBoundary() = default;
+ MeshFlatEdgeBoundary(const MeshFlatEdgeBoundary&) = default;
+ MeshFlatEdgeBoundary(MeshFlatEdgeBoundary&&) = default;
+ ~MeshFlatEdgeBoundary() = default;
+};
+
+template <size_t Dimension>
+MeshFlatEdgeBoundary<Dimension> getMeshFlatEdgeBoundary(const Mesh<Connectivity<Dimension>>& mesh,
+ const IBoundaryDescriptor& boundary_descriptor);
+
+#endif // MESH_FLAT_EDGE_BOUNDARY_HPP
diff --git a/tests/CMakeLists.txt b/tests/CMakeLists.txt
index dfd6616c3d5f0a56b363b8e246c3e771566d5343..0b2849df2e0da1394b80800c150ded28b4cf1ada 100644
--- a/tests/CMakeLists.txt
+++ b/tests/CMakeLists.txt
@@ -177,6 +177,7 @@ add_executable (mpi_unit_tests
test_ItemValueUtils.cpp
test_MeshEdgeBoundary.cpp
test_MeshFaceBoundary.cpp
+ test_MeshFlatEdgeBoundary.cpp
test_MeshFlatFaceBoundary.cpp
test_MeshFlatNodeBoundary.cpp
test_MeshLineFaceBoundary.cpp
diff --git a/tests/test_MeshFlatEdgeBoundary.cpp b/tests/test_MeshFlatEdgeBoundary.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..aebe2e01eafaacfa1d387ab3d8792f2edd4abde3
--- /dev/null
+++ b/tests/test_MeshFlatEdgeBoundary.cpp
@@ -0,0 +1,1471 @@
+#include <catch2/catch_test_macros.hpp>
+#include <catch2/matchers/catch_matchers_all.hpp>
+
+#include <MeshDataBaseForTests.hpp>
+
+#include <algebra/TinyMatrix.hpp>
+#include <mesh/Connectivity.hpp>
+#include <mesh/Mesh.hpp>
+#include <mesh/MeshFlatEdgeBoundary.hpp>
+#include <mesh/NamedBoundaryDescriptor.hpp>
+#include <mesh/NumberedBoundaryDescriptor.hpp>
+
+// clazy:excludeall=non-pod-global-static
+
+TEST_CASE("MeshFlatEdgeBoundary", "[mesh]")
+{
+ auto is_same = [](const auto& a, const auto& b) -> bool {
+ if (a.size() > 0 and b.size() > 0) {
+ return (a[0] == b[0]);
+ } else {
+ return (a.size() == b.size());
+ }
+ };
+
+ auto get_edge_list_from_tag = [](const size_t tag, const auto& connectivity) -> Array<const EdgeId> {
+ for (size_t i = 0; i < connectivity.template numberOfRefItemList<ItemType::edge>(); ++i) {
+ const auto& ref_edge_list = connectivity.template refItemList<ItemType::edge>(i);
+ const RefId ref_id = ref_edge_list.refId();
+ if (ref_id.tagNumber() == tag) {
+ return ref_edge_list.list();
+ }
+ }
+ return {};
+ };
+
+ auto get_edge_list_from_name = [](const std::string& name, const auto& connectivity) -> Array<const EdgeId> {
+ for (size_t i = 0; i < connectivity.template numberOfRefItemList<ItemType::edge>(); ++i) {
+ const auto& ref_edge_list = connectivity.template refItemList<ItemType::edge>(i);
+ const RefId ref_id = ref_edge_list.refId();
+ if (ref_id.tagName() == name) {
+ return ref_edge_list.list();
+ }
+ }
+ return {};
+ };
+
+ SECTION("aligned axis")
+ {
+ SECTION("1D")
+ {
+ static constexpr size_t Dimension = 1;
+
+ using ConnectivityType = Connectivity<Dimension>;
+ using MeshType = Mesh<ConnectivityType>;
+
+ using R1 = TinyVector<1>;
+
+ SECTION("cartesian 1d")
+ {
+ std::shared_ptr p_mesh = MeshDataBaseForTests::get().cartesian1DMesh();
+ const MeshType& mesh = *p_mesh;
+
+ const ConnectivityType& connectivity = mesh.connectivity();
+
+ {
+ const std::set<size_t> tag_set = {0, 1};
+
+ for (auto tag : tag_set) {
+ NumberedBoundaryDescriptor numbered_boundary_descriptor(tag);
+ const auto& edge_boundary = getMeshFlatEdgeBoundary(mesh, numbered_boundary_descriptor);
+
+ auto edge_list = get_edge_list_from_tag(tag, connectivity);
+ REQUIRE(is_same(edge_boundary.edgeList(), edge_list));
+
+ R1 normal = zero;
+
+ switch (tag) {
+ case 0: {
+ normal = R1{-1};
+ break;
+ }
+ case 1: {
+ normal = R1{1};
+ break;
+ }
+ default: {
+ FAIL("unexpected tag number");
+ }
+ }
+
+ REQUIRE(l2Norm(edge_boundary.outgoingNormal() - normal) == Catch::Approx(0).margin(1E-13));
+ }
+ }
+
+ {
+ const std::set<std::string> name_set = {"XMIN", "XMAX"};
+
+ for (auto name : name_set) {
+ NamedBoundaryDescriptor named_boundary_descriptor(name);
+ const auto& edge_boundary = getMeshFlatEdgeBoundary(mesh, named_boundary_descriptor);
+
+ auto edge_list = get_edge_list_from_name(name, connectivity);
+
+ REQUIRE(is_same(edge_boundary.edgeList(), edge_list));
+
+ R1 normal = zero;
+
+ if (name == "XMIN") {
+ normal = R1{-1};
+ } else if (name == "XMAX") {
+ normal = R1{1};
+ } else {
+ FAIL("unexpected name: " + name);
+ }
+
+ REQUIRE(l2Norm(edge_boundary.outgoingNormal() - normal) == Catch::Approx(0).margin(1E-13));
+ }
+ }
+ }
+
+ SECTION("unordered 1d")
+ {
+ std::shared_ptr p_mesh = MeshDataBaseForTests::get().unordered1DMesh();
+ const MeshType& mesh = *p_mesh;
+
+ const ConnectivityType& connectivity = mesh.connectivity();
+
+ {
+ const std::set<size_t> tag_set = {1, 2};
+
+ for (auto tag : tag_set) {
+ NumberedBoundaryDescriptor numbered_boundary_descriptor(tag);
+ const auto& edge_boundary = getMeshFlatEdgeBoundary(mesh, numbered_boundary_descriptor);
+
+ auto edge_list = get_edge_list_from_tag(tag, connectivity);
+ REQUIRE(is_same(edge_boundary.edgeList(), edge_list));
+
+ R1 normal = zero;
+
+ switch (tag) {
+ case 1: {
+ normal = R1{-1};
+ break;
+ }
+ case 2: {
+ normal = R1{1};
+ break;
+ }
+ default: {
+ FAIL("unexpected tag number");
+ }
+ }
+
+ REQUIRE(l2Norm(edge_boundary.outgoingNormal() - normal) == Catch::Approx(0).margin(1E-13));
+ }
+ }
+
+ {
+ const std::set<std::string> name_set = {"XMIN", "XMAX"};
+
+ for (auto name : name_set) {
+ NamedBoundaryDescriptor named_boundary_descriptor(name);
+ const auto& edge_boundary = getMeshFlatEdgeBoundary(mesh, named_boundary_descriptor);
+
+ auto edge_list = get_edge_list_from_name(name, connectivity);
+
+ REQUIRE(is_same(edge_boundary.edgeList(), edge_list));
+
+ R1 normal = zero;
+
+ if (name == "XMIN") {
+ normal = R1{-1};
+ } else if (name == "XMAX") {
+ normal = R1{1};
+ } else {
+ FAIL("unexpected name: " + name);
+ }
+
+ REQUIRE(l2Norm(edge_boundary.outgoingNormal() - normal) == Catch::Approx(0).margin(1E-13));
+ }
+ }
+ }
+ }
+
+ SECTION("2D")
+ {
+ static constexpr size_t Dimension = 2;
+
+ using ConnectivityType = Connectivity<Dimension>;
+ using MeshType = Mesh<ConnectivityType>;
+
+ using R2 = TinyVector<2>;
+
+ SECTION("cartesian 2d")
+ {
+ std::shared_ptr p_mesh = MeshDataBaseForTests::get().cartesian2DMesh();
+ const MeshType& mesh = *p_mesh;
+
+ const ConnectivityType& connectivity = mesh.connectivity();
+
+ {
+ const std::set<size_t> tag_set = {0, 1, 2, 3};
+
+ for (auto tag : tag_set) {
+ NumberedBoundaryDescriptor numbered_boundary_descriptor(tag);
+ const auto& edge_boundary = getMeshFlatEdgeBoundary(mesh, numbered_boundary_descriptor);
+
+ auto edge_list = get_edge_list_from_tag(tag, connectivity);
+ REQUIRE(is_same(edge_boundary.edgeList(), edge_list));
+
+ R2 normal = zero;
+
+ switch (tag) {
+ case 0: {
+ normal = R2{-1, 0};
+ break;
+ }
+ case 1: {
+ normal = R2{1, 0};
+ break;
+ }
+ case 2: {
+ normal = R2{0, -1};
+ break;
+ }
+ case 3: {
+ normal = R2{0, 1};
+ break;
+ }
+ default: {
+ FAIL("unexpected tag number");
+ }
+ }
+
+ REQUIRE(l2Norm(edge_boundary.outgoingNormal() - normal) == Catch::Approx(0).margin(1E-13));
+ }
+
+ {
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NumberedBoundaryDescriptor(10)),
+ "error: cannot find edge list with name \"10\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NumberedBoundaryDescriptor(11)),
+ "error: cannot find edge list with name \"11\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NumberedBoundaryDescriptor(12)),
+ "error: cannot find edge list with name \"12\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NumberedBoundaryDescriptor(13)),
+ "error: cannot find edge list with name \"13\"");
+ }
+ }
+
+ {
+ const std::set<std::string> name_set = {"XMIN", "XMAX", "YMIN", "YMAX"};
+
+ for (auto name : name_set) {
+ NamedBoundaryDescriptor named_boundary_descriptor(name);
+ const auto& edge_boundary = getMeshFlatEdgeBoundary(mesh, named_boundary_descriptor);
+
+ auto edge_list = get_edge_list_from_name(name, connectivity);
+ REQUIRE(is_same(edge_boundary.edgeList(), edge_list));
+
+ R2 normal = zero;
+
+ if (name == "XMIN") {
+ normal = R2{-1, 0};
+ } else if (name == "XMAX") {
+ normal = R2{1, 0};
+ } else if (name == "YMIN") {
+ normal = R2{0, -1};
+ } else if (name == "YMAX") {
+ normal = R2{0, 1};
+ } else {
+ FAIL("unexpected name: " + name);
+ }
+
+ REQUIRE(l2Norm(edge_boundary.outgoingNormal() - normal) == Catch::Approx(0).margin(1E-13));
+ }
+
+ {
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NamedBoundaryDescriptor("XMINYMIN")),
+ "error: cannot find edge list with name \"XMINYMIN\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NamedBoundaryDescriptor("XMINYMAX")),
+ "error: cannot find edge list with name \"XMINYMAX\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NamedBoundaryDescriptor("XMAXYMIN")),
+ "error: cannot find edge list with name \"XMAXYMIN\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NamedBoundaryDescriptor("XMAXYMAX")),
+ "error: cannot find edge list with name \"XMAXYMAX\"");
+ }
+ }
+ }
+
+ SECTION("hybrid 2d")
+ {
+ std::shared_ptr p_mesh = MeshDataBaseForTests::get().hybrid2DMesh();
+ const MeshType& mesh = *p_mesh;
+
+ const ConnectivityType& connectivity = mesh.connectivity();
+
+ {
+ const std::set<size_t> tag_set = {1, 2, 3, 4};
+
+ for (auto tag : tag_set) {
+ NumberedBoundaryDescriptor numbered_boundary_descriptor(tag);
+ const auto& edge_boundary = getMeshFlatEdgeBoundary(mesh, numbered_boundary_descriptor);
+
+ auto edge_list = get_edge_list_from_tag(tag, connectivity);
+ REQUIRE(is_same(edge_boundary.edgeList(), edge_list));
+
+ R2 normal = zero;
+
+ switch (tag) {
+ case 1: {
+ normal = R2{-1, 0};
+ break;
+ }
+ case 2: {
+ normal = R2{1, 0};
+ break;
+ }
+ case 3: {
+ normal = R2{0, 1};
+ break;
+ }
+ case 4: {
+ normal = R2{0, -1};
+ break;
+ }
+ default: {
+ FAIL("unexpected tag number");
+ }
+ }
+
+ REQUIRE(l2Norm(edge_boundary.outgoingNormal() - normal) == Catch::Approx(0).margin(1E-13));
+ }
+
+ {
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NumberedBoundaryDescriptor(8)),
+ "error: cannot find edge list with name \"8\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NumberedBoundaryDescriptor(9)),
+ "error: cannot find edge list with name \"9\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NumberedBoundaryDescriptor(10)),
+ "error: cannot find edge list with name \"10\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NumberedBoundaryDescriptor(11)),
+ "error: cannot find edge list with name \"11\"");
+ }
+ }
+
+ {
+ const std::set<std::string> name_set = {"XMIN", "XMAX", "YMIN", "YMAX"};
+
+ for (auto name : name_set) {
+ NamedBoundaryDescriptor named_boundary_descriptor(name);
+ const auto& edge_boundary = getMeshFlatEdgeBoundary(mesh, named_boundary_descriptor);
+
+ auto edge_list = get_edge_list_from_name(name, connectivity);
+
+ R2 normal = zero;
+
+ REQUIRE(is_same(edge_boundary.edgeList(), edge_list));
+
+ if (name == "XMIN") {
+ normal = R2{-1, 0};
+ } else if (name == "XMAX") {
+ normal = R2{1, 0};
+ } else if (name == "YMIN") {
+ normal = R2{0, -1};
+ } else if (name == "YMAX") {
+ normal = R2{0, 1};
+ } else {
+ FAIL("unexpected name: " + name);
+ }
+
+ REQUIRE(l2Norm(edge_boundary.outgoingNormal() - normal) == Catch::Approx(0).margin(1E-13));
+ }
+
+ {
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NamedBoundaryDescriptor("XMINYMIN")),
+ "error: cannot find edge list with name \"XMINYMIN\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NamedBoundaryDescriptor("XMINYMAX")),
+ "error: cannot find edge list with name \"XMINYMAX\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NamedBoundaryDescriptor("XMAXYMIN")),
+ "error: cannot find edge list with name \"XMAXYMIN\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NamedBoundaryDescriptor("XMAXYMAX")),
+ "error: cannot find edge list with name \"XMAXYMAX\"");
+ }
+ }
+ }
+ }
+
+ SECTION("3D")
+ {
+ static constexpr size_t Dimension = 3;
+
+ using ConnectivityType = Connectivity<Dimension>;
+ using MeshType = Mesh<ConnectivityType>;
+
+ using R3 = TinyVector<3>;
+
+ SECTION("cartesian 3d")
+ {
+ std::shared_ptr p_mesh = MeshDataBaseForTests::get().cartesian3DMesh();
+ const MeshType& mesh = *p_mesh;
+
+ const ConnectivityType& connectivity = mesh.connectivity();
+
+ {
+ 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& edge_boundary = getMeshFlatEdgeBoundary(mesh, numbered_boundary_descriptor);
+
+ auto edge_list = get_edge_list_from_tag(tag, connectivity);
+ REQUIRE(is_same(edge_boundary.edgeList(), edge_list));
+
+ R3 normal = zero;
+
+ switch (tag) {
+ case 0: {
+ normal = R3{-1, 0, 0};
+ break;
+ }
+ case 1: {
+ normal = R3{1, 0, 0};
+ break;
+ }
+ case 2: {
+ normal = R3{0, -1, 0};
+ break;
+ }
+ case 3: {
+ normal = R3{0, 1, 0};
+ break;
+ }
+ case 4: {
+ normal = R3{0, 0, -1};
+ break;
+ }
+ case 5: {
+ normal = R3{0, 0, 1};
+ break;
+ }
+ default: {
+ FAIL("unexpected tag number");
+ }
+ }
+
+ REQUIRE(l2Norm(edge_boundary.outgoingNormal() - normal) == Catch::Approx(0).margin(1E-13));
+ }
+
+ {
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NumberedBoundaryDescriptor(20)),
+ "error: cannot find surface with name \"20\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NumberedBoundaryDescriptor(21)),
+ "error: cannot find surface with name \"21\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NumberedBoundaryDescriptor(22)),
+ "error: cannot find surface with name \"22\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NumberedBoundaryDescriptor(23)),
+ "error: cannot find surface with name \"23\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NumberedBoundaryDescriptor(24)),
+ "error: cannot find surface with name \"24\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NumberedBoundaryDescriptor(25)),
+ "error: cannot find surface with name \"25\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NumberedBoundaryDescriptor(26)),
+ "error: cannot find surface with name \"26\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NumberedBoundaryDescriptor(27)),
+ "error: cannot find surface with name \"27\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NumberedBoundaryDescriptor(28)),
+ "error: cannot find surface with name \"28\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NumberedBoundaryDescriptor(29)),
+ "error: cannot find surface with name \"29\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NumberedBoundaryDescriptor(30)),
+ "error: cannot find surface with name \"30\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NumberedBoundaryDescriptor(31)),
+ "error: cannot find surface with name \"31\"");
+
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NumberedBoundaryDescriptor(10)),
+ "error: cannot find edge list with name \"10\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NumberedBoundaryDescriptor(11)),
+ "error: cannot find edge list with name \"11\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NumberedBoundaryDescriptor(12)),
+ "error: cannot find edge list with name \"12\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NumberedBoundaryDescriptor(13)),
+ "error: cannot find edge list with name \"13\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NumberedBoundaryDescriptor(14)),
+ "error: cannot find edge list with name \"14\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NumberedBoundaryDescriptor(15)),
+ "error: cannot find edge list with name \"15\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NumberedBoundaryDescriptor(16)),
+ "error: cannot find edge list with name \"16\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NumberedBoundaryDescriptor(17)),
+ "error: cannot find edge list with name \"17\"");
+ }
+ }
+
+ {
+ const std::set<std::string> name_set = {"XMIN", "XMAX", "YMIN", "YMAX", "ZMIN", "ZMAX"};
+
+ for (auto name : name_set) {
+ NamedBoundaryDescriptor named_boundary_descriptor(name);
+ const auto& edge_boundary = getMeshFlatEdgeBoundary(mesh, named_boundary_descriptor);
+
+ auto edge_list = get_edge_list_from_name(name, connectivity);
+
+ REQUIRE(is_same(edge_boundary.edgeList(), edge_list));
+
+ R3 normal = zero;
+
+ if (name == "XMIN") {
+ normal = R3{-1, 0, 0};
+ } else if (name == "XMAX") {
+ normal = R3{1, 0, 0};
+ } else if (name == "YMIN") {
+ normal = R3{0, -1, 0};
+ } else if (name == "YMAX") {
+ normal = R3{0, 1, 0};
+ } else if (name == "ZMIN") {
+ normal = R3{0, 0, -1};
+ } else if (name == "ZMAX") {
+ normal = R3{0, 0, 1};
+ } else {
+ FAIL("unexpected name: " + name);
+ }
+
+ REQUIRE(l2Norm(edge_boundary.outgoingNormal() - normal) == Catch::Approx(0).margin(1E-13));
+ }
+
+ {
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NamedBoundaryDescriptor("XMINYMIN")),
+ "error: cannot find surface with name \"XMINYMIN\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NamedBoundaryDescriptor("XMINYMAX")),
+ "error: cannot find surface with name \"XMINYMAX\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NamedBoundaryDescriptor("XMAXYMAX")),
+ "error: cannot find surface with name \"XMAXYMAX\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NamedBoundaryDescriptor("XMINYMAX")),
+ "error: cannot find surface with name \"XMINYMAX\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NamedBoundaryDescriptor("XMINZMIN")),
+ "error: cannot find surface with name \"XMINZMIN\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NamedBoundaryDescriptor("XMINZMAX")),
+ "error: cannot find surface with name \"XMINZMAX\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NamedBoundaryDescriptor("XMAXZMAX")),
+ "error: cannot find surface with name \"XMAXZMAX\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NamedBoundaryDescriptor("XMINZMAX")),
+ "error: cannot find surface with name \"XMINZMAX\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NamedBoundaryDescriptor("YMINZMIN")),
+ "error: cannot find surface with name \"YMINZMIN\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NamedBoundaryDescriptor("YMINZMAX")),
+ "error: cannot find surface with name \"YMINZMAX\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NamedBoundaryDescriptor("YMAXZMAX")),
+ "error: cannot find surface with name \"YMAXZMAX\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NamedBoundaryDescriptor("YMINZMAX")),
+ "error: cannot find surface with name \"YMINZMAX\"");
+
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NamedBoundaryDescriptor("XMINYMINZMIN")),
+ "error: cannot find edge list with name \"XMINYMINZMIN\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NamedBoundaryDescriptor("XMAXYMINZMIN")),
+ "error: cannot find edge list with name \"XMAXYMINZMIN\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NamedBoundaryDescriptor("XMAXYMAXZMIN")),
+ "error: cannot find edge list with name \"XMAXYMAXZMIN\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NamedBoundaryDescriptor("XMINYMAXZMIN")),
+ "error: cannot find edge list with name \"XMINYMAXZMIN\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NamedBoundaryDescriptor("XMINYMINZMAX")),
+ "error: cannot find edge list with name \"XMINYMINZMAX\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NamedBoundaryDescriptor("XMAXYMINZMAX")),
+ "error: cannot find edge list with name \"XMAXYMINZMAX\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NamedBoundaryDescriptor("XMAXYMAXZMAX")),
+ "error: cannot find edge list with name \"XMAXYMAXZMAX\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NamedBoundaryDescriptor("XMINYMAXZMAX")),
+ "error: cannot find edge list with name \"XMINYMAXZMAX\"");
+ }
+ }
+ }
+
+ SECTION("hybrid 3d")
+ {
+ std::shared_ptr p_mesh = MeshDataBaseForTests::get().hybrid3DMesh();
+ const MeshType& mesh = *p_mesh;
+
+ const ConnectivityType& connectivity = mesh.connectivity();
+
+ {
+ const std::set<size_t> tag_set = {22, 23, 24, 25, 26, 27};
+
+ for (auto tag : tag_set) {
+ NumberedBoundaryDescriptor numbered_boundary_descriptor(tag);
+ const auto& edge_boundary = getMeshFlatEdgeBoundary(mesh, numbered_boundary_descriptor);
+
+ auto edge_list = get_edge_list_from_tag(tag, connectivity);
+ REQUIRE(is_same(edge_boundary.edgeList(), edge_list));
+
+ R3 normal = zero;
+
+ switch (tag) {
+ case 22: {
+ normal = R3{-1, 0, 0};
+ break;
+ }
+ case 23: {
+ normal = R3{1, 0, 0};
+ break;
+ }
+ case 24: {
+ normal = R3{0, 0, 1};
+ break;
+ }
+ case 25: {
+ normal = R3{0, 0, -1};
+ break;
+ }
+ case 26: {
+ normal = R3{0, 1, 0};
+ break;
+ }
+ case 27: {
+ normal = R3{0, -1, 0};
+ break;
+ }
+ default: {
+ FAIL("unexpected tag number");
+ }
+ }
+
+ REQUIRE(l2Norm(edge_boundary.outgoingNormal() - normal) == Catch::Approx(0).margin(1E-13));
+ }
+
+ {
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NumberedBoundaryDescriptor(28)),
+ "error: cannot find surface with name \"28\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NumberedBoundaryDescriptor(29)),
+ "error: cannot find surface with name \"29\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NumberedBoundaryDescriptor(30)),
+ "error: cannot find surface with name \"30\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NumberedBoundaryDescriptor(31)),
+ "error: cannot find surface with name \"31\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NumberedBoundaryDescriptor(32)),
+ "error: cannot find surface with name \"32\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NumberedBoundaryDescriptor(33)),
+ "error: cannot find surface with name \"33\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NumberedBoundaryDescriptor(34)),
+ "error: cannot find surface with name \"34\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NumberedBoundaryDescriptor(35)),
+ "error: cannot find surface with name \"35\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NumberedBoundaryDescriptor(36)),
+ "error: cannot find surface with name \"36\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NumberedBoundaryDescriptor(37)),
+ "error: cannot find surface with name \"37\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NumberedBoundaryDescriptor(38)),
+ "error: cannot find surface with name \"38\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NumberedBoundaryDescriptor(39)),
+ "error: cannot find surface with name \"39\"");
+
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NumberedBoundaryDescriptor(40)),
+ "error: cannot find edge list with name \"40\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NumberedBoundaryDescriptor(41)),
+ "error: cannot find edge list with name \"41\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NumberedBoundaryDescriptor(42)),
+ "error: cannot find edge list with name \"42\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NumberedBoundaryDescriptor(43)),
+ "error: cannot find edge list with name \"43\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NumberedBoundaryDescriptor(44)),
+ "error: cannot find edge list with name \"44\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NumberedBoundaryDescriptor(45)),
+ "error: cannot find edge list with name \"45\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NumberedBoundaryDescriptor(47)),
+ "error: cannot find edge list with name \"47\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NumberedBoundaryDescriptor(51)),
+ "error: cannot find edge list with name \"51\"");
+ }
+ }
+
+ {
+ const std::set<std::string> name_set = {"XMIN", "XMAX", "YMIN", "YMAX", "ZMIN", "ZMAX"};
+
+ for (auto name : name_set) {
+ NamedBoundaryDescriptor named_boundary_descriptor(name);
+ const auto& edge_boundary = getMeshFlatEdgeBoundary(mesh, named_boundary_descriptor);
+
+ auto edge_list = get_edge_list_from_name(name, connectivity);
+
+ REQUIRE(is_same(edge_boundary.edgeList(), edge_list));
+
+ R3 normal = zero;
+
+ if (name == "XMIN") {
+ normal = R3{-1, 0, 0};
+ } else if (name == "XMAX") {
+ normal = R3{1, 0, 0};
+ } else if (name == "YMIN") {
+ normal = R3{0, -1, 0};
+ } else if (name == "YMAX") {
+ normal = R3{0, 1, 0};
+ } else if (name == "ZMIN") {
+ normal = R3{0, 0, -1};
+ } else if (name == "ZMAX") {
+ normal = R3{0, 0, 1};
+ } else {
+ FAIL("unexpected name: " + name);
+ }
+
+ REQUIRE(l2Norm(edge_boundary.outgoingNormal() - normal) == Catch::Approx(0).margin(1E-13));
+ }
+
+ {
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NamedBoundaryDescriptor("XMINYMIN")),
+ "error: cannot find surface with name \"XMINYMIN\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NamedBoundaryDescriptor("XMINYMAX")),
+ "error: cannot find surface with name \"XMINYMAX\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NamedBoundaryDescriptor("XMAXYMAX")),
+ "error: cannot find surface with name \"XMAXYMAX\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NamedBoundaryDescriptor("XMINYMAX")),
+ "error: cannot find surface with name \"XMINYMAX\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NamedBoundaryDescriptor("XMINZMIN")),
+ "error: cannot find surface with name \"XMINZMIN\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NamedBoundaryDescriptor("XMINZMAX")),
+ "error: cannot find surface with name \"XMINZMAX\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NamedBoundaryDescriptor("XMAXZMAX")),
+ "error: cannot find surface with name \"XMAXZMAX\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NamedBoundaryDescriptor("XMINZMAX")),
+ "error: cannot find surface with name \"XMINZMAX\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NamedBoundaryDescriptor("YMINZMIN")),
+ "error: cannot find surface with name \"YMINZMIN\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NamedBoundaryDescriptor("YMINZMAX")),
+ "error: cannot find surface with name \"YMINZMAX\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NamedBoundaryDescriptor("YMAXZMAX")),
+ "error: cannot find surface with name \"YMAXZMAX\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NamedBoundaryDescriptor("YMINZMAX")),
+ "error: cannot find surface with name \"YMINZMAX\"");
+
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NamedBoundaryDescriptor("XMINYMINZMIN")),
+ "error: cannot find edge list with name \"XMINYMINZMIN\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NamedBoundaryDescriptor("XMAXYMINZMIN")),
+ "error: cannot find edge list with name \"XMAXYMINZMIN\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NamedBoundaryDescriptor("XMAXYMAXZMIN")),
+ "error: cannot find edge list with name \"XMAXYMAXZMIN\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NamedBoundaryDescriptor("XMINYMAXZMIN")),
+ "error: cannot find edge list with name \"XMINYMAXZMIN\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NamedBoundaryDescriptor("XMINYMINZMAX")),
+ "error: cannot find edge list with name \"XMINYMINZMAX\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NamedBoundaryDescriptor("XMAXYMINZMAX")),
+ "error: cannot find edge list with name \"XMAXYMINZMAX\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NamedBoundaryDescriptor("XMAXYMAXZMAX")),
+ "error: cannot find edge list with name \"XMAXYMAXZMAX\"");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, NamedBoundaryDescriptor("XMINYMAXZMAX")),
+ "error: cannot find edge list with name \"XMINYMAXZMAX\"");
+ }
+ }
+ }
+ }
+ }
+
+ SECTION("rotated axis")
+ {
+ SECTION("2D")
+ {
+ static constexpr size_t Dimension = 2;
+
+ using ConnectivityType = Connectivity<Dimension>;
+ using MeshType = Mesh<ConnectivityType>;
+
+ using R2 = TinyVector<2>;
+
+ const double theta = 0.3;
+ const TinyMatrix<2> R{std::cos(theta), -std::sin(theta), //
+ std::sin(theta), std::cos(theta)};
+
+ SECTION("cartesian 2d")
+ {
+ std::shared_ptr p_mesh = MeshDataBaseForTests::get().cartesian2DMesh();
+
+ const ConnectivityType& connectivity = p_mesh->connectivity();
+
+ auto xr = p_mesh->xr();
+
+ NodeValue<R2> rotated_xr{connectivity};
+
+ parallel_for(
+ connectivity.numberOfNodes(), PUGS_LAMBDA(const NodeId node_id) { rotated_xr[node_id] = R * xr[node_id]; });
+
+ MeshType mesh{p_mesh->shared_connectivity(), rotated_xr};
+
+ {
+ const std::set<size_t> tag_set = {0, 1, 2, 3};
+
+ for (auto tag : tag_set) {
+ NumberedBoundaryDescriptor numbered_boundary_descriptor(tag);
+ const auto& edge_boundary = getMeshFlatEdgeBoundary(mesh, numbered_boundary_descriptor);
+
+ auto edge_list = get_edge_list_from_tag(tag, connectivity);
+ REQUIRE(is_same(edge_boundary.edgeList(), edge_list));
+
+ R2 normal = zero;
+
+ switch (tag) {
+ case 0: {
+ normal = R * R2{-1, 0};
+ break;
+ }
+ case 1: {
+ normal = R * R2{1, 0};
+ break;
+ }
+ case 2: {
+ normal = R * R2{0, -1};
+ break;
+ }
+ case 3: {
+ normal = R * R2{0, 1};
+ break;
+ }
+ default: {
+ FAIL("unexpected tag number");
+ }
+ }
+
+ REQUIRE(l2Norm(edge_boundary.outgoingNormal() - normal) == Catch::Approx(0).margin(1E-13));
+ }
+ }
+
+ {
+ const std::set<std::string> name_set = {"XMIN", "XMAX", "YMIN", "YMAX"};
+
+ for (auto name : name_set) {
+ NamedBoundaryDescriptor named_boundary_descriptor(name);
+ const auto& edge_boundary = getMeshFlatEdgeBoundary(mesh, named_boundary_descriptor);
+
+ auto edge_list = get_edge_list_from_name(name, connectivity);
+ REQUIRE(is_same(edge_boundary.edgeList(), edge_list));
+
+ R2 normal = zero;
+
+ if (name == "XMIN") {
+ normal = R * R2{-1, 0};
+ } else if (name == "XMAX") {
+ normal = R * R2{1, 0};
+ } else if (name == "YMIN") {
+ normal = R * R2{0, -1};
+ } else if (name == "YMAX") {
+ normal = R * R2{0, 1};
+ } else {
+ FAIL("unexpected name: " + name);
+ }
+
+ REQUIRE(l2Norm(edge_boundary.outgoingNormal() - normal) == Catch::Approx(0).margin(1E-13));
+ }
+ }
+ }
+
+ SECTION("hybrid 2d")
+ {
+ std::shared_ptr p_mesh = MeshDataBaseForTests::get().hybrid2DMesh();
+
+ const ConnectivityType& connectivity = p_mesh->connectivity();
+
+ auto xr = p_mesh->xr();
+
+ NodeValue<R2> rotated_xr{connectivity};
+ parallel_for(
+ connectivity.numberOfNodes(), PUGS_LAMBDA(const NodeId node_id) { rotated_xr[node_id] = R * xr[node_id]; });
+
+ MeshType mesh{p_mesh->shared_connectivity(), rotated_xr};
+
+ {
+ const std::set<size_t> tag_set = {1, 2, 3, 4};
+
+ for (auto tag : tag_set) {
+ NumberedBoundaryDescriptor numbered_boundary_descriptor(tag);
+ const auto& edge_boundary = getMeshFlatEdgeBoundary(mesh, numbered_boundary_descriptor);
+
+ auto edge_list = get_edge_list_from_tag(tag, connectivity);
+ REQUIRE(is_same(edge_boundary.edgeList(), edge_list));
+
+ R2 normal = zero;
+
+ switch (tag) {
+ case 1: {
+ normal = R * R2{-1, 0};
+ break;
+ }
+ case 2: {
+ normal = R * R2{1, 0};
+ break;
+ }
+ case 3: {
+ normal = R * R2{0, 1};
+ break;
+ }
+ case 4: {
+ normal = R * R2{0, -1};
+ break;
+ }
+ default: {
+ FAIL("unexpected tag number");
+ }
+ }
+
+ REQUIRE(l2Norm(edge_boundary.outgoingNormal() - normal) == Catch::Approx(0).margin(1E-13));
+ }
+ }
+
+ {
+ const std::set<std::string> name_set = {"XMIN", "XMAX", "YMIN", "YMAX"};
+
+ for (auto name : name_set) {
+ NamedBoundaryDescriptor named_boundary_descriptor(name);
+ const auto& edge_boundary = getMeshFlatEdgeBoundary(mesh, named_boundary_descriptor);
+
+ auto edge_list = get_edge_list_from_name(name, connectivity);
+
+ R2 normal = zero;
+
+ REQUIRE(is_same(edge_boundary.edgeList(), edge_list));
+
+ if (name == "XMIN") {
+ normal = R * R2{-1, 0};
+ } else if (name == "XMAX") {
+ normal = R * R2{1, 0};
+ } else if (name == "YMIN") {
+ normal = R * R2{0, -1};
+ } else if (name == "YMAX") {
+ normal = R * R2{0, 1};
+ } else {
+ FAIL("unexpected name: " + name);
+ }
+
+ REQUIRE(l2Norm(edge_boundary.outgoingNormal() - normal) == Catch::Approx(0).margin(1E-13));
+ }
+ }
+ }
+ }
+
+ SECTION("3D")
+ {
+ static constexpr size_t Dimension = 3;
+
+ using ConnectivityType = Connectivity<Dimension>;
+ using MeshType = Mesh<ConnectivityType>;
+
+ using R3 = TinyVector<3>;
+
+ const double theta = 0.3;
+ const double phi = 0.4;
+ const TinyMatrix<3> R =
+ TinyMatrix<3>{std::cos(theta), -std::sin(theta), 0, std::sin(theta), std::cos(theta), 0, 0, 0, 1} *
+ TinyMatrix<3>{0, std::cos(phi), -std::sin(phi), 0, std::sin(phi), std::cos(phi), 1, 0, 0};
+
+ SECTION("cartesian 3d")
+ {
+ std::shared_ptr p_mesh = MeshDataBaseForTests::get().cartesian3DMesh();
+
+ const ConnectivityType& connectivity = p_mesh->connectivity();
+
+ auto xr = p_mesh->xr();
+
+ NodeValue<R3> rotated_xr{connectivity};
+
+ parallel_for(
+ connectivity.numberOfNodes(), PUGS_LAMBDA(const NodeId node_id) { rotated_xr[node_id] = R * xr[node_id]; });
+
+ MeshType mesh{p_mesh->shared_connectivity(), rotated_xr};
+
+ {
+ 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& edge_boundary = getMeshFlatEdgeBoundary(mesh, numbered_boundary_descriptor);
+
+ auto edge_list = get_edge_list_from_tag(tag, connectivity);
+ REQUIRE(is_same(edge_boundary.edgeList(), edge_list));
+
+ R3 normal = zero;
+
+ switch (tag) {
+ case 0: {
+ normal = R * R3{-1, 0, 0};
+ break;
+ }
+ case 1: {
+ normal = R * R3{1, 0, 0};
+ break;
+ }
+ case 2: {
+ normal = R * R3{0, -1, 0};
+ break;
+ }
+ case 3: {
+ normal = R * R3{0, 1, 0};
+ break;
+ }
+ case 4: {
+ normal = R * R3{0, 0, -1};
+ break;
+ }
+ case 5: {
+ normal = R * R3{0, 0, 1};
+ break;
+ }
+ default: {
+ FAIL("unexpected tag number");
+ }
+ }
+
+ REQUIRE(l2Norm(edge_boundary.outgoingNormal() - normal) == Catch::Approx(0).margin(1E-13));
+ }
+ }
+
+ {
+ const std::set<std::string> name_set = {"XMIN", "XMAX", "YMIN", "YMAX", "ZMIN", "ZMAX"};
+
+ for (auto name : name_set) {
+ NamedBoundaryDescriptor named_boundary_descriptor(name);
+ const auto& edge_boundary = getMeshFlatEdgeBoundary(mesh, named_boundary_descriptor);
+
+ auto edge_list = get_edge_list_from_name(name, connectivity);
+
+ REQUIRE(is_same(edge_boundary.edgeList(), edge_list));
+
+ R3 normal = zero;
+
+ if (name == "XMIN") {
+ normal = R * R3{-1, 0, 0};
+ } else if (name == "XMAX") {
+ normal = R * R3{1, 0, 0};
+ } else if (name == "YMIN") {
+ normal = R * R3{0, -1, 0};
+ } else if (name == "YMAX") {
+ normal = R * R3{0, 1, 0};
+ } else if (name == "ZMIN") {
+ normal = R * R3{0, 0, -1};
+ } else if (name == "ZMAX") {
+ normal = R * R3{0, 0, 1};
+ } else {
+ FAIL("unexpected name: " + name);
+ }
+
+ REQUIRE(l2Norm(edge_boundary.outgoingNormal() - normal) == Catch::Approx(0).margin(1E-13));
+ }
+ }
+ }
+
+ SECTION("hybrid 3d")
+ {
+ std::shared_ptr p_mesh = MeshDataBaseForTests::get().hybrid3DMesh();
+
+ const ConnectivityType& connectivity = p_mesh->connectivity();
+
+ auto xr = p_mesh->xr();
+
+ NodeValue<R3> rotated_xr{connectivity};
+
+ parallel_for(
+ connectivity.numberOfNodes(), PUGS_LAMBDA(const NodeId node_id) { rotated_xr[node_id] = R * xr[node_id]; });
+
+ MeshType mesh{p_mesh->shared_connectivity(), rotated_xr};
+
+ {
+ const std::set<size_t> tag_set = {22, 23, 24, 25, 26, 27};
+
+ for (auto tag : tag_set) {
+ NumberedBoundaryDescriptor numbered_boundary_descriptor(tag);
+ const auto& edge_boundary = getMeshFlatEdgeBoundary(mesh, numbered_boundary_descriptor);
+
+ auto edge_list = get_edge_list_from_tag(tag, connectivity);
+ REQUIRE(is_same(edge_boundary.edgeList(), edge_list));
+
+ R3 normal = zero;
+
+ switch (tag) {
+ case 22: {
+ normal = R * R3{-1, 0, 0};
+ break;
+ }
+ case 23: {
+ normal = R * R3{1, 0, 0};
+ break;
+ }
+ case 24: {
+ normal = R * R3{0, 0, 1};
+ break;
+ }
+ case 25: {
+ normal = R * R3{0, 0, -1};
+ break;
+ }
+ case 26: {
+ normal = R * R3{0, 1, 0};
+ break;
+ }
+ case 27: {
+ normal = R * R3{0, -1, 0};
+ break;
+ }
+ default: {
+ FAIL("unexpected tag number");
+ }
+ }
+
+ REQUIRE(l2Norm(edge_boundary.outgoingNormal() - normal) == Catch::Approx(0).margin(1E-13));
+ }
+ }
+
+ {
+ const std::set<std::string> name_set = {"XMIN", "XMAX", "YMIN", "YMAX", "ZMIN", "ZMAX"};
+
+ for (auto name : name_set) {
+ NamedBoundaryDescriptor named_boundary_descriptor(name);
+ const auto& edge_boundary = getMeshFlatEdgeBoundary(mesh, named_boundary_descriptor);
+
+ auto edge_list = get_edge_list_from_name(name, connectivity);
+
+ REQUIRE(is_same(edge_boundary.edgeList(), edge_list));
+
+ R3 normal = zero;
+
+ if (name == "XMIN") {
+ normal = R * R3{-1, 0, 0};
+ } else if (name == "XMAX") {
+ normal = R * R3{1, 0, 0};
+ } else if (name == "YMIN") {
+ normal = R * R3{0, -1, 0};
+ } else if (name == "YMAX") {
+ normal = R * R3{0, 1, 0};
+ } else if (name == "ZMIN") {
+ normal = R * R3{0, 0, -1};
+ } else if (name == "ZMAX") {
+ normal = R * R3{0, 0, 1};
+ } else {
+ FAIL("unexpected name: " + name);
+ }
+
+ REQUIRE(l2Norm(edge_boundary.outgoingNormal() - normal) == Catch::Approx(0).margin(1E-13));
+ }
+ }
+ }
+ }
+ }
+
+ SECTION("curved mesh")
+ {
+ SECTION("2D")
+ {
+ static constexpr size_t Dimension = 2;
+
+ using ConnectivityType = Connectivity<Dimension>;
+ using MeshType = Mesh<ConnectivityType>;
+
+ using R2 = TinyVector<2>;
+
+ auto curve = [](const R2& X) -> R2 { return R2{X[0], (1 + X[0] * X[0]) * X[1]}; };
+
+ SECTION("hybrid 2d")
+ {
+ std::shared_ptr p_mesh = MeshDataBaseForTests::get().hybrid2DMesh();
+
+ const ConnectivityType& connectivity = p_mesh->connectivity();
+
+ auto xr = p_mesh->xr();
+
+ NodeValue<TinyVector<2>> curved_xr{connectivity};
+ parallel_for(
+ connectivity.numberOfNodes(), PUGS_LAMBDA(const NodeId node_id) { curved_xr[node_id] = curve(xr[node_id]); });
+
+ MeshType mesh{p_mesh->shared_connectivity(), curved_xr};
+
+ {
+ const std::set<size_t> tag_set = {1, 2, 4};
+
+ for (auto tag : tag_set) {
+ NumberedBoundaryDescriptor numbered_boundary_descriptor(tag);
+ const auto& edge_boundary = getMeshFlatEdgeBoundary(mesh, numbered_boundary_descriptor);
+
+ auto edge_list = get_edge_list_from_tag(tag, connectivity);
+ REQUIRE(is_same(edge_boundary.edgeList(), edge_list));
+
+ R2 normal = zero;
+
+ switch (tag) {
+ case 1: {
+ normal = R2{-1, 0};
+ break;
+ }
+ case 2: {
+ normal = R2{1, 0};
+ break;
+ }
+ case 4: {
+ normal = R2{0, -1};
+ break;
+ }
+ default: {
+ FAIL("unexpected tag number");
+ }
+ }
+
+ REQUIRE(l2Norm(edge_boundary.outgoingNormal() - normal) == Catch::Approx(0).margin(1E-13));
+ }
+
+ {
+ NumberedBoundaryDescriptor numbered_boundary_descriptor(3);
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, numbered_boundary_descriptor),
+ "error: invalid boundary \"YMAX\": boundary is not flat!");
+ }
+ }
+
+ {
+ const std::set<std::string> name_set = {"XMIN", "XMAX", "YMIN"};
+
+ for (auto name : name_set) {
+ NamedBoundaryDescriptor named_boundary_descriptor(name);
+ const auto& edge_boundary = getMeshFlatEdgeBoundary(mesh, named_boundary_descriptor);
+
+ auto edge_list = get_edge_list_from_name(name, connectivity);
+
+ R2 normal = zero;
+
+ REQUIRE(is_same(edge_boundary.edgeList(), edge_list));
+
+ if (name == "XMIN") {
+ normal = R2{-1, 0};
+ } else if (name == "XMAX") {
+ normal = R2{1, 0};
+ } else if (name == "YMIN") {
+ normal = R2{0, -1};
+ } else {
+ FAIL("unexpected name: " + name);
+ }
+
+ REQUIRE(l2Norm(edge_boundary.outgoingNormal() - normal) == Catch::Approx(0).margin(1E-13));
+ }
+
+ {
+ NamedBoundaryDescriptor named_boundary_descriptor("YMAX");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, named_boundary_descriptor),
+ "error: invalid boundary \"YMAX\": boundary is not flat!");
+ }
+ }
+ }
+ }
+
+ SECTION("3D")
+ {
+ static constexpr size_t Dimension = 3;
+
+ using ConnectivityType = Connectivity<Dimension>;
+ using MeshType = Mesh<ConnectivityType>;
+
+ using R3 = TinyVector<3>;
+
+ SECTION("cartesian 3d")
+ {
+ std::shared_ptr p_mesh = MeshDataBaseForTests::get().cartesian3DMesh();
+
+ const ConnectivityType& connectivity = p_mesh->connectivity();
+
+ auto xr = p_mesh->xr();
+
+ auto curve = [](const R3& X) -> R3 {
+ return R3{X[0], (1 + X[0] * X[0]) * (X[1] + 1), (1 - 0.2 * X[0] * X[0]) * X[2]};
+ };
+
+ NodeValue<R3> curved_xr{connectivity};
+
+ parallel_for(
+ connectivity.numberOfNodes(), PUGS_LAMBDA(const NodeId node_id) { curved_xr[node_id] = curve(xr[node_id]); });
+
+ MeshType mesh{p_mesh->shared_connectivity(), curved_xr};
+
+ {
+ const std::set<size_t> tag_set = {0, 1, 2, 4};
+
+ for (auto tag : tag_set) {
+ auto edge_list = get_edge_list_from_tag(tag, connectivity);
+
+ NumberedBoundaryDescriptor numbered_boundary_descriptor(tag);
+ const auto& edge_boundary = getMeshFlatEdgeBoundary(mesh, numbered_boundary_descriptor);
+
+ REQUIRE(is_same(edge_boundary.edgeList(), edge_list));
+
+ R3 normal = zero;
+
+ switch (tag) {
+ case 0: {
+ normal = R3{-1, 0, 0};
+ break;
+ }
+ case 1: {
+ normal = R3{1, 0, 0};
+ break;
+ }
+ case 2: {
+ normal = R3{0, -1, 0};
+ break;
+ }
+ case 4: {
+ normal = R3{0, 0, -1};
+ break;
+ }
+ default: {
+ FAIL("unexpected tag number");
+ }
+ }
+
+ REQUIRE(l2Norm(edge_boundary.outgoingNormal() - normal) == Catch::Approx(0).margin(1E-13));
+ }
+
+ {
+ NumberedBoundaryDescriptor numbered_boundary_descriptor(3);
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, numbered_boundary_descriptor),
+ "error: invalid boundary \"YMAX\": boundary is not flat!");
+ }
+
+ {
+ NumberedBoundaryDescriptor numbered_boundary_descriptor(5);
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, numbered_boundary_descriptor),
+ "error: invalid boundary \"ZMAX\": boundary is not flat!");
+ }
+ }
+
+ {
+ const std::set<std::string> name_set = {"XMIN", "XMAX", "YMIN", "ZMIN"};
+
+ for (auto name : name_set) {
+ NamedBoundaryDescriptor named_boundary_descriptor(name);
+ const auto& edge_boundary = getMeshFlatEdgeBoundary(mesh, named_boundary_descriptor);
+
+ auto edge_list = get_edge_list_from_name(name, connectivity);
+
+ REQUIRE(is_same(edge_boundary.edgeList(), edge_list));
+
+ R3 normal = zero;
+
+ if (name == "XMIN") {
+ normal = R3{-1, 0, 0};
+ } else if (name == "XMAX") {
+ normal = R3{1, 0, 0};
+ } else if (name == "YMIN") {
+ normal = R3{0, -1, 0};
+ } else if (name == "ZMIN") {
+ normal = R3{0, 0, -1};
+ } else {
+ FAIL("unexpected name: " + name);
+ }
+
+ REQUIRE(l2Norm(edge_boundary.outgoingNormal() - normal) == Catch::Approx(0).margin(1E-13));
+ }
+
+ {
+ NamedBoundaryDescriptor named_boundary_descriptor("YMAX");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, named_boundary_descriptor),
+ "error: invalid boundary \"YMAX\": boundary is not flat!");
+ }
+
+ {
+ NamedBoundaryDescriptor named_boundary_descriptor("ZMAX");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, named_boundary_descriptor),
+ "error: invalid boundary \"ZMAX\": boundary is not flat!");
+ }
+ }
+ }
+
+ SECTION("hybrid 3d")
+ {
+ std::shared_ptr p_mesh = MeshDataBaseForTests::get().hybrid3DMesh();
+
+ const ConnectivityType& connectivity = p_mesh->connectivity();
+
+ auto xr = p_mesh->xr();
+
+ auto curve = [](const R3& X) -> R3 {
+ return R3{X[0], (1 + X[0] * X[0]) * X[1], (1 - 0.2 * X[0] * X[0]) * X[2]};
+ };
+
+ NodeValue<R3> curved_xr{connectivity};
+
+ parallel_for(
+ connectivity.numberOfNodes(), PUGS_LAMBDA(const NodeId node_id) { curved_xr[node_id] = curve(xr[node_id]); });
+
+ MeshType mesh{p_mesh->shared_connectivity(), curved_xr};
+
+ {
+ const std::set<size_t> tag_set = {22, 23, 25, 27};
+
+ for (auto tag : tag_set) {
+ NumberedBoundaryDescriptor numbered_boundary_descriptor(tag);
+ const auto& edge_boundary = getMeshFlatEdgeBoundary(mesh, numbered_boundary_descriptor);
+
+ auto edge_list = get_edge_list_from_tag(tag, connectivity);
+ REQUIRE(is_same(edge_boundary.edgeList(), edge_list));
+
+ R3 normal = zero;
+
+ switch (tag) {
+ case 22: {
+ normal = R3{-1, 0, 0};
+ break;
+ }
+ case 23: {
+ normal = R3{1, 0, 0};
+ break;
+ }
+ case 25: {
+ normal = R3{0, 0, -1};
+ break;
+ }
+ case 27: {
+ normal = R3{0, -1, 0};
+ break;
+ }
+ default: {
+ FAIL("unexpected tag number");
+ }
+ }
+
+ REQUIRE(l2Norm(edge_boundary.outgoingNormal() - normal) == Catch::Approx(0).margin(1E-13));
+ }
+
+ {
+ NumberedBoundaryDescriptor numbered_boundary_descriptor(24);
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, numbered_boundary_descriptor),
+ "error: invalid boundary \"ZMAX\": boundary is not flat!");
+ }
+
+ {
+ NumberedBoundaryDescriptor numbered_boundary_descriptor(26);
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, numbered_boundary_descriptor),
+ "error: invalid boundary \"YMAX\": boundary is not flat!");
+ }
+ }
+
+ {
+ const std::set<std::string> name_set = {"XMIN", "XMAX", "YMIN", "ZMIN"};
+
+ for (auto name : name_set) {
+ NamedBoundaryDescriptor named_boundary_descriptor(name);
+ const auto& edge_boundary = getMeshFlatEdgeBoundary(mesh, named_boundary_descriptor);
+
+ auto edge_list = get_edge_list_from_name(name, connectivity);
+
+ REQUIRE(is_same(edge_boundary.edgeList(), edge_list));
+
+ R3 normal = zero;
+
+ if (name == "XMIN") {
+ normal = R3{-1, 0, 0};
+ } else if (name == "XMAX") {
+ normal = R3{1, 0, 0};
+ } else if (name == "YMIN") {
+ normal = R3{0, -1, 0};
+ } else if (name == "ZMIN") {
+ normal = R3{0, 0, -1};
+ } else {
+ FAIL("unexpected name: " + name);
+ }
+
+ REQUIRE(l2Norm(edge_boundary.outgoingNormal() - normal) == Catch::Approx(0).margin(1E-13));
+ }
+
+ {
+ NamedBoundaryDescriptor named_boundary_descriptor("YMAX");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, named_boundary_descriptor),
+ "error: invalid boundary \"YMAX\": boundary is not flat!");
+ }
+
+ {
+ NamedBoundaryDescriptor named_boundary_descriptor("ZMAX");
+ REQUIRE_THROWS_WITH(getMeshFlatEdgeBoundary(mesh, named_boundary_descriptor),
+ "error: invalid boundary \"ZMAX\": boundary is not flat!");
+ }
+ }
+ }
+ }
+ }
+}