Skip to content
Snippets Groups Projects
Commit c779ff9f authored by Stéphane Del Pino's avatar Stéphane Del Pino
Browse files

Add tests for arbitrary number of ghost layers construction

parent 9e507a50
No related branches found
No related tags found
1 merge request!205High-order polynomial reconstruction
......@@ -10,6 +10,9 @@
class GlobalVariableManager
{
private:
// Give some special access for testing
friend class NbGhostLayersTester;
size_t m_connectivity_id = 0;
size_t m_mesh_id = 0;
......
......@@ -161,6 +161,7 @@ add_executable (unit_tests
add_executable (mpi_unit_tests
mpi_test_main.cpp
test_Connectivity.cpp
test_ConnectivityDispatcher.cpp
test_DiscreteFunctionDPk.cpp
test_DiscreteFunctionDPkVector.cpp
test_DiscreteFunctionIntegrator.cpp
......
#include <catch2/catch_test_macros.hpp>
#include <catch2/matchers/catch_matchers_all.hpp>
#include <mesh/CartesianMeshBuilder.hpp>
#include <mesh/Connectivity.hpp>
#include <mesh/GmshReader.hpp>
#include <mesh/Mesh.hpp>
#include <mesh/MeshVariant.hpp>
#include <utils/Messenger.hpp>
#include <MeshDataBaseForTests.hpp>
#include <filesystem>
// clazy:excludeall=non-pod-global-static
class NbGhostLayersTester
{
private:
const size_t m_original_number_of_ghost_layers;
public:
NbGhostLayersTester(const size_t number_of_ghost_layers)
: m_original_number_of_ghost_layers{GlobalVariableManager::instance().getNumberOfGhostLayers()}
{
GlobalVariableManager::instance().m_number_of_ghost_layers = number_of_ghost_layers;
}
~NbGhostLayersTester()
{
GlobalVariableManager::instance().m_number_of_ghost_layers = m_original_number_of_ghost_layers;
}
};
TEST_CASE("ConnectivityDispatcher", "[mesh]")
{
auto check_number_of_ghost_layers = [](const auto& connectivity, const size_t number_of_layers) {
// We assume that the specified number of layers can be built
// (there are enough non owned layer of cells in the connectivity)
const auto cell_is_owned = connectivity.cellIsOwned();
CellValue<size_t> cell_layer{connectivity};
cell_layer.fill(number_of_layers + 1);
NodeValue<size_t> node_layer{connectivity};
node_layer.fill(number_of_layers + 1);
auto node_to_cell_matrix = connectivity.nodeToCellMatrix();
auto cell_to_node_matrix = connectivity.cellToNodeMatrix();
parallel_for(
connectivity.numberOfCells(), PUGS_LAMBDA(const CellId cell_id) {
if (cell_is_owned[cell_id]) {
cell_layer[cell_id] = 0;
}
});
for (size_t i_layer = 0; i_layer < number_of_layers; ++i_layer) {
parallel_for(
connectivity.numberOfNodes(), PUGS_LAMBDA(const NodeId node_id) {
auto node_cell_list = node_to_cell_matrix[node_id];
size_t min_layer = cell_layer[node_cell_list[0]];
for (size_t i_cell = 1; i_cell < node_cell_list.size(); ++i_cell) {
min_layer = std::min(min_layer, cell_layer[node_cell_list[i_cell]]);
}
if (min_layer < number_of_layers + 1) {
node_layer[node_id] = min_layer;
}
});
parallel_for(
connectivity.numberOfCells(), PUGS_LAMBDA(const CellId cell_id) {
auto cell_node_list = cell_to_node_matrix[cell_id];
size_t min_layer = node_layer[cell_node_list[0]];
size_t max_layer = min_layer;
for (size_t i_node = 1; i_node < cell_node_list.size(); ++i_node) {
min_layer = std::min(min_layer, node_layer[cell_node_list[i_node]]);
max_layer = std::max(max_layer, node_layer[cell_node_list[i_node]]);
}
if ((min_layer != max_layer) or
((min_layer < number_of_layers + 1) and (cell_layer[cell_id] == number_of_layers + 1))) {
cell_layer[cell_id] = min_layer + 1;
}
});
}
auto is_boundary_face = connectivity.isBoundaryFace();
auto face_to_cell_matrix = connectivity.faceToCellMatrix();
bool has_required_number_of_ghost_layers = true;
for (FaceId face_id = 0; face_id < connectivity.numberOfFaces(); ++face_id) {
auto face_cell_list = face_to_cell_matrix[face_id];
if ((face_cell_list.size() == 1) and (not is_boundary_face[face_id])) {
const CellId face_cell_id = face_cell_list[0];
has_required_number_of_ghost_layers &= (cell_layer[face_cell_id] == number_of_layers);
}
}
REQUIRE(parallel::allReduceAnd(has_required_number_of_ghost_layers));
bool first_ghost_layer_is_1 = true;
for (FaceId face_id = 0; face_id < connectivity.numberOfFaces(); ++face_id) {
auto face_cell_list = face_to_cell_matrix[face_id];
if (face_cell_list.size() == 2) {
const CellId face_cell0_id = face_cell_list[0];
const CellId face_cell1_id = face_cell_list[1];
if (cell_is_owned[face_cell0_id] xor cell_is_owned[face_cell1_id]) {
for (size_t i_cell = 0; i_cell < face_cell_list.size(); ++i_cell) {
const CellId face_cell_id = face_cell_list[i_cell];
if (not cell_is_owned[face_cell_id]) {
first_ghost_layer_is_1 &= (cell_layer[face_cell_id] == 1);
}
}
}
}
}
REQUIRE(parallel::allReduceAnd(first_ghost_layer_is_1));
};
SECTION("1 layer meshes")
{
for (auto named_mesh : MeshDataBaseForTests::get().all1DMeshes()) {
SECTION(named_mesh.name())
{
const std::shared_ptr p_mesh = named_mesh.mesh()->get<const Mesh<1>>();
check_number_of_ghost_layers(p_mesh->connectivity(), 1);
}
}
for (auto named_mesh : MeshDataBaseForTests::get().all2DMeshes()) {
SECTION(named_mesh.name())
{
const std::shared_ptr p_mesh = named_mesh.mesh()->get<const Mesh<2>>();
check_number_of_ghost_layers(p_mesh->connectivity(), 1);
}
}
for (auto named_mesh : MeshDataBaseForTests::get().all3DMeshes()) {
SECTION(named_mesh.name())
{
const std::shared_ptr p_mesh = named_mesh.mesh()->get<const Mesh<3>>();
check_number_of_ghost_layers(p_mesh->connectivity(), 1);
}
}
}
for (size_t nb_ghost_layers = 2; nb_ghost_layers < 5; ++nb_ghost_layers) {
std::stringstream os;
os << nb_ghost_layers << " layer meshes";
SECTION(os.str())
{
REQUIRE(GlobalVariableManager::instance().getNumberOfGhostLayers() == 1);
NbGhostLayersTester nb_ghost_layers_tester(nb_ghost_layers);
REQUIRE(GlobalVariableManager::instance().getNumberOfGhostLayers() == nb_ghost_layers);
SECTION("Cartesian 1D mesh")
{
auto cartesian_1d_mesh =
CartesianMeshBuilder{TinyVector<1>{-1}, TinyVector<1>{3}, TinyVector<1, size_t>{23}}.mesh();
const std::shared_ptr p_mesh = cartesian_1d_mesh->get<const Mesh<1>>();
check_number_of_ghost_layers(p_mesh->connectivity(), nb_ghost_layers);
}
SECTION("Cartesian 2D mesh")
{
auto cartesian_2d_mesh =
CartesianMeshBuilder{TinyVector<2>{0, -1}, TinyVector<2>{3, 2}, TinyVector<2, size_t>{6, 7}}.mesh();
const std::shared_ptr p_mesh = cartesian_2d_mesh->get<const Mesh<2>>();
check_number_of_ghost_layers(p_mesh->connectivity(), nb_ghost_layers);
}
SECTION("Cartesian 3D mesh")
{
auto cartesian_3d_mesh =
CartesianMeshBuilder{TinyVector<3>{0, 1, 0}, TinyVector<3>{2, -1, 3}, TinyVector<3, size_t>{6, 7, 4}}.mesh();
const std::shared_ptr p_mesh = cartesian_3d_mesh->get<const Mesh<3>>();
check_number_of_ghost_layers(p_mesh->connectivity(), nb_ghost_layers);
}
SECTION("unordered 1d mesh")
{
const std::string filename = std::filesystem::path{PUGS_BINARY_DIR}.append("tests").append("unordered-1d.msh");
auto mesh_v = GmshReader{filename}.mesh();
const std::shared_ptr p_mesh = mesh_v->get<const Mesh<1>>();
check_number_of_ghost_layers(p_mesh->connectivity(), nb_ghost_layers);
}
SECTION("hybrid 2d mesh")
{
const std::string filename = std::filesystem::path{PUGS_BINARY_DIR}.append("tests").append("hybrid-2d.msh");
auto mesh_v = GmshReader{filename}.mesh();
const std::shared_ptr p_mesh = mesh_v->get<const Mesh<2>>();
check_number_of_ghost_layers(p_mesh->connectivity(), nb_ghost_layers);
}
SECTION("hybrid 3d mesh")
{
const std::string filename = std::filesystem::path{PUGS_BINARY_DIR}.append("tests").append("hybrid-3d.msh");
auto mesh_v = GmshReader{filename}.mesh();
const std::shared_ptr p_mesh = mesh_v->get<const Mesh<3>>();
check_number_of_ghost_layers(p_mesh->connectivity(), nb_ghost_layers);
}
}
REQUIRE(GlobalVariableManager::instance().getNumberOfGhostLayers() == 1);
}
}
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment