diff --git a/src/mesh/MeshNodeBoundary.hpp b/src/mesh/MeshNodeBoundary.hpp index 45173f0203be1bde19b7736a094e3f9cde072055..dbb3f6fec539dbb1f5e0e8b9cffa3889ef2f5183 100644 --- a/src/mesh/MeshNodeBoundary.hpp +++ b/src/mesh/MeshNodeBoundary.hpp @@ -374,23 +374,36 @@ _getOutgoingNormal(const MeshType& mesh) const R normal = this->_getNormal(mesh); - const NodeValue<const R>& xr = mesh.xr(); - const auto& cell_to_node_matrix - = mesh.connectivity().cellToNodeMatrix(); + double max_height = 0; - const auto& node_to_cell_matrix - = mesh.connectivity().nodeToCellMatrix(); + if (m_node_list.size()>0) { + const NodeValue<const R>& xr = mesh.xr(); + const auto& cell_to_node_matrix + = mesh.connectivity().cellToNodeMatrix(); - const NodeId r0 = m_node_list[0]; - const CellId j0 = node_to_cell_matrix[r0][0]; - const auto& j0_nodes = cell_to_node_matrix[j0]; - double max_height = 0; - for (size_t r=0; r<j0_nodes.size(); ++r) { - const double height = (xr[j0_nodes[r]]-xr[r0], normal); + const auto& node_to_cell_matrix + = mesh.connectivity().nodeToCellMatrix(); + + const NodeId r0 = m_node_list[0]; + const CellId j0 = node_to_cell_matrix[r0][0]; + const auto& j0_nodes = cell_to_node_matrix[j0]; + + for (size_t r=0; r<j0_nodes.size(); ++r) { + const double height = (xr[j0_nodes[r]]-xr[r0], normal); + if (std::abs(height) > std::abs(max_height)) { + max_height = height; + } + } + } + + Array<double> max_height_array = parallel::allGather(max_height); + for (size_t i=0; i<max_height_array.size(); ++i) { + const double height = max_height_array[i]; if (std::abs(height) > std::abs(max_height)) { max_height = height; } } + if (max_height > 0) { return -normal; } else {