diff --git a/src/mesh/MeshRandomizer.hpp b/src/mesh/MeshRandomizer.hpp index 36b3b169d37d27d473e9691242bb7882f993bb62..36f38e779620cf7844303756a3a54b081da61c9a 100644 --- a/src/mesh/MeshRandomizer.hpp +++ b/src/mesh/MeshRandomizer.hpp @@ -143,29 +143,34 @@ class MeshRandomizer const ConnectivityType& connectivity = m_given_mesh.connectivity(); NodeValue<const Rd> given_xr = m_given_mesh.xr(); - auto edge_to_node_matrix = connectivity.edgeToNodeMatrix(); - EdgeValue<double> edge_length{connectivity}; - parallel_for( - connectivity.numberOfEdges(), PUGS_LAMBDA(const EdgeId edge_id) { - const auto& edge_nodes = edge_to_node_matrix[edge_id]; - const NodeId node_0 = edge_nodes[0]; - const NodeId node_1 = edge_nodes[1]; - const Rd delta_x = given_xr[node_0] - given_xr[node_1]; - - edge_length[edge_id] = std::sqrt(dot(delta_x, delta_x)); - }); + auto node_to_cell_matrix = connectivity.nodeToCellMatrix(); + auto cell_to_node_matrix = connectivity.cellToNodeMatrix(); + auto node_number_in_their_cells = connectivity.nodeLocalNumbersInTheirCells(); - auto node_to_edge_matrix = connectivity.nodeToEdgeMatrix(); NodeValue<double> max_delta_xr{connectivity}; parallel_for( connectivity.numberOfNodes(), PUGS_LAMBDA(const NodeId node_id) { - const auto& node_edges = node_to_edge_matrix[node_id]; - double max_delta = edge_length[node_edges[0]]; + const Rd& x0 = given_xr[node_id]; + + const auto& node_cell_list = node_to_cell_matrix[node_id]; + double min_distance_2 = std::numeric_limits<double>::max(); - for (size_t i_edge = 1; i_edge < node_edges.size(); ++i_edge) { - const EdgeId edge_id = node_edges[i_edge]; - max_delta = std::min(max_delta, edge_length[edge_id]); + for (size_t i_cell = 0; i_cell < node_cell_list.size(); ++i_cell) { + const size_t i_cell_node = node_number_in_their_cells(node_id, i_cell); + + const CellId cell_id = node_cell_list[i_cell]; + const auto& cell_node_list = cell_to_node_matrix[cell_id]; + + for (size_t i_node = 0; i_node < cell_node_list.size(); ++i_node) { + if (i_node != i_cell_node) { + const NodeId cell_node_id = cell_node_list[i_node]; + const Rd delta = x0 - given_xr[cell_node_id]; + min_distance_2 = std::min(min_distance_2, dot(delta, delta)); + } + } } + double max_delta = std::sqrt(min_distance_2); + max_delta_xr[node_id] = max_delta; });