From 53231621365d5d00dcaffbd904d32c2611ae80fe Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?St=C3=A9phane=20Del=20Pino?= <stephane.delpino44@gmail.com> Date: Mon, 15 May 2023 12:31:25 +0200 Subject: [PATCH] Add zone or indicator function for smoothed zone --- src/language/modules/SchemeModule.cpp | 28 +++++ src/mesh/MeshSmoother.cpp | 150 ++++++++++++++++++++++++-- src/mesh/MeshSmoother.hpp | 12 +++ 3 files changed, 181 insertions(+), 9 deletions(-) diff --git a/src/language/modules/SchemeModule.cpp b/src/language/modules/SchemeModule.cpp index cf183409d..1eec3fb02 100644 --- a/src/language/modules/SchemeModule.cpp +++ b/src/language/modules/SchemeModule.cpp @@ -264,6 +264,34 @@ SchemeModule::SchemeModule() )); + this->_addBuiltinFunction("smoothMesh", + std::function( + + [](std::shared_ptr<const IMesh> p_mesh, + const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>& + bc_descriptor_list, + const std::vector<std::shared_ptr<const IZoneDescriptor>>& smoothing_zone_list) + -> std::shared_ptr<const IMesh> { + MeshSmootherHandler handler; + return handler.getSmoothedMesh(*p_mesh, bc_descriptor_list, smoothing_zone_list); + } + + )); + + this->_addBuiltinFunction("smoothMesh", std::function( + + [](std::shared_ptr<const IMesh> p_mesh, + const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>& + bc_descriptor_list, + const std::vector<std::shared_ptr<const DiscreteFunctionVariant>>& + discrete_function_variant_list) -> std::shared_ptr<const IMesh> { + MeshSmootherHandler handler; + return handler.getSmoothedMesh(*p_mesh, bc_descriptor_list, + discrete_function_variant_list); + } + + )); + this->_addBuiltinFunction("fixed", std::function( [](std::shared_ptr<const IBoundaryDescriptor> boundary) diff --git a/src/mesh/MeshSmoother.cpp b/src/mesh/MeshSmoother.cpp index bb6bf7289..26185119c 100644 --- a/src/mesh/MeshSmoother.cpp +++ b/src/mesh/MeshSmoother.cpp @@ -6,10 +6,12 @@ #include <mesh/Connectivity.hpp> #include <mesh/ItemValueUtils.hpp> #include <mesh/Mesh.hpp> +#include <mesh/MeshCellZone.hpp> #include <mesh/MeshFlatNodeBoundary.hpp> #include <mesh/MeshLineNodeBoundary.hpp> #include <mesh/MeshNodeBoundary.hpp> #include <scheme/AxisBoundaryConditionDescriptor.hpp> +#include <scheme/DiscreteFunctionVariant.hpp> #include <scheme/FixedBoundaryConditionDescriptor.hpp> #include <scheme/SymmetryBoundaryConditionDescriptor.hpp> #include <utils/RandomEngine.hpp> @@ -183,20 +185,11 @@ class MeshSmootherHandler::MeshSmoother const NodeId cell_node_id = cell_node_list[i_node]; mean_position += given_xr[cell_node_id]; number_of_neighbours++; - // std::cout << node_id << " position " << given_xr[node_id] << " position voisin " << - // given_xr[cell_node_id] - // << " mean_position " << mean_position << "\n"; } } } mean_position = 1. / number_of_neighbours * mean_position; shift_r[node_id] = mean_position - given_xr[node_id]; - // std::cout << " mean position " << mean_position << given_xr[node_id] << "\n"; - - // double nshift = sqrt(dot(shift_r[node_id], shift_r[node_id])); - // if (nshift > max_delta_xr[node_id]) { - // shift_r[node_id] = max_delta_xr[node_id] / nshift * shift_r[node_id]; - // } }); this->_applyBC(shift_r); @@ -237,6 +230,83 @@ class MeshSmootherHandler::MeshSmoother return std::make_shared<MeshType>(m_given_mesh.shared_connectivity(), xr); } + std::shared_ptr<const IMesh> + getSmoothedMesh(const std::vector<std::shared_ptr<const IZoneDescriptor>>& zone_descriptor_list) const + { + NodeValue<const Rd> given_xr = m_given_mesh.xr(); + + auto node_to_cell_matrix = m_given_mesh.connectivity().nodeToCellMatrix(); + + NodeValue<bool> is_displaced{m_given_mesh.connectivity()}; + is_displaced.fill(false); + + for (size_t i_zone = 0; i_zone < zone_descriptor_list.size(); ++i_zone) { + MeshCellZone<Dimension> cell_zone = getMeshCellZone(m_given_mesh, *zone_descriptor_list[i_zone]); + const auto cell_list = cell_zone.cellList(); + CellValue<bool> is_zone_cell{m_given_mesh.connectivity()}; + is_zone_cell.fill(false); + parallel_for( + cell_list.size(), PUGS_LAMBDA(const size_t i_cell) { is_zone_cell[cell_list[i_cell]] = true; }); + parallel_for( + m_given_mesh.numberOfNodes(), PUGS_LAMBDA(const NodeId node_id) { + auto node_cell_list = node_to_cell_matrix[node_id]; + bool displace = true; + for (size_t i_node_cell = 0; i_node_cell < node_cell_list.size(); ++i_node_cell) { + const CellId cell_id = node_cell_list[i_node_cell]; + displace &= is_zone_cell[cell_id]; + } + if (displace) { + is_displaced[node_id] = true; + } + }); + } + synchronize(is_displaced); + NodeValue<Rd> xr = this->_getDisplacement(); + + parallel_for( + m_given_mesh.numberOfNodes(), + PUGS_LAMBDA(const NodeId node_id) { xr[node_id] = is_displaced[node_id] * xr[node_id] + given_xr[node_id]; }); + + return std::make_shared<MeshType>(m_given_mesh.shared_connectivity(), xr); + } + + std::shared_ptr<const IMesh> + getSmoothedMesh( + const std::vector<std::shared_ptr<const DiscreteFunctionVariant>>& discrete_function_variant_list) const + { + NodeValue<const Rd> given_xr = m_given_mesh.xr(); + + auto node_to_cell_matrix = m_given_mesh.connectivity().nodeToCellMatrix(); + + NodeValue<bool> is_displaced{m_given_mesh.connectivity()}; + is_displaced.fill(false); + + for (size_t i_zone = 0; i_zone < discrete_function_variant_list.size(); ++i_zone) { + auto is_zone_cell = + discrete_function_variant_list[i_zone]->get<const DiscreteFunctionP0<Dimension, const double>>(); + parallel_for( + m_given_mesh.numberOfNodes(), PUGS_LAMBDA(const NodeId node_id) { + auto node_cell_list = node_to_cell_matrix[node_id]; + bool displace = true; + for (size_t i_node_cell = 0; i_node_cell < node_cell_list.size(); ++i_node_cell) { + const CellId cell_id = node_cell_list[i_node_cell]; + displace &= (is_zone_cell[cell_id] != 0); + } + if (displace) { + is_displaced[node_id] = true; + } + }); + } + synchronize(is_displaced); + NodeValue<Rd> xr = this->_getDisplacement(); + + parallel_for( + m_given_mesh.numberOfNodes(), + PUGS_LAMBDA(const NodeId node_id) { xr[node_id] = is_displaced[node_id] * xr[node_id] + given_xr[node_id]; }); + + return std::make_shared<MeshType>(m_given_mesh.shared_connectivity(), xr); + } + MeshSmoother(const MeshSmoother&) = delete; MeshSmoother(MeshSmoother&&) = delete; @@ -396,3 +466,65 @@ MeshSmootherHandler::getSmoothedMesh( } } } + +std::shared_ptr<const IMesh> +MeshSmootherHandler::getSmoothedMesh( + const IMesh& mesh, + const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>& bc_descriptor_list, + const std::vector<std::shared_ptr<const IZoneDescriptor>>& smoothing_zone_list) const +{ + switch (mesh.dimension()) { + case 1: { + constexpr size_t Dimension = 1; + using MeshType = Mesh<Connectivity<Dimension>>; + MeshSmoother smoother(dynamic_cast<const MeshType&>(mesh), bc_descriptor_list); + return smoother.getSmoothedMesh(smoothing_zone_list); + } + case 2: { + constexpr size_t Dimension = 2; + using MeshType = Mesh<Connectivity<Dimension>>; + MeshSmoother smoother(dynamic_cast<const MeshType&>(mesh), bc_descriptor_list); + return smoother.getSmoothedMesh(smoothing_zone_list); + } + case 3: { + constexpr size_t Dimension = 3; + using MeshType = Mesh<Connectivity<Dimension>>; + MeshSmoother smoother(dynamic_cast<const MeshType&>(mesh), bc_descriptor_list); + return smoother.getSmoothedMesh(smoothing_zone_list); + } + default: { + throw UnexpectedError("invalid mesh dimension"); + } + } +} + +std::shared_ptr<const IMesh> +MeshSmootherHandler::getSmoothedMesh( + const IMesh& mesh, + const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>& bc_descriptor_list, + const std::vector<std::shared_ptr<const DiscreteFunctionVariant>>& discrete_function_variant_list) const +{ + switch (mesh.dimension()) { + case 1: { + constexpr size_t Dimension = 1; + using MeshType = Mesh<Connectivity<Dimension>>; + MeshSmoother smoother(dynamic_cast<const MeshType&>(mesh), bc_descriptor_list); + return smoother.getSmoothedMesh(discrete_function_variant_list); + } + case 2: { + constexpr size_t Dimension = 2; + using MeshType = Mesh<Connectivity<Dimension>>; + MeshSmoother smoother(dynamic_cast<const MeshType&>(mesh), bc_descriptor_list); + return smoother.getSmoothedMesh(discrete_function_variant_list); + } + case 3: { + constexpr size_t Dimension = 3; + using MeshType = Mesh<Connectivity<Dimension>>; + MeshSmoother smoother(dynamic_cast<const MeshType&>(mesh), bc_descriptor_list); + return smoother.getSmoothedMesh(discrete_function_variant_list); + } + default: { + throw UnexpectedError("invalid mesh dimension"); + } + } +} diff --git a/src/mesh/MeshSmoother.hpp b/src/mesh/MeshSmoother.hpp index 32ba5c9e6..74fde4855 100644 --- a/src/mesh/MeshSmoother.hpp +++ b/src/mesh/MeshSmoother.hpp @@ -8,6 +8,8 @@ #include <vector> class FunctionSymbolId; +class IZoneDescriptor; +class DiscreteFunctionVariant; class MeshSmootherHandler { @@ -25,6 +27,16 @@ class MeshSmootherHandler const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>& bc_descriptor_list, const FunctionSymbolId& function_symbol_id) const; + std::shared_ptr<const IMesh> getSmoothedMesh( + const IMesh& mesh, + const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>& bc_descriptor_list, + const std::vector<std::shared_ptr<const IZoneDescriptor>>& smoothing_zone_list) const; + + std::shared_ptr<const IMesh> getSmoothedMesh( + const IMesh& mesh, + const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>& bc_descriptor_list, + const std::vector<std::shared_ptr<const DiscreteFunctionVariant>>& smoothing_zone_list) const; + MeshSmootherHandler() = default; MeshSmootherHandler(MeshSmootherHandler&&) = default; ~MeshSmootherHandler() = default; -- GitLab