diff --git a/src/mesh/MeshLineNodeBoundary.cpp b/src/mesh/MeshLineNodeBoundary.cpp
index c7e3bb97ae7d786294dbd3d24bd8209b600404b0..71ff3e65a96e122a01f8a9a246837c07be7e6a14 100644
--- a/src/mesh/MeshLineNodeBoundary.cpp
+++ b/src/mesh/MeshLineNodeBoundary.cpp
@@ -15,21 +15,21 @@ MeshLineNodeBoundary<Dimension>::_checkBoundaryIsLine(const TinyVector<Dimension
const NodeValue<const Rd>& xr = mesh.xr();
- Rdxd P = Rdxd{identity} - tensorProduct(direction, direction);
+ const Rdxd P = Rdxd{identity} - tensorProduct(direction, direction);
bool is_bad = false;
- parallel_for(this->m_node_list.size(), [=, &is_bad](int node_id) {
- const Rd& x = xr[this->m_node_list[node_id]];
- const Rd delta = x - origin;
- if (dot(P * delta, direction) > 1E-13 * length) {
+ parallel_for(this->m_node_list.size(), [=, &is_bad](int i_node) {
+ const Rd& x = xr[this->m_node_list[i_node]];
+ const Rd delta = P * (x - origin);
+ if (dot(delta, delta) > 1E-13 * length) {
is_bad = true;
}
});
if (parallel::allReduceOr(is_bad)) {
std::ostringstream ost;
- ost << "invalid boundary " << rang::fgB::yellow << this->m_boundary_name << rang::style::reset
- << ": boundary is not a line!";
+ ost << "invalid boundary \"" << rang::fgB::yellow << this->m_boundary_name << rang::style::reset
+ << "\": boundary is not a line!";
throw NormalError(ost.str());
}
}
@@ -48,8 +48,8 @@ MeshLineNodeBoundary<2>::_getDirection(const Mesh<Connectivity<2>>& mesh)
if (xmin == xmax) {
std::ostringstream ost;
- ost << "invalid boundary " << rang::fgB::yellow << this->m_boundary_name << rang::style::reset
- << ": unable to compute direction";
+ ost << "invalid boundary \"" << rang::fgB::yellow << this->m_boundary_name << rang::style::reset
+ << "\": unable to compute direction";
throw NormalError(ost.str());
}
@@ -92,6 +92,12 @@ MeshLineNodeBoundary<3>::_getDirection(const Mesh<Connectivity<3>>& mesh)
}
const double length = l2Norm(direction);
+ if (length == 0) {
+ std::ostringstream ost;
+ ost << "invalid boundary \"" << rang::fgB::yellow << this->m_boundary_name << rang::style::reset
+ << "\": unable to compute direction";
+ throw NormalError(ost.str());
+ }
direction *= 1. / length;
this->_checkBoundaryIsLine(direction, xmin, length, mesh);
diff --git a/tests/CMakeLists.txt b/tests/CMakeLists.txt
index 85c99ead032dccef8ad4bb68269495ac3da88951..6db69f436ba9fcae7ae5957a4edddf14be12cf05 100644
--- a/tests/CMakeLists.txt
+++ b/tests/CMakeLists.txt
@@ -178,6 +178,7 @@ add_executable (mpi_unit_tests
test_MeshFaceBoundary.cpp
test_MeshFlatFaceBoundary.cpp
test_MeshFlatNodeBoundary.cpp
+ test_MeshLineNodeBoundary.cpp
test_MeshNodeBoundary.cpp
test_Messenger.cpp
test_OFStream.cpp
diff --git a/tests/test_MeshLineNodeBoundary.cpp b/tests/test_MeshLineNodeBoundary.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..af3bb4e66ecbbb6f107a2ee3a97a08f3917db2d0
--- /dev/null
+++ b/tests/test_MeshLineNodeBoundary.cpp
@@ -0,0 +1,1258 @@
+#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/MeshLineNodeBoundary.hpp>
+#include <mesh/NamedBoundaryDescriptor.hpp>
+#include <mesh/NumberedBoundaryDescriptor.hpp>
+
+// clazy:excludeall=non-pod-global-static
+
+TEST_CASE("MeshLineNodeBoundary", "[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_node_list_from_tag = [](const size_t tag, const auto& connectivity) -> Array<const NodeId> {
+ for (size_t i = 0; i < connectivity.template numberOfRefItemList<ItemType::node>(); ++i) {
+ const auto& ref_node_list = connectivity.template refItemList<ItemType::node>(i);
+ const RefId ref_id = ref_node_list.refId();
+ if (ref_id.tagNumber() == tag) {
+ return ref_node_list.list();
+ }
+ }
+ return {};
+ };
+
+ auto get_node_list_from_name = [](const std::string& name, const auto& connectivity) -> Array<const NodeId> {
+ for (size_t i = 0; i < connectivity.template numberOfRefItemList<ItemType::node>(); ++i) {
+ const auto& ref_node_list = connectivity.template refItemList<ItemType::node>(i);
+ const RefId ref_id = ref_node_list.refId();
+ if (ref_id.tagName() == name) {
+ return ref_node_list.list();
+ }
+ }
+ return {};
+ };
+
+ SECTION("aligned axis")
+ {
+ 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& node_boundary = getMeshLineNodeBoundary(mesh, numbered_boundary_descriptor);
+
+ auto node_list = get_node_list_from_tag(tag, connectivity);
+ REQUIRE(is_same(node_boundary.nodeList(), node_list));
+
+ R2 direction = zero;
+
+ switch (tag) {
+ case 0: {
+ direction = R2{0, 1};
+ break;
+ }
+ case 1: {
+ direction = R2{0, 1};
+ break;
+ }
+ case 2: {
+ direction = R2{1, 0};
+ break;
+ }
+ case 3: {
+ direction = R2{1, 0};
+ break;
+ }
+ default: {
+ FAIL("unexpected tag number");
+ }
+ }
+
+ REQUIRE(l2Norm(node_boundary.direction() - direction) == Catch::Approx(0).margin(1E-13));
+ }
+
+ {
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NumberedBoundaryDescriptor(10)),
+ "error: invalid boundary \"XMINYMIN\": unable to compute direction");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NumberedBoundaryDescriptor(11)),
+ "error: invalid boundary \"XMAXYMIN\": unable to compute direction");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NumberedBoundaryDescriptor(12)),
+ "error: invalid boundary \"XMAXYMAX\": unable to compute direction");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NumberedBoundaryDescriptor(13)),
+ "error: invalid boundary \"XMINYMAX\": unable to compute direction");
+ }
+ }
+
+ {
+ const std::set<std::string> name_set = {"XMIN", "XMAX", "YMIN", "YMAX"};
+
+ for (auto name : name_set) {
+ NamedBoundaryDescriptor named_boundary_descriptor(name);
+ const auto& node_boundary = getMeshLineNodeBoundary(mesh, named_boundary_descriptor);
+
+ auto node_list = get_node_list_from_name(name, connectivity);
+ REQUIRE(is_same(node_boundary.nodeList(), node_list));
+
+ R2 direction = zero;
+
+ if (name == "XMIN") {
+ direction = R2{0, 1};
+ } else if (name == "XMAX") {
+ direction = R2{0, 1};
+ } else if (name == "YMIN") {
+ direction = R2{1, 0};
+ } else if (name == "YMAX") {
+ direction = R2{1, 0};
+ } else {
+ FAIL("unexpected name: " + name);
+ }
+
+ REQUIRE(l2Norm(node_boundary.direction() - direction) == Catch::Approx(0).margin(1E-13));
+ }
+
+ {
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NamedBoundaryDescriptor("XMINYMIN")),
+ "error: invalid boundary \"XMINYMIN\": unable to compute direction");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NamedBoundaryDescriptor("XMINYMAX")),
+ "error: invalid boundary \"XMINYMAX\": unable to compute direction");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NamedBoundaryDescriptor("XMAXYMIN")),
+ "error: invalid boundary \"XMAXYMIN\": unable to compute direction");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NamedBoundaryDescriptor("XMAXYMAX")),
+ "error: invalid boundary \"XMAXYMAX\": unable to compute direction");
+ }
+ }
+ }
+
+ 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& node_boundary = getMeshLineNodeBoundary(mesh, numbered_boundary_descriptor);
+
+ auto node_list = get_node_list_from_tag(tag, connectivity);
+ REQUIRE(is_same(node_boundary.nodeList(), node_list));
+
+ R2 direction = zero;
+
+ switch (tag) {
+ case 1: {
+ direction = R2{0, 1};
+ break;
+ }
+ case 2: {
+ direction = R2{0, 1};
+ break;
+ }
+ case 3: {
+ direction = R2{1, 0};
+ break;
+ }
+ case 4: {
+ direction = R2{1, 0};
+ break;
+ }
+ default: {
+ FAIL("unexpected tag number");
+ }
+ }
+
+ REQUIRE(l2Norm(node_boundary.direction() - direction) == Catch::Approx(0).margin(1E-13));
+ }
+
+ {
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NumberedBoundaryDescriptor(8)),
+ "error: invalid boundary \"XMINYMIN\": unable to compute direction");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NumberedBoundaryDescriptor(9)),
+ "error: invalid boundary \"XMINYMAX\": unable to compute direction");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NumberedBoundaryDescriptor(10)),
+ "error: invalid boundary \"XMAXYMIN\": unable to compute direction");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NumberedBoundaryDescriptor(11)),
+ "error: invalid boundary \"XMAXYMAX\": unable to compute direction");
+ }
+ }
+
+ {
+ const std::set<std::string> name_set = {"XMIN", "XMAX", "YMIN", "YMAX"};
+
+ for (auto name : name_set) {
+ NamedBoundaryDescriptor named_boundary_descriptor(name);
+ const auto& node_boundary = getMeshLineNodeBoundary(mesh, named_boundary_descriptor);
+
+ auto node_list = get_node_list_from_name(name, connectivity);
+
+ R2 direction = zero;
+
+ REQUIRE(is_same(node_boundary.nodeList(), node_list));
+
+ if (name == "XMIN") {
+ direction = R2{0, 1};
+ } else if (name == "XMAX") {
+ direction = R2{0, 1};
+ } else if (name == "YMIN") {
+ direction = R2{1, 0};
+ } else if (name == "YMAX") {
+ direction = R2{1, 0};
+ } else {
+ FAIL("unexpected name: " + name);
+ }
+
+ REQUIRE(l2Norm(node_boundary.direction() - direction) == Catch::Approx(0).margin(1E-13));
+ }
+
+ {
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NamedBoundaryDescriptor("XMINYMIN")),
+ "error: invalid boundary \"XMINYMIN\": unable to compute direction");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NamedBoundaryDescriptor("XMINYMAX")),
+ "error: invalid boundary \"XMINYMAX\": unable to compute direction");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NamedBoundaryDescriptor("XMAXYMIN")),
+ "error: invalid boundary \"XMAXYMIN\": unable to compute direction");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NamedBoundaryDescriptor("XMAXYMAX")),
+ "error: invalid boundary \"XMAXYMAX\": unable to compute direction");
+ }
+ }
+ }
+ }
+
+ 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 = {20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31};
+
+ for (auto tag : tag_set) {
+ NumberedBoundaryDescriptor numbered_boundary_descriptor(tag);
+ const auto& node_boundary = getMeshLineNodeBoundary(mesh, numbered_boundary_descriptor);
+
+ auto node_list = get_node_list_from_tag(tag, connectivity);
+ REQUIRE(is_same(node_boundary.nodeList(), node_list));
+
+ R3 direction = zero;
+
+ switch (tag) {
+ case 20:
+ case 21:
+ case 22:
+ case 23: {
+ direction = R3{0, 0, -1};
+ break;
+ }
+ case 24:
+ case 25:
+ case 26:
+ case 27: {
+ direction = R3{0, -1, 0};
+ break;
+ }
+ case 28:
+ case 29:
+ case 30:
+ case 31: {
+ direction = R3{1, 0, 0};
+ break;
+ }
+ default: {
+ FAIL("unexpected tag number");
+ }
+ }
+
+ REQUIRE(l2Norm(node_boundary.direction() - direction) == Catch::Approx(0).margin(1E-13));
+ }
+
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NumberedBoundaryDescriptor(0)),
+ "error: invalid boundary \"XMIN\": boundary is not a line!");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NumberedBoundaryDescriptor(1)),
+ "error: invalid boundary \"XMAX\": boundary is not a line!");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NumberedBoundaryDescriptor(2)),
+ "error: invalid boundary \"YMIN\": boundary is not a line!");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NumberedBoundaryDescriptor(3)),
+ "error: invalid boundary \"YMAX\": boundary is not a line!");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NumberedBoundaryDescriptor(4)),
+ "error: invalid boundary \"ZMIN\": boundary is not a line!");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NumberedBoundaryDescriptor(5)),
+ "error: invalid boundary \"ZMAX\": boundary is not a line!");
+
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NumberedBoundaryDescriptor(10)),
+ "error: invalid boundary \"XMINYMINZMIN\": unable to compute direction");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NumberedBoundaryDescriptor(11)),
+ "error: invalid boundary \"XMAXYMINZMIN\": unable to compute direction");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NumberedBoundaryDescriptor(12)),
+ "error: invalid boundary \"XMAXYMAXZMIN\": unable to compute direction");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NumberedBoundaryDescriptor(13)),
+ "error: invalid boundary \"XMINYMAXZMIN\": unable to compute direction");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NumberedBoundaryDescriptor(14)),
+ "error: invalid boundary \"XMINYMINZMAX\": unable to compute direction");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NumberedBoundaryDescriptor(15)),
+ "error: invalid boundary \"XMAXYMINZMAX\": unable to compute direction");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NumberedBoundaryDescriptor(16)),
+ "error: invalid boundary \"XMAXYMAXZMAX\": unable to compute direction");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NumberedBoundaryDescriptor(17)),
+ "error: invalid boundary \"XMINYMAXZMAX\": unable to compute direction");
+ }
+
+ {
+ const std::set<std::string> name_set = {"XMINYMIN", "XMINYMAX", "XMAXYMIN", "XMAXYMAX",
+ "XMINZMIN", "XMINZMAX", "XMAXZMAX", "XMINZMAX",
+ "YMINZMIN", "YMINZMAX", "YMAXZMAX", "YMAXZMIN"};
+
+ for (const auto& name : name_set) {
+ NamedBoundaryDescriptor named_boundary_descriptor(name);
+ const auto& node_boundary = getMeshLineNodeBoundary(mesh, named_boundary_descriptor);
+
+ auto node_list = get_node_list_from_name(name, connectivity);
+
+ REQUIRE(is_same(node_boundary.nodeList(), node_list));
+
+ R3 direction = zero;
+
+ if ((name == "XMINYMIN") or (name == "XMINYMAX") or (name == "XMAXYMIN") or (name == "XMAXYMAX")) {
+ direction = R3{0, 0, -1};
+ } else if ((name == "XMINZMIN") or (name == "XMINZMAX") or (name == "XMAXZMIN") or (name == "XMAXZMAX")) {
+ direction = R3{0, -1, 0};
+ } else if ((name == "YMINZMIN") or (name == "YMINZMAX") or (name == "YMAXZMIN") or (name == "YMAXZMAX")) {
+ direction = R3{1, 0, 0};
+ } else {
+ FAIL("unexpected name: " + name);
+ }
+
+ REQUIRE(l2Norm(node_boundary.direction() - direction) == Catch::Approx(0).margin(1E-13));
+ }
+
+ {
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NamedBoundaryDescriptor("XMIN")),
+ "error: invalid boundary \"XMIN\": boundary is not a line!");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NamedBoundaryDescriptor("XMAX")),
+ "error: invalid boundary \"XMAX\": boundary is not a line!");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NamedBoundaryDescriptor("YMIN")),
+ "error: invalid boundary \"YMIN\": boundary is not a line!");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NamedBoundaryDescriptor("YMAX")),
+ "error: invalid boundary \"YMAX\": boundary is not a line!");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NamedBoundaryDescriptor("ZMIN")),
+ "error: invalid boundary \"ZMIN\": boundary is not a line!");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NamedBoundaryDescriptor("ZMAX")),
+ "error: invalid boundary \"ZMAX\": boundary is not a line!");
+
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NamedBoundaryDescriptor("XMINYMINZMIN")),
+ "error: invalid boundary \"XMINYMINZMIN\": unable to compute direction");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NamedBoundaryDescriptor("XMAXYMINZMIN")),
+ "error: invalid boundary \"XMAXYMINZMIN\": unable to compute direction");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NamedBoundaryDescriptor("XMAXYMAXZMIN")),
+ "error: invalid boundary \"XMAXYMAXZMIN\": unable to compute direction");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NamedBoundaryDescriptor("XMINYMAXZMIN")),
+ "error: invalid boundary \"XMINYMAXZMIN\": unable to compute direction");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NamedBoundaryDescriptor("XMINYMINZMAX")),
+ "error: invalid boundary \"XMINYMINZMAX\": unable to compute direction");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NamedBoundaryDescriptor("XMAXYMINZMAX")),
+ "error: invalid boundary \"XMAXYMINZMAX\": unable to compute direction");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NamedBoundaryDescriptor("XMAXYMAXZMAX")),
+ "error: invalid boundary \"XMAXYMAXZMAX\": unable to compute direction");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NamedBoundaryDescriptor("XMINYMAXZMAX")),
+ "error: invalid boundary \"XMINYMAXZMAX\": unable to compute direction");
+ }
+ }
+ }
+
+ 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 = {28, 29, 30, 31, 32, 33, 34, 35, 36, 37, 38, 39};
+
+ for (auto tag : tag_set) {
+ NumberedBoundaryDescriptor numbered_boundary_descriptor(tag);
+ const auto& node_boundary = getMeshLineNodeBoundary(mesh, numbered_boundary_descriptor);
+
+ auto node_list = get_node_list_from_tag(tag, connectivity);
+ REQUIRE(is_same(node_boundary.nodeList(), node_list));
+
+ R3 direction = zero;
+
+ switch (tag) {
+ case 30:
+ case 31:
+ case 34:
+ case 35: {
+ direction = R3{0, 0, -1};
+ break;
+ }
+ case 28:
+ case 29:
+ case 32:
+ case 33: {
+ direction = R3{0, -1, 0};
+ break;
+ }
+ case 36:
+ case 37:
+ case 38:
+ case 39: {
+ direction = R3{1, 0, 0};
+ break;
+ }
+ default: {
+ FAIL("unexpected tag number");
+ }
+ }
+
+ REQUIRE(l2Norm(node_boundary.direction() - direction) == Catch::Approx(0).margin(1E-13));
+ }
+
+ {
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NumberedBoundaryDescriptor(22)),
+ "error: invalid boundary \"XMIN\": boundary is not a line!");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NumberedBoundaryDescriptor(23)),
+ "error: invalid boundary \"XMAX\": boundary is not a line!");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NumberedBoundaryDescriptor(24)),
+ "error: invalid boundary \"ZMAX\": boundary is not a line!");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NumberedBoundaryDescriptor(25)),
+ "error: invalid boundary \"ZMIN\": boundary is not a line!");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NumberedBoundaryDescriptor(26)),
+ "error: invalid boundary \"YMAX\": boundary is not a line!");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NumberedBoundaryDescriptor(27)),
+ "error: invalid boundary \"YMIN\": boundary is not a line!");
+
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NumberedBoundaryDescriptor(40)),
+ "error: invalid boundary \"XMINYMINZMIN\": unable to compute direction");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NumberedBoundaryDescriptor(41)),
+ "error: invalid boundary \"XMAXYMINZMIN\": unable to compute direction");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NumberedBoundaryDescriptor(42)),
+ "error: invalid boundary \"XMINYMAXZMIN\": unable to compute direction");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NumberedBoundaryDescriptor(43)),
+ "error: invalid boundary \"XMINYMAXZMAX\": unable to compute direction");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NumberedBoundaryDescriptor(44)),
+ "error: invalid boundary \"XMINYMINZMAX\": unable to compute direction");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NumberedBoundaryDescriptor(45)),
+ "error: invalid boundary \"XMAXYMINZMAX\": unable to compute direction");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NumberedBoundaryDescriptor(47)),
+ "error: invalid boundary \"XMAXYMAXZMAX\": unable to compute direction");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NumberedBoundaryDescriptor(51)),
+ "error: invalid boundary \"XMAXYMAXZMIN\": unable to compute direction");
+ }
+ }
+
+ {
+ const std::set<std::string> name_set = {"XMINYMIN", "XMINYMAX", "XMAXYMIN", "XMAXYMAX",
+ "XMINZMIN", "XMINZMAX", "XMAXZMAX", "XMINZMAX",
+ "YMINZMIN", "YMINZMAX", "YMAXZMAX", "YMAXZMIN"};
+
+ for (auto name : name_set) {
+ NamedBoundaryDescriptor named_boundary_descriptor(name);
+ const auto& node_boundary = getMeshLineNodeBoundary(mesh, named_boundary_descriptor);
+
+ auto node_list = get_node_list_from_name(name, connectivity);
+
+ REQUIRE(is_same(node_boundary.nodeList(), node_list));
+
+ R3 direction = zero;
+
+ if ((name == "XMINYMIN") or (name == "XMINYMAX") or (name == "XMAXYMIN") or (name == "XMAXYMAX")) {
+ direction = R3{0, 0, -1};
+ } else if ((name == "XMINZMIN") or (name == "XMINZMAX") or (name == "XMAXZMIN") or (name == "XMAXZMAX")) {
+ direction = R3{0, -1, 0};
+ } else if ((name == "YMINZMIN") or (name == "YMINZMAX") or (name == "YMAXZMIN") or (name == "YMAXZMAX")) {
+ direction = R3{1, 0, 0};
+ } else {
+ FAIL("unexpected name: " + name);
+ }
+
+ REQUIRE(l2Norm(node_boundary.direction() - direction) == Catch::Approx(0).margin(1E-13));
+ }
+
+ {
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NamedBoundaryDescriptor("XMIN")),
+ "error: invalid boundary \"XMIN\": boundary is not a line!");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NamedBoundaryDescriptor("XMAX")),
+ "error: invalid boundary \"XMAX\": boundary is not a line!");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NamedBoundaryDescriptor("YMIN")),
+ "error: invalid boundary \"YMIN\": boundary is not a line!");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NamedBoundaryDescriptor("YMAX")),
+ "error: invalid boundary \"YMAX\": boundary is not a line!");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NamedBoundaryDescriptor("ZMIN")),
+ "error: invalid boundary \"ZMIN\": boundary is not a line!");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NamedBoundaryDescriptor("ZMAX")),
+ "error: invalid boundary \"ZMAX\": boundary is not a line!");
+
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NamedBoundaryDescriptor("XMINYMINZMIN")),
+ "error: invalid boundary \"XMINYMINZMIN\": unable to compute direction");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NamedBoundaryDescriptor("XMAXYMINZMIN")),
+ "error: invalid boundary \"XMAXYMINZMIN\": unable to compute direction");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NamedBoundaryDescriptor("XMAXYMAXZMIN")),
+ "error: invalid boundary \"XMAXYMAXZMIN\": unable to compute direction");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NamedBoundaryDescriptor("XMINYMAXZMIN")),
+ "error: invalid boundary \"XMINYMAXZMIN\": unable to compute direction");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NamedBoundaryDescriptor("XMINYMINZMAX")),
+ "error: invalid boundary \"XMINYMINZMAX\": unable to compute direction");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NamedBoundaryDescriptor("XMAXYMINZMAX")),
+ "error: invalid boundary \"XMAXYMINZMAX\": unable to compute direction");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NamedBoundaryDescriptor("XMAXYMAXZMAX")),
+ "error: invalid boundary \"XMAXYMAXZMAX\": unable to compute direction");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NamedBoundaryDescriptor("XMINYMAXZMAX")),
+ "error: invalid boundary \"XMINYMAXZMAX\": unable to compute direction");
+ }
+ }
+ }
+ }
+ }
+
+ 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& node_boundary = getMeshLineNodeBoundary(mesh, numbered_boundary_descriptor);
+
+ auto node_list = get_node_list_from_tag(tag, connectivity);
+ REQUIRE(is_same(node_boundary.nodeList(), node_list));
+
+ R2 direction = zero;
+
+ switch (tag) {
+ case 0: {
+ direction = R * R2{0, -1};
+ break;
+ }
+ case 1: {
+ direction = R * R2{0, -1};
+ break;
+ }
+ case 2: {
+ direction = R * R2{1, 0};
+ break;
+ }
+ case 3: {
+ direction = R * R2{1, 0};
+ break;
+ }
+ default: {
+ FAIL("unexpected tag number");
+ }
+ }
+
+ REQUIRE(l2Norm(node_boundary.direction() - direction) == 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& node_boundary = getMeshLineNodeBoundary(mesh, named_boundary_descriptor);
+
+ auto node_list = get_node_list_from_name(name, connectivity);
+ REQUIRE(is_same(node_boundary.nodeList(), node_list));
+
+ R2 direction = zero;
+
+ if (name == "XMIN") {
+ direction = R * R2{0, -1};
+ } else if (name == "XMAX") {
+ direction = R * R2{0, -1};
+ } else if (name == "YMIN") {
+ direction = R * R2{1, 0};
+ } else if (name == "YMAX") {
+ direction = R * R2{1, 0};
+ } else {
+ FAIL("unexpected name: " + name);
+ }
+
+ REQUIRE(l2Norm(node_boundary.direction() - direction) == 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& node_boundary = getMeshLineNodeBoundary(mesh, numbered_boundary_descriptor);
+
+ auto node_list = get_node_list_from_tag(tag, connectivity);
+ REQUIRE(is_same(node_boundary.nodeList(), node_list));
+
+ R2 direction = zero;
+
+ switch (tag) {
+ case 1: {
+ direction = R * R2{0, -1};
+ break;
+ }
+ case 2: {
+ direction = R * R2{0, -1};
+ break;
+ }
+ case 3: {
+ direction = R * R2{1, 0};
+ break;
+ }
+ case 4: {
+ direction = R * R2{1, 0};
+ break;
+ }
+ default: {
+ FAIL("unexpected tag number");
+ }
+ }
+
+ REQUIRE(l2Norm(node_boundary.direction() - direction) == 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& node_boundary = getMeshLineNodeBoundary(mesh, named_boundary_descriptor);
+
+ auto node_list = get_node_list_from_name(name, connectivity);
+
+ R2 direction = zero;
+
+ REQUIRE(is_same(node_boundary.nodeList(), node_list));
+
+ if (name == "XMIN") {
+ direction = R * R2{0, -1};
+ } else if (name == "XMAX") {
+ direction = R * R2{0, -1};
+ } else if (name == "YMIN") {
+ direction = R * R2{1, 0};
+ } else if (name == "YMAX") {
+ direction = R * R2{1, 0};
+ } else {
+ FAIL("unexpected name: " + name);
+ }
+
+ REQUIRE(l2Norm(node_boundary.direction() - direction) == 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 = {20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31};
+
+ for (auto tag : tag_set) {
+ NumberedBoundaryDescriptor numbered_boundary_descriptor(tag);
+ const auto& node_boundary = getMeshLineNodeBoundary(mesh, numbered_boundary_descriptor);
+
+ auto node_list = get_node_list_from_tag(tag, connectivity);
+ REQUIRE(is_same(node_boundary.nodeList(), node_list));
+
+ R3 direction = zero;
+
+ switch (tag) {
+ case 20:
+ case 21:
+ case 22:
+ case 23: {
+ direction = R * R3{0, 0, -1};
+ break;
+ }
+ case 24:
+ case 25:
+ case 26:
+ case 27: {
+ direction = R * R3{0, 1, 0};
+ break;
+ }
+ case 28:
+ case 29:
+ case 30:
+ case 31: {
+ direction = R * R3{-1, 0, 0};
+ break;
+ }
+ default: {
+ FAIL("unexpected tag number");
+ }
+ }
+
+ REQUIRE(l2Norm(node_boundary.direction() - direction) == Catch::Approx(0).margin(1E-13));
+ }
+ }
+
+ {
+ const std::set<std::string> name_set = {"XMINYMIN", "XMINYMAX", "XMAXYMIN", "XMAXYMAX",
+ "XMINZMIN", "XMINZMAX", "XMAXZMAX", "XMINZMAX",
+ "YMINZMIN", "YMINZMAX", "YMAXZMAX", "YMAXZMIN"};
+
+ for (auto name : name_set) {
+ NamedBoundaryDescriptor named_boundary_descriptor(name);
+ const auto& node_boundary = getMeshLineNodeBoundary(mesh, named_boundary_descriptor);
+
+ auto node_list = get_node_list_from_name(name, connectivity);
+
+ REQUIRE(is_same(node_boundary.nodeList(), node_list));
+
+ R3 direction = zero;
+
+ if ((name == "XMINYMIN") or (name == "XMINYMAX") or (name == "XMAXYMIN") or (name == "XMAXYMAX")) {
+ direction = R * R3{0, 0, -1};
+ } else if ((name == "XMINZMIN") or (name == "XMINZMAX") or (name == "XMAXZMIN") or (name == "XMAXZMAX")) {
+ direction = R * R3{0, 1, 0};
+ } else if ((name == "YMINZMIN") or (name == "YMINZMAX") or (name == "YMAXZMIN") or (name == "YMAXZMAX")) {
+ direction = R * R3{-1, 0, 0};
+ } else {
+ FAIL("unexpected name: " + name);
+ }
+
+ REQUIRE(l2Norm(node_boundary.direction() - direction) == 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 = {28, 29, 30, 31, 32, 33, 34, 35, 36, 37, 38, 39};
+
+ for (auto tag : tag_set) {
+ NumberedBoundaryDescriptor numbered_boundary_descriptor(tag);
+ const auto& node_boundary = getMeshLineNodeBoundary(mesh, numbered_boundary_descriptor);
+
+ auto node_list = get_node_list_from_tag(tag, connectivity);
+ REQUIRE(is_same(node_boundary.nodeList(), node_list));
+
+ R3 direction = zero;
+
+ switch (tag) {
+ case 30:
+ case 31:
+ case 34:
+ case 35: {
+ direction = R * R3{0, 0, -1};
+ break;
+ }
+ case 28:
+ case 29:
+ case 32:
+ case 33: {
+ direction = R * R3{0, 1, 0};
+ break;
+ }
+ case 36:
+ case 37:
+ case 38:
+ case 39: {
+ direction = R * R3{-1, 0, 0};
+ break;
+ }
+ default: {
+ FAIL("unexpected tag number");
+ }
+ }
+
+ REQUIRE(l2Norm(node_boundary.direction() - direction) == Catch::Approx(0).margin(1E-13));
+ }
+ }
+
+ {
+ const std::set<std::string> name_set = {"XMINYMIN", "XMINYMAX", "XMAXYMIN", "XMAXYMAX",
+ "XMINZMIN", "XMINZMAX", "XMAXZMAX", "XMINZMAX",
+ "YMINZMIN", "YMINZMAX", "YMAXZMAX", "YMAXZMIN"};
+
+ for (auto name : name_set) {
+ NamedBoundaryDescriptor named_boundary_descriptor(name);
+ const auto& node_boundary = getMeshLineNodeBoundary(mesh, named_boundary_descriptor);
+
+ auto node_list = get_node_list_from_name(name, connectivity);
+
+ REQUIRE(is_same(node_boundary.nodeList(), node_list));
+
+ R3 direction = zero;
+
+ if ((name == "XMINYMIN") or (name == "XMINYMAX") or (name == "XMAXYMIN") or (name == "XMAXYMAX")) {
+ direction = R * R3{0, 0, -1};
+ } else if ((name == "XMINZMIN") or (name == "XMINZMAX") or (name == "XMAXZMIN") or (name == "XMAXZMAX")) {
+ direction = R * R3{0, 1, 0};
+ } else if ((name == "YMINZMIN") or (name == "YMINZMAX") or (name == "YMAXZMIN") or (name == "YMAXZMAX")) {
+ direction = R * R3{-1, 0, 0};
+ } else {
+ FAIL("unexpected name: " + name);
+ }
+
+ REQUIRE(l2Norm(node_boundary.direction() - direction) == 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& node_boundary = getMeshLineNodeBoundary(mesh, numbered_boundary_descriptor);
+
+ auto node_list = get_node_list_from_tag(tag, connectivity);
+ REQUIRE(is_same(node_boundary.nodeList(), node_list));
+
+ R2 direction = zero;
+
+ switch (tag) {
+ case 1: {
+ direction = R2{0, 1};
+ break;
+ }
+ case 2: {
+ direction = R2{0, 1};
+ break;
+ }
+ case 4: {
+ direction = R2{1, 0};
+ break;
+ }
+ default: {
+ FAIL("unexpected tag number");
+ }
+ }
+
+ REQUIRE(l2Norm(node_boundary.direction() - direction) == Catch::Approx(0).margin(1E-13));
+ }
+
+ {
+ NumberedBoundaryDescriptor numbered_boundary_descriptor(3);
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, numbered_boundary_descriptor),
+ "error: invalid boundary \"YMAX\": boundary is not a line!");
+ }
+ }
+
+ {
+ const std::set<std::string> name_set = {"XMIN", "XMAX", "YMIN"};
+
+ for (auto name : name_set) {
+ NamedBoundaryDescriptor named_boundary_descriptor(name);
+ const auto& node_boundary = getMeshLineNodeBoundary(mesh, named_boundary_descriptor);
+
+ auto node_list = get_node_list_from_name(name, connectivity);
+
+ R2 direction = zero;
+
+ REQUIRE(is_same(node_boundary.nodeList(), node_list));
+
+ if (name == "XMIN") {
+ direction = R2{0, 1};
+ } else if (name == "XMAX") {
+ direction = R2{0, 1};
+ } else if (name == "YMIN") {
+ direction = R2{1, 0};
+ } else {
+ FAIL("unexpected name: " + name);
+ }
+
+ REQUIRE(l2Norm(node_boundary.direction() - direction) == Catch::Approx(0).margin(1E-13));
+ }
+
+ {
+ NamedBoundaryDescriptor named_boundary_descriptor("YMAX");
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, named_boundary_descriptor),
+ "error: invalid boundary \"YMAX\": boundary is not a line!");
+ }
+ }
+ }
+ }
+
+ 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 = {20, 21, 22, 23, 24, 25, 26, 27, 28};
+
+ for (auto tag : tag_set) {
+ auto node_list = get_node_list_from_tag(tag, connectivity);
+
+ NumberedBoundaryDescriptor numbered_boundary_descriptor(tag);
+ const auto& node_boundary = getMeshLineNodeBoundary(mesh, numbered_boundary_descriptor);
+
+ REQUIRE(is_same(node_boundary.nodeList(), node_list));
+
+ R3 direction = zero;
+
+ switch (tag) {
+ case 20:
+ case 21:
+ case 22:
+ case 23: {
+ direction = R3{0, 0, -1};
+ break;
+ }
+ case 24:
+ case 25:
+ case 26:
+ case 27: {
+ direction = R3{0, -1, 0};
+ break;
+ }
+ case 28: {
+ direction = R3{1, 0, 0};
+ break;
+ }
+ default: {
+ FAIL("unexpected tag number");
+ }
+ }
+
+ REQUIRE(l2Norm(node_boundary.direction() - direction) == Catch::Approx(0).margin(1E-13));
+ }
+
+ {
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NumberedBoundaryDescriptor(29)),
+ "error: invalid boundary \"YMINZMAX\": boundary is not a line!");
+
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NumberedBoundaryDescriptor(30)),
+ "error: invalid boundary \"YMAXZMAX\": boundary is not a line!");
+
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NumberedBoundaryDescriptor(31)),
+ "error: invalid boundary \"YMAXZMIN\": boundary is not a line!");
+ }
+ }
+
+ {
+ const std::set<std::string> name_set = {"XMINYMIN", "XMINYMAX", "XMAXYMIN", "XMAXYMAX", "XMINZMIN",
+ "XMINZMAX", "XMAXZMAX", "XMINZMAX", "YMINZMIN"};
+
+ for (auto name : name_set) {
+ NamedBoundaryDescriptor named_boundary_descriptor(name);
+ const auto& node_boundary = getMeshLineNodeBoundary(mesh, named_boundary_descriptor);
+
+ auto node_list = get_node_list_from_name(name, connectivity);
+
+ REQUIRE(is_same(node_boundary.nodeList(), node_list));
+
+ R3 direction = zero;
+
+ if ((name == "XMINYMIN") or (name == "XMINYMAX") or (name == "XMAXYMIN") or (name == "XMAXYMAX")) {
+ direction = R3{0, 0, -1};
+ } else if ((name == "XMINZMIN") or (name == "XMINZMAX") or (name == "XMAXZMIN") or (name == "XMAXZMAX")) {
+ direction = R3{0, -1, 0};
+ } else if ((name == "YMINZMIN")) {
+ direction = R3{1, 0, 0};
+ } else {
+ FAIL("unexpected name: " + name);
+ }
+
+ REQUIRE(l2Norm(node_boundary.direction() - direction) == Catch::Approx(0).margin(1E-13));
+ }
+
+ {
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NamedBoundaryDescriptor("YMAXZMAX")),
+ "error: invalid boundary \"YMAXZMAX\": boundary is not a line!");
+
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NamedBoundaryDescriptor("YMINZMAX")),
+ "error: invalid boundary \"YMINZMAX\": boundary is not a line!");
+
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NamedBoundaryDescriptor("YMAXZMIN")),
+ "error: invalid boundary \"YMAXZMIN\": boundary is not a line!");
+ }
+ }
+ }
+
+ 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 = {28, 29, 30, 31, 32, 33, 34, 35, 36};
+
+ for (auto tag : tag_set) {
+ NumberedBoundaryDescriptor numbered_boundary_descriptor(tag);
+ const auto& node_boundary = getMeshLineNodeBoundary(mesh, numbered_boundary_descriptor);
+
+ auto node_list = get_node_list_from_tag(tag, connectivity);
+ REQUIRE(is_same(node_boundary.nodeList(), node_list));
+
+ R3 direction = zero;
+
+ switch (tag) {
+ case 30:
+ case 31:
+ case 34:
+ case 35: {
+ direction = R3{0, 0, -1};
+ break;
+ }
+ case 28:
+ case 29:
+ case 32:
+ case 33: {
+ direction = R3{0, -1, 0};
+ break;
+ }
+ case 36:
+ // case 37:
+ // case 38:
+ case 39: {
+ direction = R3{1, 0, 0};
+ break;
+ }
+ default: {
+ FAIL("unexpected tag number");
+ }
+ }
+
+ REQUIRE(l2Norm(node_boundary.direction() - direction) == Catch::Approx(0).margin(1E-13));
+ }
+
+ {
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NumberedBoundaryDescriptor(37)),
+ "error: invalid boundary \"YMINZMAX\": boundary is not a line!");
+
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NumberedBoundaryDescriptor(38)),
+ "error: invalid boundary \"YMAXZMIN\": boundary is not a line!");
+
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NumberedBoundaryDescriptor(39)),
+ "error: invalid boundary \"YMAXZMAX\": boundary is not a line!");
+ }
+ }
+
+ {
+ const std::set<std::string> name_set = {"XMINYMIN", "XMINYMAX", "XMAXYMIN", "XMAXYMAX", "XMINZMIN",
+ "XMINZMAX", "XMAXZMAX", "XMINZMAX", "YMINZMIN"};
+
+ for (auto name : name_set) {
+ NamedBoundaryDescriptor named_boundary_descriptor(name);
+ const auto& node_boundary = getMeshLineNodeBoundary(mesh, named_boundary_descriptor);
+
+ auto node_list = get_node_list_from_name(name, connectivity);
+
+ REQUIRE(is_same(node_boundary.nodeList(), node_list));
+
+ R3 direction = zero;
+
+ if ((name == "XMINYMIN") or (name == "XMINYMAX") or (name == "XMAXYMIN") or (name == "XMAXYMAX")) {
+ direction = R3{0, 0, -1};
+ } else if ((name == "XMINZMIN") or (name == "XMINZMAX") or (name == "XMAXZMIN") or (name == "XMAXZMAX")) {
+ direction = R3{0, -1, 0};
+ } else if ((name == "YMINZMIN")) {
+ direction = R3{1, 0, 0};
+ } else {
+ FAIL("unexpected name: " + name);
+ }
+
+ REQUIRE(l2Norm(node_boundary.direction() - direction) == Catch::Approx(0).margin(1E-13));
+ }
+
+ {
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NamedBoundaryDescriptor("YMAXZMAX")),
+ "error: invalid boundary \"YMAXZMAX\": boundary is not a line!");
+
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NamedBoundaryDescriptor("YMINZMAX")),
+ "error: invalid boundary \"YMINZMAX\": boundary is not a line!");
+
+ REQUIRE_THROWS_WITH(getMeshLineNodeBoundary(mesh, NamedBoundaryDescriptor("YMAXZMIN")),
+ "error: invalid boundary \"YMAXZMIN\": boundary is not a line!");
+ }
+ }
+ }
+ }
+ }
+}