From b417f2c678a5cf8476fef4c12877fee9631ee409 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?St=C3=A9phane=20Del=20Pino?= <stephane.delpino44@gmail.com>
Date: Fri, 11 Oct 2024 19:38:09 +0200
Subject: [PATCH] Add a (first?) stencil descriptor and prepare other kinds of
stencils
---
src/mesh/StencilArray.hpp | 2 +
src/mesh/StencilBuilder.cpp | 14 +++---
src/mesh/StencilBuilder.hpp | 3 +-
src/mesh/StencilDescriptor.hpp | 62 +++++++++++++++++++++++++
src/mesh/StencilManager.cpp | 60 +++++++++++++++---------
src/mesh/StencilManager.hpp | 31 ++++++++++---
src/scheme/PolynomialReconstruction.cpp | 5 +-
tests/test_StencilBuilder.cpp | 48 +++++++++++++++----
8 files changed, 177 insertions(+), 48 deletions(-)
create mode 100644 src/mesh/StencilDescriptor.hpp
diff --git a/src/mesh/StencilArray.hpp b/src/mesh/StencilArray.hpp
index d1b35def0..273e38445 100644
--- a/src/mesh/StencilArray.hpp
+++ b/src/mesh/StencilArray.hpp
@@ -99,5 +99,7 @@ class StencilArray
};
using CellToCellStencilArray = StencilArray<ItemType::cell, ItemType::cell>;
+using CellToFaceStencilArray = StencilArray<ItemType::cell, ItemType::face>;
+using NodeToCellStencilArray = StencilArray<ItemType::cell, ItemType::node>;
#endif // STENCIL_ARRAY_HPP
diff --git a/src/mesh/StencilBuilder.cpp b/src/mesh/StencilBuilder.cpp
index 920dad9e2..350492419 100644
--- a/src/mesh/StencilBuilder.cpp
+++ b/src/mesh/StencilBuilder.cpp
@@ -335,21 +335,21 @@ StencilBuilder::_build(const ConnectivityType& connectivity,
CellToCellStencilArray
StencilBuilder::build(const IConnectivity& connectivity,
- size_t number_of_layers,
+ const StencilDescriptor& stencil_descriptor,
const BoundaryDescriptorList& symmetry_boundary_descriptor_list) const
{
switch (connectivity.dimension()) {
case 1: {
- return StencilBuilder::_build(dynamic_cast<const Connectivity<1>&>(connectivity), number_of_layers,
- symmetry_boundary_descriptor_list);
+ return StencilBuilder::_build(dynamic_cast<const Connectivity<1>&>(connectivity),
+ stencil_descriptor.numberOfLayers(), symmetry_boundary_descriptor_list);
}
case 2: {
- return StencilBuilder::_build(dynamic_cast<const Connectivity<2>&>(connectivity), number_of_layers,
- symmetry_boundary_descriptor_list);
+ return StencilBuilder::_build(dynamic_cast<const Connectivity<2>&>(connectivity),
+ stencil_descriptor.numberOfLayers(), symmetry_boundary_descriptor_list);
}
case 3: {
- return StencilBuilder::_build(dynamic_cast<const Connectivity<3>&>(connectivity), number_of_layers,
- symmetry_boundary_descriptor_list);
+ return StencilBuilder::_build(dynamic_cast<const Connectivity<3>&>(connectivity),
+ stencil_descriptor.numberOfLayers(), symmetry_boundary_descriptor_list);
}
default: {
throw UnexpectedError("invalid connectivity dimension");
diff --git a/src/mesh/StencilBuilder.hpp b/src/mesh/StencilBuilder.hpp
index aaf5076f1..8b0cc09d1 100644
--- a/src/mesh/StencilBuilder.hpp
+++ b/src/mesh/StencilBuilder.hpp
@@ -3,6 +3,7 @@
#include <mesh/IBoundaryDescriptor.hpp>
#include <mesh/StencilArray.hpp>
+#include <mesh/StencilDescriptor.hpp>
#include <string>
#include <vector>
@@ -28,7 +29,7 @@ class StencilBuilder
friend class StencilManager;
CellToCellStencilArray build(const IConnectivity& connectivity,
- size_t number_of_layers,
+ const StencilDescriptor& stencil_descriptor,
const BoundaryDescriptorList& symmetry_boundary_descriptor_list) const;
public:
diff --git a/src/mesh/StencilDescriptor.hpp b/src/mesh/StencilDescriptor.hpp
new file mode 100644
index 000000000..8dfffd3e5
--- /dev/null
+++ b/src/mesh/StencilDescriptor.hpp
@@ -0,0 +1,62 @@
+#ifndef STENCIL_DESCRIPTOR_HPP
+#define STENCIL_DESCRIPTOR_HPP
+
+#include <utils/PugsMacros.hpp>
+
+#include <cstddef>
+
+class StencilDescriptor
+{
+ public:
+ enum class Type : int
+ {
+ by_nodes = 1,
+ by_edges = 2,
+ by_faces = 3
+ };
+
+ private:
+ size_t m_number_of_layers;
+ Type m_type;
+
+ public:
+ PUGS_INLINE
+ const size_t&
+ numberOfLayers() const
+ {
+ return m_number_of_layers;
+ };
+
+ PUGS_INLINE
+ const Type&
+ type() const
+ {
+ return m_type;
+ };
+
+ PUGS_INLINE
+ friend bool
+ operator==(const StencilDescriptor& sd0, const StencilDescriptor& sd1)
+ {
+ return (sd0.m_number_of_layers == sd1.m_number_of_layers) and (sd0.m_type == sd1.m_type);
+ }
+
+ PUGS_INLINE
+ friend bool
+ operator!=(const StencilDescriptor& sd0, const StencilDescriptor& sd1)
+ {
+ return not(sd0 == sd1);
+ }
+
+ StencilDescriptor(const size_t number_of_layers, const Type type) : m_number_of_layers{number_of_layers}, m_type{type}
+ {}
+
+ // StencilDescriptor() : m_number_of_layers(1), m_type{Type::by_nodes} {};
+
+ StencilDescriptor(const StencilDescriptor&) = default;
+ StencilDescriptor(StencilDescriptor&&) = default;
+
+ ~StencilDescriptor() = default;
+};
+
+#endif // STENCIL_DESCRIPTOR_HPP
diff --git a/src/mesh/StencilManager.cpp b/src/mesh/StencilManager.cpp
index ed1f9ef76..6de72c08a 100644
--- a/src/mesh/StencilManager.cpp
+++ b/src/mesh/StencilManager.cpp
@@ -18,46 +18,60 @@ StencilManager::destroy()
Assert(m_instance != nullptr, "StencilManager was not created");
// LCOV_EXCL_START
- if (m_instance->m_stored_cell_to_cell_stencil_map.size() > 0) {
- std::stringstream error;
- error << ": some connectivities are still registered\n";
- for (const auto& [key, stencil_p] : m_instance->m_stored_cell_to_cell_stencil_map) {
- error << " - connectivity of id " << rang::fgB::magenta << key.connectivity_id << rang::style::reset << '\n';
+ auto check_still_registered = [](auto& stored_stencil_map) {
+ if (stored_stencil_map.size() > 0) {
+ std::stringstream error;
+ error << ": some connectivities are still registered\n";
+ for (const auto& [key, stencil_p] : stored_stencil_map) {
+ error << " - connectivity of id " << rang::fgB::magenta << key.connectivity_id << rang::style::reset << '\n';
+ }
+ throw UnexpectedError(error.str());
}
- throw UnexpectedError(error.str());
- // LCOV_EXCL_STOP
- }
+ };
+
+ check_still_registered(m_instance->m_stored_cell_to_cell_stencil_map);
+ check_still_registered(m_instance->m_stored_cell_to_face_stencil_map);
+ check_still_registered(m_instance->m_stored_node_to_cell_stencil_map);
+ // LCOV_EXCL_STOP
+
delete m_instance;
m_instance = nullptr;
}
const CellToCellStencilArray&
StencilManager::getCellToCellStencilArray(const IConnectivity& connectivity,
- size_t degree,
+ const StencilDescriptor& stencil_descriptor,
const BoundaryDescriptorList& symmetry_boundary_descriptor_list)
{
if (not m_stored_cell_to_cell_stencil_map.contains(
- Key{connectivity.id(), degree, symmetry_boundary_descriptor_list})) {
- m_stored_cell_to_cell_stencil_map[Key{connectivity.id(), degree, symmetry_boundary_descriptor_list}] =
+ Key{connectivity.id(), stencil_descriptor, symmetry_boundary_descriptor_list})) {
+ m_stored_cell_to_cell_stencil_map[Key{connectivity.id(), stencil_descriptor, symmetry_boundary_descriptor_list}] =
std::make_shared<CellToCellStencilArray>(
- StencilBuilder{}.build(connectivity, degree, symmetry_boundary_descriptor_list));
+ StencilBuilder{}.build(connectivity, stencil_descriptor, symmetry_boundary_descriptor_list));
}
- return *m_stored_cell_to_cell_stencil_map.at(Key{connectivity.id(), degree, symmetry_boundary_descriptor_list});
+ return *m_stored_cell_to_cell_stencil_map.at(
+ Key{connectivity.id(), stencil_descriptor, symmetry_boundary_descriptor_list});
}
void
StencilManager::deleteConnectivity(const size_t connectivity_id)
{
- bool has_removed = false;
- do {
- has_removed = false;
- for (const auto& [key, p_stencil] : m_stored_cell_to_cell_stencil_map) {
- if (connectivity_id == key.connectivity_id) {
- m_stored_cell_to_cell_stencil_map.erase(key);
- has_removed = true;
- break;
+ auto delete_connectivity_stencil = [&connectivity_id](auto& stored_stencil_map) {
+ bool has_removed = false;
+ do {
+ has_removed = false;
+ for (const auto& [key, p_stencil] : stored_stencil_map) {
+ if (connectivity_id == key.connectivity_id) {
+ stored_stencil_map.erase(key);
+ has_removed = true;
+ break;
+ }
}
- }
- } while (has_removed);
+ } while (has_removed);
+ };
+
+ delete_connectivity_stencil(m_stored_cell_to_cell_stencil_map);
+ delete_connectivity_stencil(m_stored_cell_to_face_stencil_map);
+ delete_connectivity_stencil(m_stored_node_to_cell_stencil_map);
}
diff --git a/src/mesh/StencilManager.hpp b/src/mesh/StencilManager.hpp
index de0007631..6deb86536 100644
--- a/src/mesh/StencilManager.hpp
+++ b/src/mesh/StencilManager.hpp
@@ -4,6 +4,7 @@
#include <mesh/IBoundaryDescriptor.hpp>
#include <mesh/IConnectivity.hpp>
#include <mesh/StencilArray.hpp>
+#include <mesh/StencilDescriptor.hpp>
#include <memory>
#include <set>
@@ -24,13 +25,13 @@ class StencilManager
struct Key
{
size_t connectivity_id;
- size_t degree;
+ StencilDescriptor stencil_descriptor;
BoundaryDescriptorList symmetry_boundary_descriptor_list;
PUGS_INLINE bool
operator==(const Key& k) const
{
- if ((connectivity_id != k.connectivity_id) or (degree != k.degree) or
+ if ((connectivity_id != k.connectivity_id) or (stencil_descriptor != k.stencil_descriptor) or
(symmetry_boundary_descriptor_list.size() != k.symmetry_boundary_descriptor_list.size())) {
return false;
}
@@ -57,12 +58,20 @@ class StencilManager
operator()(const Key& key) const
{
// We do not use the symmetry boundary set by now
- return (std::hash<decltype(Key::connectivity_id)>()(key.connectivity_id)) ^
- (std::hash<decltype(Key::degree)>()(key.degree) << 1);
+ return ((std::hash<decltype(Key::connectivity_id)>()(key.connectivity_id)) ^
+ (std::hash<std::decay_t<decltype(Key::stencil_descriptor.numberOfLayers())>>()(
+ key.stencil_descriptor.numberOfLayers())
+ << 1));
}
};
- std::unordered_map<Key, std::shared_ptr<const CellToCellStencilArray>, HashKey> m_stored_cell_to_cell_stencil_map;
+ template <ItemType source_item_type, ItemType target_item_type>
+ using StoredStencilTMap =
+ std::unordered_map<Key, std::shared_ptr<const StencilArray<source_item_type, target_item_type>>, HashKey>;
+
+ StoredStencilTMap<ItemType::cell, ItemType::cell> m_stored_cell_to_cell_stencil_map;
+ StoredStencilTMap<ItemType::cell, ItemType::face> m_stored_cell_to_face_stencil_map;
+ StoredStencilTMap<ItemType::node, ItemType::cell> m_stored_node_to_cell_stencil_map;
public:
static void create();
@@ -80,7 +89,17 @@ class StencilManager
const CellToCellStencilArray& getCellToCellStencilArray(
const IConnectivity& i_connectivity,
- size_t degree = 1,
+ const StencilDescriptor& stencil_descriptor,
+ const BoundaryDescriptorList& symmetry_boundary_descriptor_list = {});
+
+ const CellToFaceStencilArray& getCellToFaceStencilArray(
+ const IConnectivity& i_connectivity,
+ const StencilDescriptor& stencil_descriptor,
+ const BoundaryDescriptorList& symmetry_boundary_descriptor_list = {});
+
+ const NodeToCellStencilArray& getNodeToCellStencilArray(
+ const IConnectivity& i_connectivity,
+ const StencilDescriptor& stencil_descriptor,
const BoundaryDescriptorList& symmetry_boundary_descriptor_list = {});
StencilManager(const StencilManager&) = delete;
diff --git a/src/scheme/PolynomialReconstruction.cpp b/src/scheme/PolynomialReconstruction.cpp
index 82f5081a7..76d27b500 100644
--- a/src/scheme/PolynomialReconstruction.cpp
+++ b/src/scheme/PolynomialReconstruction.cpp
@@ -19,6 +19,7 @@
#include <mesh/MeshDataManager.hpp>
#include <mesh/MeshFlatFaceBoundary.hpp>
#include <mesh/NamedBoundaryDescriptor.hpp>
+#include <mesh/StencilDescriptor.hpp>
#include <mesh/StencilManager.hpp>
#include <scheme/DiscreteFunctionDPkVariant.hpp>
#include <scheme/DiscreteFunctionUtils.hpp>
@@ -770,7 +771,9 @@ PolynomialReconstruction::_build(
DiscreteFunctionDPk<MeshType::Dimension, double>::BasisViewType::dimensionFromDegree(m_descriptor.degree());
const auto& stencil_array =
- StencilManager::instance().getCellToCellStencilArray(mesh.connectivity(), m_descriptor.degree(),
+ StencilManager::instance().getCellToCellStencilArray(mesh.connectivity(),
+ StencilDescriptor{m_descriptor.degree(),
+ StencilDescriptor::Type::by_nodes},
m_descriptor.symmetryBoundaryDescriptorList());
auto xr = mesh.xr();
diff --git a/tests/test_StencilBuilder.cpp b/tests/test_StencilBuilder.cpp
index b68f0883d..b8c6387e1 100644
--- a/tests/test_StencilBuilder.cpp
+++ b/tests/test_StencilBuilder.cpp
@@ -63,7 +63,9 @@ TEST_CASE("StencilBuilder", "[mesh]")
const Connectivity<1>& connectivity = mesh.connectivity();
- auto stencil_array = StencilManager::instance().getCellToCellStencilArray(connectivity);
+ auto stencil_array =
+ StencilManager::instance().getCellToCellStencilArray(connectivity,
+ StencilDescriptor{1, StencilDescriptor::Type::by_nodes});
REQUIRE(is_valid(connectivity, stencil_array));
}
@@ -73,7 +75,9 @@ TEST_CASE("StencilBuilder", "[mesh]")
const auto& mesh = *MeshDataBaseForTests::get().unordered1DMesh()->get<Mesh<1>>();
const Connectivity<1>& connectivity = mesh.connectivity();
- auto stencil_array = StencilManager::instance().getCellToCellStencilArray(connectivity);
+ auto stencil_array =
+ StencilManager::instance().getCellToCellStencilArray(connectivity,
+ StencilDescriptor{1, StencilDescriptor::Type::by_nodes});
REQUIRE(is_valid(connectivity, stencil_array));
}
@@ -86,7 +90,10 @@ TEST_CASE("StencilBuilder", "[mesh]")
const auto& mesh = *MeshDataBaseForTests::get().cartesian2DMesh()->get<Mesh<2>>();
const Connectivity<2>& connectivity = mesh.connectivity();
- REQUIRE(is_valid(connectivity, StencilManager::instance().getCellToCellStencilArray(connectivity)));
+ REQUIRE(is_valid(connectivity,
+ StencilManager::instance()
+ .getCellToCellStencilArray(connectivity,
+ StencilDescriptor{1, StencilDescriptor::Type::by_nodes})));
}
SECTION("hybrid")
@@ -94,7 +101,10 @@ TEST_CASE("StencilBuilder", "[mesh]")
const auto& mesh = *MeshDataBaseForTests::get().hybrid2DMesh()->get<Mesh<2>>();
const Connectivity<2>& connectivity = mesh.connectivity();
- REQUIRE(is_valid(connectivity, StencilManager::instance().getCellToCellStencilArray(connectivity)));
+ REQUIRE(is_valid(connectivity,
+ StencilManager::instance()
+ .getCellToCellStencilArray(connectivity,
+ StencilDescriptor{1, StencilDescriptor::Type::by_nodes})));
}
}
@@ -105,7 +115,10 @@ TEST_CASE("StencilBuilder", "[mesh]")
const auto& mesh = *MeshDataBaseForTests::get().cartesian3DMesh()->get<Mesh<3>>();
const Connectivity<3>& connectivity = mesh.connectivity();
- REQUIRE(is_valid(connectivity, StencilManager::instance().getCellToCellStencilArray(connectivity)));
+ REQUIRE(is_valid(connectivity,
+ StencilManager::instance()
+ .getCellToCellStencilArray(connectivity,
+ StencilDescriptor{1, StencilDescriptor::Type::by_nodes})));
}
SECTION("hybrid")
@@ -113,7 +126,10 @@ TEST_CASE("StencilBuilder", "[mesh]")
const auto& mesh = *MeshDataBaseForTests::get().hybrid3DMesh()->get<Mesh<3>>();
const Connectivity<3>& connectivity = mesh.connectivity();
- REQUIRE(is_valid(connectivity, StencilManager::instance().getCellToCellStencilArray(connectivity)));
+ REQUIRE(is_valid(connectivity,
+ StencilManager::instance()
+ .getCellToCellStencilArray(connectivity,
+ StencilDescriptor{1, StencilDescriptor::Type::by_nodes})));
}
}
}
@@ -221,7 +237,10 @@ TEST_CASE("StencilBuilder", "[mesh]")
const Connectivity<2>& connectivity = mesh.connectivity();
- auto stencil_array = StencilManager::instance().getCellToCellStencilArray(connectivity, 1, list);
+ auto stencil_array =
+ StencilManager::instance().getCellToCellStencilArray(connectivity,
+ StencilDescriptor{1, StencilDescriptor::Type::by_nodes},
+ list);
REQUIRE(stencil_array.symmetryBoundaryStencilArrayList().size() == 4);
REQUIRE(check_ghost_cells_have_empty_stencils(stencil_array, connectivity));
@@ -233,7 +252,10 @@ TEST_CASE("StencilBuilder", "[mesh]")
const auto& mesh = *MeshDataBaseForTests::get().hybrid2DMesh()->get<Mesh<2>>();
const Connectivity<2>& connectivity = mesh.connectivity();
- auto stencil = StencilManager::instance().getCellToCellStencilArray(connectivity, 1, list);
+ auto stencil =
+ StencilManager::instance().getCellToCellStencilArray(connectivity,
+ StencilDescriptor{1, StencilDescriptor::Type::by_nodes},
+ list);
REQUIRE(stencil.symmetryBoundaryStencilArrayList().size() == 4);
REQUIRE(check_ghost_cells_have_empty_stencils(stencil, connectivity));
@@ -256,7 +278,10 @@ TEST_CASE("StencilBuilder", "[mesh]")
const auto& mesh = *MeshDataBaseForTests::get().cartesian3DMesh()->get<Mesh<3>>();
const Connectivity<3>& connectivity = mesh.connectivity();
- auto stencil = StencilManager::instance().getCellToCellStencilArray(connectivity, 1, list);
+ auto stencil =
+ StencilManager::instance().getCellToCellStencilArray(connectivity,
+ StencilDescriptor{1, StencilDescriptor::Type::by_nodes},
+ list);
REQUIRE(stencil.symmetryBoundaryStencilArrayList().size() == 6);
REQUIRE(check_ghost_cells_have_empty_stencils(stencil, connectivity));
@@ -268,7 +293,10 @@ TEST_CASE("StencilBuilder", "[mesh]")
const auto& mesh = *MeshDataBaseForTests::get().hybrid3DMesh()->get<Mesh<3>>();
const Connectivity<3>& connectivity = mesh.connectivity();
- auto stencil = StencilManager::instance().getCellToCellStencilArray(connectivity, 1, list);
+ auto stencil =
+ StencilManager::instance().getCellToCellStencilArray(connectivity,
+ StencilDescriptor{1, StencilDescriptor::Type::by_nodes},
+ list);
REQUIRE(stencil.symmetryBoundaryStencilArrayList().size() == 6);
REQUIRE(check_ghost_cells_have_empty_stencils(stencil, connectivity));
--
GitLab