diff --git a/src/mesh/MeshNodeBoundary.hpp b/src/mesh/MeshNodeBoundary.hpp index 3d470179ffd541130593e1a308ee7a0907f2738b..45173f0203be1bde19b7736a094e3f9cde072055 100644 --- a/src/mesh/MeshNodeBoundary.hpp +++ b/src/mesh/MeshNodeBoundary.hpp @@ -183,7 +183,7 @@ _getNormal(const MeshType&) std::exit(1); } - return R(1); + return R{1}; } template <> @@ -288,6 +288,39 @@ _getNormal(const MeshType& mesh) zmax = x; } } +#warning re work this part to avoir parallelism dependance + Array<R3> xmin_array = parallel::allGather(xmin); + Array<R3> xmax_array = parallel::allGather(xmax); + Array<R3> ymin_array = parallel::allGather(ymin); + Array<R3> ymax_array = parallel::allGather(ymax); + Array<R3> zmin_array = parallel::allGather(zmin); + Array<R3> zmax_array = parallel::allGather(zmax); + + for (size_t i=0; i<xmin_array.size(); ++i) { + const R3& x = xmin_array[i]; + if (x[0] < xmin[0]) { xmin = x; } + } + for (size_t i=0; i<ymin_array.size(); ++i) { + const R3& x = ymin_array[i]; + if (x[1] < ymin[1]) { ymin = x; } + } + for (size_t i=0; i<zmin_array.size(); ++i) { + const R3& x = zmin_array[i]; + if (x[2] < zmin[2]) { zmin = x; } + } + for (size_t i=0; i<xmax_array.size(); ++i) { + const R3& x = xmax_array[i]; + if (x[0] > xmax[0]) { xmax = x; } + } + for (size_t i=0; i<ymax_array.size(); ++i) { + const R3& x = ymax_array[i]; + if (x[1] > ymax[1]) { ymax = x; } + } + for (size_t i=0; i<zmax_array.size(); ++i) { + const R3& x = zmax_array[i]; + if (x[2] > zmax[2]) { zmax = x; } + } + const R3 u = xmax-xmin; const R3 v = ymax-ymin; @@ -425,23 +458,37 @@ _getOutgoingNormal(const MeshType& mesh) const R3 normal = this->_getNormal(mesh); - const NodeValue<const R3>& 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 R3>& 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 {