diff --git a/src/mesh/CMakeLists.txt b/src/mesh/CMakeLists.txt index 3569306a44f5fce39435bd24a457a25a6e660f8e..9c232ea60db552e1f032a8b4a1b45cf199c33d4a 100644 --- a/src/mesh/CMakeLists.txt +++ b/src/mesh/CMakeLists.txt @@ -8,6 +8,7 @@ add_library( Connectivity.cpp ConnectivityComputer.cpp GmshReader.cpp + MeshDispatcher.cpp SynchronizerManager.cpp) include_directories(${PASTIS_BINARY_DIR}/src/utils) diff --git a/src/mesh/GmshReader.cpp b/src/mesh/GmshReader.cpp index 23e74c0fd5769d6268781fee690cef9f677f88a7..29b6e4ebe81d6e9b14c4042232f11fdc3f9dcb22 100644 --- a/src/mesh/GmshReader.cpp +++ b/src/mesh/GmshReader.cpp @@ -14,11 +14,12 @@ #include <RefFaceList.hpp> #include <Messenger.hpp> -#include <Partitioner.hpp> #include <ArrayUtils.hpp> #include <ItemValueUtils.hpp> +#include <MeshDispatcher.hpp> + #include <unordered_map> #include <map> #include <regex> @@ -130,6 +131,26 @@ ErrorHandler(const std::string& filename, ; } +template <int Dimension> +void GmshReader::_dispatch() +{ + using ConnectivityType = Connectivity<Dimension>; + using Rd = TinyVector<Dimension>; + using MeshType = Mesh<ConnectivityType>; + + if (not m_mesh) { + ConnectivityDescriptor descriptor; + std::shared_ptr connectivity = ConnectivityType::build(descriptor); + NodeValue<Rd> xr; + m_mesh = std::make_shared<MeshType>(connectivity, xr); + } + const MeshType& mesh = static_cast<const MeshType&>(*m_mesh); + + MeshDispatcher<Dimension> dispatcher(mesh); + m_mesh = dispatcher.dispatchedMesh(); +} + + template <size_t Dimension> class ConnectivityFace; @@ -328,864 +349,6 @@ class ConnectivityFace<3> ~ConnectivityFace() = default; }; -template <int Dimension> -class MeshDispatcher -{ - public: - using MeshType = Mesh<Connectivity<Dimension>>; - - private: - const MeshType& m_mesh; - CellValue<const int> m_cell_new_owner; - FaceValue<const int> m_face_new_owner; - NodeValue<const int> m_node_new_owner; - - using CellListToSendByProc = std::vector<Array<const CellId>>; - const CellListToSendByProc m_cell_list_to_send_by_proc; - - Array<int> m_nb_cell_to_send_by_proc; - Array<int> m_nb_cell_to_recv_by_proc; - - CellValue<int> _getCellNewOwner() - { - CSRGraph mesh_graph = m_mesh.cellToCellGraph(); - Partitioner P; - - CellValue<int> cell_new_owner(m_mesh.connectivity()); - cell_new_owner = P.partition(mesh_graph); - return cell_new_owner; - } - - FaceValue<int> _getFaceNewOwner() - { - const auto& face_to_cell_matrix - = m_mesh.connectivity().faceToCellMatrix(); - const auto& cell_number = m_mesh.connectivity().cellNumber(); - -#warning could use a better policy - FaceValue<int> face_new_owner(m_mesh.connectivity()); - parallel_for(m_mesh.numberOfFaces(), PASTIS_LAMBDA(const FaceId& l) { - const auto& face_to_cell = face_to_cell_matrix[l]; - CellId Jmin = face_to_cell[0]; - - for (size_t j=1; j<face_to_cell.size(); ++j) { - const CellId J = face_to_cell[j]; - if (cell_number[J] < cell_number[Jmin]) { - Jmin=J; - } - } - face_new_owner[l] = m_cell_new_owner[Jmin]; - }); - - synchronize(face_new_owner); - return face_new_owner; - } - - NodeValue<int> _getNodeNewOwner() - { - const auto& node_to_cell_matrix - = m_mesh.connectivity().nodeToCellMatrix(); - const auto& cell_number = m_mesh.connectivity().cellNumber(); - -#warning could use a better policy - NodeValue<int> node_new_owner(m_mesh.connectivity()); - parallel_for(m_mesh.numberOfNodes(), PASTIS_LAMBDA(const NodeId& r) { - const auto& node_to_cell = node_to_cell_matrix[r]; - CellId Jmin = node_to_cell[0]; - - for (size_t j=1; j<node_to_cell.size(); ++j) { - const CellId J = node_to_cell[j]; - if (cell_number[J] < cell_number[Jmin]) { - Jmin=J; - } - } - node_new_owner[r] = m_cell_new_owner[Jmin]; - }); - - synchronize(node_new_owner); - return node_new_owner; - } - - const CellListToSendByProc - _buildCellListToSend() const - { - const auto& node_to_cell_matrix - = m_mesh.connectivity().nodeToCellMatrix(); - const auto& cell_to_node_matrix - = m_mesh.connectivity().cellToNodeMatrix(); - - std::vector<std::vector<CellId>> cell_vector_to_send_by_proc(parallel::size()); - Array<bool> send_to_rank(parallel::size()); - for (CellId j=0; j<m_mesh.numberOfCells(); ++j) { - send_to_rank.fill(false); - const auto& cell_to_node = cell_to_node_matrix[j]; - - for (size_t R=0; R<cell_to_node.size(); ++R) { - const NodeId& r = cell_to_node[R]; - const auto& node_to_cell = node_to_cell_matrix[r]; - for (size_t K=0; K<node_to_cell.size(); ++K) { - const CellId& k = node_to_cell[K]; - send_to_rank[m_cell_new_owner[k]] = true; - } - } - - for (size_t k=0; k<send_to_rank.size(); ++k) { - if (send_to_rank[k]) { - cell_vector_to_send_by_proc[k].push_back(j); - } - } - } - - std::vector<Array<const CellId>> cell_list_to_send_by_proc(parallel::size()); - for (size_t i=0; i<parallel::size(); ++i) { - cell_list_to_send_by_proc[i] = convert_to_array(cell_vector_to_send_by_proc[i]); - } - - return cell_list_to_send_by_proc; - } - - Array<int> _buildNbCellToSend() - { - Array<int> nb_cell_to_send_by_proc(parallel::size()); - for (size_t i=0; i<parallel::size(); ++i) { - nb_cell_to_send_by_proc[i] = m_cell_list_to_send_by_proc[i].size(); - } - return nb_cell_to_send_by_proc; - } - - public: - - PASTIS_INLINE - const CellValue<const int>& cellNewOwner() const - { - return m_cell_new_owner; - } - - PASTIS_INLINE - const FaceValue<const int>& faceNewOwner() const - { - return m_face_new_owner; - } - - PASTIS_INLINE - const NodeValue<const int>& nodeNewOwner() const - { - return m_node_new_owner; - } - - template<typename DataType, typename ConnectivityPtr> - std::vector<Array<std::remove_const_t<DataType>>> - exchange(ItemValue<DataType, ItemType::cell, ConnectivityPtr> cell_value) const - { - using MutableDataType = std::remove_const_t<DataType>; - std::vector<Array<DataType>> cell_value_to_send_by_proc(parallel::size()); - for (size_t i=0; i<parallel::size(); ++i) { - const Array<const CellId>& cell_list = m_cell_list_to_send_by_proc[i]; - Array<MutableDataType> cell_value_list(cell_list.size()); - parallel_for (cell_list.size(), PASTIS_LAMBDA(const CellId& cell_id) { - cell_value_list[cell_id] = cell_value[cell_list[cell_id]]; - }); - cell_value_to_send_by_proc[i] = cell_value_list; - } - - std::vector<Array<MutableDataType>> recv_cell_value_by_proc(parallel::size()); - for (size_t i=0; i<parallel::size(); ++i) { - recv_cell_value_by_proc[i] = Array<MutableDataType>(m_nb_cell_to_recv_by_proc[i]); - } - - parallel::exchange(cell_value_to_send_by_proc, recv_cell_value_by_proc); - return recv_cell_value_by_proc; - } - - [[deprecated]] - const CellListToSendByProc& cell_list_to_send_by_proc() const - { - return m_cell_list_to_send_by_proc; - } - - MeshDispatcher(const MeshType& mesh) - : m_mesh(mesh), - m_cell_new_owner(_getCellNewOwner()), - m_face_new_owner(_getFaceNewOwner()), - m_node_new_owner(_getNodeNewOwner()), - m_cell_list_to_send_by_proc(_buildCellListToSend()), - m_nb_cell_to_send_by_proc(_buildNbCellToSend()), - m_nb_cell_to_recv_by_proc(parallel::allToAll(m_nb_cell_to_send_by_proc)) - { - ; - } - - MeshDispatcher(const MeshDispatcher&) = delete; - ~MeshDispatcher() = default; -}; - - - -template <int Dimension> -void GmshReader::_dispatch() -{ - using ConnectivityType = Connectivity<Dimension>; - using Rd = TinyVector<Dimension>; - using MeshType = Mesh<ConnectivityType>; - - if (not m_mesh) { - ConnectivityDescriptor descriptor; - std::shared_ptr connectivity = ConnectivityType::build(descriptor); - NodeValue<Rd> xr; - m_mesh = std::make_shared<MeshType>(connectivity, xr); - } - const MeshType& mesh = static_cast<const MeshType&>(*m_mesh); - - MeshDispatcher<Dimension> dispatcher(mesh); - - std::vector<Array<int>> recv_cell_number_by_proc - = dispatcher.exchange(mesh.connectivity().cellNumber()); - - std::vector<Array<CellType>> recv_cell_type_by_proc - = dispatcher.exchange(mesh.connectivity().cellType()); - - std::vector<Array<int>> recv_cell_new_owner_by_proc - = dispatcher.exchange(dispatcher.cellNewOwner()); - - CellValue<int> number_of_node_per_cell(mesh.connectivity()); - - const auto& cell_to_node_matrix = mesh.connectivity().cellToNodeMatrix(); - parallel_for(mesh.numberOfCells(), PASTIS_LAMBDA(const CellId& j){ - number_of_node_per_cell[j] = cell_to_node_matrix[j].size(); - }); - - std::vector<Array<int>> recv_number_of_node_per_cell_by_proc - = dispatcher.exchange(number_of_node_per_cell); - - const auto& cell_list_to_send_by_proc = dispatcher.cell_list_to_send_by_proc(); - - const NodeValue<const int>& node_number = mesh.connectivity().nodeNumber(); - std::vector<Array<int>> cell_node_number_to_send_by_proc(parallel::size()); - for (size_t i_rank=0; i_rank < parallel::size(); ++i_rank) { - std::vector<int> node_number_by_cell_vector; - for (size_t j=0; j<cell_list_to_send_by_proc[i_rank].size(); ++j) { - const CellId& cell_id = cell_list_to_send_by_proc[i_rank][j]; - const auto& cell_node_list = cell_to_node_matrix[cell_id]; - for (size_t r=0; r<cell_node_list.size(); ++r) { - const NodeId& node_id = cell_node_list[r]; - node_number_by_cell_vector.push_back(node_number[node_id]); - } - } - cell_node_number_to_send_by_proc[i_rank] = convert_to_array(node_number_by_cell_vector); - } - - std::vector<Array<int>> recv_cell_node_number_by_proc(parallel::size()); - for (size_t i_rank=0; i_rank < parallel::size(); ++i_rank) { - recv_cell_node_number_by_proc[i_rank] - = Array<int>(sum(recv_number_of_node_per_cell_by_proc[i_rank])); - } - - parallel::exchange(cell_node_number_to_send_by_proc, recv_cell_node_number_by_proc); - - const std::unordered_map<int, int> node_number_id_map - = [&] () { - std::unordered_map<int, int> node_number_id_map; - for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) { - int cpt=0; - for (size_t i=0; i<recv_cell_node_number_by_proc[i_rank].size(); ++i) { - int node_number = recv_cell_node_number_by_proc[i_rank][i]; - auto [iterator, inserted] = node_number_id_map.insert(std::make_pair(node_number, cpt)); - if (inserted) cpt++; - } - } - return node_number_id_map; - } (); - - - ConnectivityDescriptor new_descriptor; - - new_descriptor.node_number_vector.resize(node_number_id_map.size()); - for (const auto& [number, id] : node_number_id_map) { - new_descriptor.node_number_vector[id] = number; - } - - for (size_t i_rank=0; i_rank < parallel::size(); ++i_rank) { - int l=0; - for (size_t i=0; i<recv_cell_number_by_proc[i_rank].size(); ++i) { - std::vector<unsigned int> node_vector; - for (int k=0; k<recv_number_of_node_per_cell_by_proc[i_rank][i]; ++k) { - const auto& searched_node_id = node_number_id_map.find(recv_cell_node_number_by_proc[i_rank][l++]); - Assert(searched_node_id != node_number_id_map.end()); - node_vector.push_back(searched_node_id->second); - } - new_descriptor.cell_by_node_vector.emplace_back(node_vector); - } - } - - const std::unordered_map<int, int> cell_number_id_map - = [&] () { - std::unordered_map<int, int> cell_number_id_map; - for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) { - int cpt=0; - for (size_t i=0; i<recv_cell_number_by_proc[i_rank].size(); ++i) { - int cell_number = recv_cell_number_by_proc[i_rank][i]; - auto [iterator, inserted] = cell_number_id_map.insert(std::make_pair(cell_number, cpt)); - if (inserted) cpt++; - } - } - return cell_number_id_map; - } (); - - new_descriptor.cell_number_vector.resize(cell_number_id_map.size()); - for (const auto& [number, id] : cell_number_id_map) { - new_descriptor.cell_number_vector[id] = number; - } - - new_descriptor.cell_type_vector.resize(cell_number_id_map.size()); - for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) { - for (size_t i=0; i<recv_cell_number_by_proc[i_rank].size(); ++i) { - int cell_number = recv_cell_number_by_proc[i_rank][i]; - const auto& searched_cell_id = cell_number_id_map.find(cell_number); - Assert(searched_cell_id != cell_number_id_map.end()); - new_descriptor.cell_type_vector[searched_cell_id->second] = recv_cell_type_by_proc[i_rank][i]; - } - } - - std::vector<Array<const NodeId>> send_node_id_by_proc(parallel::size()); - for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) { - Array<bool> tag(mesh.numberOfNodes()); - tag.fill(false); - std::vector<NodeId> node_id_vector; - for (size_t j=0; j<cell_list_to_send_by_proc[i_rank].size(); ++j) { - const CellId& cell_id = cell_list_to_send_by_proc[i_rank][j]; - const auto& cell_node_list = cell_to_node_matrix[cell_id]; - for (size_t r=0; r<cell_node_list.size(); ++r) { - const NodeId& node_id = cell_node_list[r]; - if (not tag[node_id]) { - node_id_vector.push_back(node_id); - tag[node_id] = true; - } - } - } - send_node_id_by_proc[i_rank] = convert_to_array(node_id_vector); - } - - std::vector<Array<const int>> send_node_number_by_proc(parallel::size()); - for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) { - Array<int> send_node_number(send_node_id_by_proc[i_rank].size()); - const Array<const NodeId> send_node_id = send_node_id_by_proc[i_rank]; - parallel_for(send_node_number.size(), PASTIS_LAMBDA(const size_t& j){ - send_node_number[j] = node_number[send_node_id[j]]; - }); - send_node_number_by_proc[i_rank] = send_node_number; - } - - Array<unsigned int> nb_node_to_send_by_proc(parallel::size()); - for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) { - nb_node_to_send_by_proc[i_rank] = send_node_id_by_proc[i_rank].size(); - } - Array<const unsigned int> nb_node_to_recv_by_proc - = parallel::allToAll(nb_node_to_send_by_proc); - - std::vector<Array<int>> recv_node_number_by_proc(parallel::size()); - for (size_t i_rank=0; i_rank < parallel::size(); ++i_rank) { - recv_node_number_by_proc[i_rank] = Array<int>(nb_node_to_recv_by_proc[i_rank]); - } - parallel::exchange(send_node_number_by_proc, recv_node_number_by_proc); - - std::vector<Array<const NodeId>> recv_node_id_correspondance_by_proc(parallel::size()); - for (size_t i_rank=0; i_rank<nb_node_to_recv_by_proc.size(); ++i_rank) { - Array<NodeId> node_id_correspondace(nb_node_to_recv_by_proc[i_rank]); - for (size_t l=0; l<nb_node_to_recv_by_proc[i_rank]; ++l) { - const int& node_number = recv_node_number_by_proc[i_rank][l]; - const auto& searched_node_id = node_number_id_map.find(node_number); - Assert(searched_node_id != node_number_id_map.end()); - node_id_correspondace[l] = searched_node_id->second; - } - recv_node_id_correspondance_by_proc[i_rank] = node_id_correspondace; - } - - new_descriptor.cell_owner_vector.resize(cell_number_id_map.size()); - for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) { - for (size_t i=0; i<recv_cell_number_by_proc[i_rank].size(); ++i) { - int cell_number = recv_cell_number_by_proc[i_rank][i]; - const auto& searched_cell_id = cell_number_id_map.find(cell_number); - Assert(searched_cell_id != cell_number_id_map.end()); - new_descriptor.cell_owner_vector[searched_cell_id->second] = recv_cell_new_owner_by_proc[i_rank][i]; - } - } - - const NodeValue<const int>& new_node_owner = dispatcher.nodeNewOwner(); - std::vector<Array<const int>> send_node_owner_by_proc(parallel::size()); - for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) { - Array<int> send_node_owner(nb_node_to_send_by_proc[i_rank]); - const Array<const NodeId> send_node_id = send_node_id_by_proc[i_rank]; - parallel_for(send_node_id.size(), PASTIS_LAMBDA(const size_t& r) { - const NodeId& node_id = send_node_id[r]; - send_node_owner[r] = new_node_owner[node_id]; - }); - send_node_owner_by_proc[i_rank] = send_node_owner; - } - - std::vector<Array<int>> recv_node_owner_by_proc(parallel::size()); - for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) { - recv_node_owner_by_proc[i_rank] = Array<int>(nb_node_to_recv_by_proc[i_rank]); - } - parallel::exchange(send_node_owner_by_proc, recv_node_owner_by_proc); - - new_descriptor.node_owner_vector.resize(new_descriptor.node_number_vector.size()); - for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) { - for (size_t r=0; r<recv_node_owner_by_proc[i_rank].size(); ++r) { - const NodeId& node_id = recv_node_id_correspondance_by_proc[i_rank][r]; - new_descriptor.node_owner_vector[node_id] = recv_node_owner_by_proc[i_rank][r]; - } - } - - if constexpr(Dimension>1) { // Faces - CellValue<int> number_of_face_per_cell(mesh.connectivity()); - const auto& cell_to_face_matrix = mesh.connectivity().cellToFaceMatrix(); - parallel_for(mesh.numberOfCells(), PASTIS_LAMBDA(const CellId& j){ - number_of_face_per_cell[j] = cell_to_face_matrix[j].size(); - }); - - std::vector<Array<int>> recv_number_of_face_per_cell_by_proc - = dispatcher.exchange(number_of_face_per_cell); - - - const FaceValue<const int>& face_number = mesh.connectivity().faceNumber(); - std::vector<Array<int>> cell_face_number_to_send_by_proc(parallel::size()); - for (size_t i_rank=0; i_rank < parallel::size(); ++i_rank) { - std::vector<int> face_number_by_cell_vector; - for (size_t j=0; j<cell_list_to_send_by_proc[i_rank].size(); ++j) { - const CellId& cell_id = cell_list_to_send_by_proc[i_rank][j]; - const auto& cell_face_list = cell_to_face_matrix[cell_id]; - for (size_t r=0; r<cell_face_list.size(); ++r) { - const FaceId& face_id = cell_face_list[r]; - face_number_by_cell_vector.push_back(face_number[face_id]); - } - } - cell_face_number_to_send_by_proc[i_rank] = convert_to_array(face_number_by_cell_vector); - } - - - std::vector<Array<int>> recv_cell_face_number_by_proc(parallel::size()); - for (size_t i_rank=0; i_rank < parallel::size(); ++i_rank) { - recv_cell_face_number_by_proc[i_rank] - = Array<int>(sum(recv_number_of_face_per_cell_by_proc[i_rank])); - } - - parallel::exchange(cell_face_number_to_send_by_proc, recv_cell_face_number_by_proc); - - const std::unordered_map<int, int> face_number_id_map - = [&] () { - std::unordered_map<int, int> face_number_id_map; - for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) { - int cpt=0; - for (size_t i=0; i<recv_cell_face_number_by_proc[i_rank].size(); ++i) { - int face_number = recv_cell_face_number_by_proc[i_rank][i]; - auto [iterator, inserted] = face_number_id_map.insert(std::make_pair(face_number, cpt)); - if (inserted) cpt++; - } - } - return face_number_id_map; - } (); - - - new_descriptor.face_number_vector.resize(face_number_id_map.size()); - for (const auto& [number, id] : face_number_id_map) { - new_descriptor.face_number_vector[id] = number; - } - - for (size_t i_rank=0; i_rank < parallel::size(); ++i_rank) { - int l=0; - for (size_t i=0; i<recv_cell_number_by_proc[i_rank].size(); ++i) { - std::vector<unsigned int> face_vector; - for (int k=0; k<recv_number_of_face_per_cell_by_proc[i_rank][i]; ++k) { - const auto& searched_face_id = face_number_id_map.find(recv_cell_face_number_by_proc[i_rank][l++]); - Assert(searched_face_id != face_number_id_map.end()); - face_vector.push_back(searched_face_id->second); - } - new_descriptor.cell_to_face_vector.emplace_back(face_vector); - } - } - - std::vector<Array<bool>> cell_face_is_reversed_to_send_by_proc(parallel::size()); - { - const auto& cell_face_is_reversed = mesh.connectivity().cellFaceIsReversed(); - for (size_t i_rank=0; i_rank < parallel::size(); ++i_rank) { - std::vector<bool> face_is_reversed_by_cell_vector; - for (size_t j=0; j<cell_list_to_send_by_proc[i_rank].size(); ++j) { - const CellId& cell_id = cell_list_to_send_by_proc[i_rank][j]; - const auto& face_is_reversed = cell_face_is_reversed.itemValues(cell_id); - for (size_t L=0; L<face_is_reversed.size(); ++L) { - face_is_reversed_by_cell_vector.push_back(face_is_reversed[L]); - } - } - cell_face_is_reversed_to_send_by_proc[i_rank] = convert_to_array(face_is_reversed_by_cell_vector); - } - } - - std::vector<Array<bool>> recv_cell_face_is_reversed_by_proc(parallel::size()); - for (size_t i_rank=0; i_rank < parallel::size(); ++i_rank) { - recv_cell_face_is_reversed_by_proc[i_rank] - = Array<bool>(sum(recv_number_of_face_per_cell_by_proc[i_rank])); - } - - parallel::exchange(cell_face_is_reversed_to_send_by_proc, recv_cell_face_is_reversed_by_proc); - - for (size_t i_rank=0; i_rank < parallel::size(); ++i_rank) { - int l=0; - for (size_t i=0; i<recv_cell_number_by_proc[i_rank].size(); ++i) { - std::vector<bool> face_is_reversed_vector; - for (int k=0; k<recv_number_of_face_per_cell_by_proc[i_rank][i]; ++k) { - face_is_reversed_vector.push_back(recv_cell_face_is_reversed_by_proc[i_rank][l++]); - } - new_descriptor.cell_face_is_reversed_vector.emplace_back(face_is_reversed_vector); - } - } - - std::vector<Array<const FaceId>> send_face_id_by_proc(parallel::size()); - for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) { - Array<bool> tag(mesh.numberOfFaces()); - tag.fill(false); - std::vector<FaceId> face_id_vector; - for (size_t j=0; j<cell_list_to_send_by_proc[i_rank].size(); ++j) { - const CellId& cell_id = cell_list_to_send_by_proc[i_rank][j]; - const auto& cell_face_list = cell_to_face_matrix[cell_id]; - for (size_t l=0; l<cell_face_list.size(); ++l) { - const FaceId& face_id = cell_face_list[l]; - if (not tag[face_id]) { - face_id_vector.push_back(face_id); - tag[face_id] = true; - } - } - } - send_face_id_by_proc[i_rank] = convert_to_array(face_id_vector); - } - - std::vector<Array<const int>> send_face_number_by_proc(parallel::size()); - for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) { - Array<int> send_face_number(send_face_id_by_proc[i_rank].size()); - const Array<const FaceId> send_face_id = send_face_id_by_proc[i_rank]; - parallel_for(send_face_number.size(), PASTIS_LAMBDA(const size_t& j){ - send_face_number[j] = face_number[send_face_id[j]]; - }); - send_face_number_by_proc[i_rank] = send_face_number; - } - - Array<unsigned int> nb_face_to_send_by_proc(parallel::size()); - for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) { - nb_face_to_send_by_proc[i_rank] = send_face_id_by_proc[i_rank].size(); - } - Array<const unsigned int> nb_face_to_recv_by_proc - = parallel::allToAll(nb_face_to_send_by_proc); - - std::vector<Array<int>> recv_face_number_by_proc(parallel::size()); - for (size_t i_rank=0; i_rank < parallel::size(); ++i_rank) { - recv_face_number_by_proc[i_rank] = Array<int>(nb_face_to_recv_by_proc[i_rank]); - } - parallel::exchange(send_face_number_by_proc, recv_face_number_by_proc); - - std::vector<Array<const FaceId>> recv_face_id_correspondance_by_proc(parallel::size()); - for (size_t i_rank=0; i_rank<nb_face_to_recv_by_proc.size(); ++i_rank) { - Array<FaceId> face_id_correspondace(nb_face_to_recv_by_proc[i_rank]); - for (size_t l=0; l<nb_face_to_recv_by_proc[i_rank]; ++l) { - const int& face_number = recv_face_number_by_proc[i_rank][l]; - const auto& searched_face_id = face_number_id_map.find(face_number); - Assert(searched_face_id != face_number_id_map.end()); - face_id_correspondace[l] = searched_face_id->second; - } - recv_face_id_correspondance_by_proc[i_rank] = face_id_correspondace; - } - - - const FaceValue<const int>& new_face_owner = dispatcher.faceNewOwner(); - std::vector<Array<const int>> send_face_owner_by_proc(parallel::size()); - for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) { - Array<int> send_face_owner(nb_face_to_send_by_proc[i_rank]); - const Array<const FaceId> send_face_id = send_face_id_by_proc[i_rank]; - parallel_for(send_face_id.size(), PASTIS_LAMBDA(const size_t& l) { - const FaceId& face_id = send_face_id[l]; - send_face_owner[l] = new_face_owner[face_id]; - }); - send_face_owner_by_proc[i_rank] = send_face_owner; - } - - std::vector<Array<int>> recv_face_owner_by_proc(parallel::size()); - for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) { - recv_face_owner_by_proc[i_rank] = Array<int>(nb_face_to_recv_by_proc[i_rank]); - } - parallel::exchange(send_face_owner_by_proc, recv_face_owner_by_proc); - - new_descriptor.face_owner_vector.resize(new_descriptor.face_number_vector.size()); - for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) { - for (size_t r=0; r<recv_face_owner_by_proc[i_rank].size(); ++r) { - const FaceId& face_id = recv_face_id_correspondance_by_proc[i_rank][r]; - new_descriptor.face_owner_vector[face_id] = recv_face_owner_by_proc[i_rank][r]; - } - } - - new_descriptor.face_owner_vector.resize(new_descriptor.face_number_vector.size()); - for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) { - for (size_t r=0; r<recv_face_owner_by_proc[i_rank].size(); ++r) { - const FaceId& face_id = recv_face_id_correspondance_by_proc[i_rank][r]; - new_descriptor.face_owner_vector[face_id] = recv_face_owner_by_proc[i_rank][r]; - } - } - - FaceValue<int> number_of_node_per_face(mesh.connectivity()); - const auto& face_to_node_matrix = mesh.connectivity().faceToNodeMatrix(); - parallel_for(mesh.numberOfFaces(), PASTIS_LAMBDA(const FaceId& j){ - number_of_node_per_face[j] = face_to_node_matrix[j].size(); - }); - - std::vector<Array<const int>> send_face_number_of_node_by_proc(parallel::size()); - for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) { - Array<int> send_face_number_of_node(nb_face_to_send_by_proc[i_rank]); - const Array<const FaceId> send_face_id = send_face_id_by_proc[i_rank]; - parallel_for(send_face_id.size(), PASTIS_LAMBDA(const size_t& l) { - const FaceId& face_id = send_face_id[l]; - send_face_number_of_node[l] = number_of_node_per_face[face_id]; - }); - send_face_number_of_node_by_proc[i_rank] = send_face_number_of_node; - } - - std::vector<Array<int>> recv_face_number_of_node_by_proc(parallel::size()); - for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) { - recv_face_number_of_node_by_proc[i_rank] = Array<int>(nb_face_to_recv_by_proc[i_rank]); - } - parallel::exchange(send_face_number_of_node_by_proc, recv_face_number_of_node_by_proc); - - std::vector<Array<int>> face_node_number_to_send_by_proc(parallel::size()); - for (size_t i_rank=0; i_rank < parallel::size(); ++i_rank) { - std::vector<int> node_number_by_face_vector; - for (size_t l=0; l<send_face_id_by_proc[i_rank].size(); ++l) { - const FaceId& face_id = send_face_id_by_proc[i_rank][l]; - const auto& face_node_list = face_to_node_matrix[face_id]; - for (size_t r=0; r<face_node_list.size(); ++r) { - const NodeId& node_id = face_node_list[r]; - node_number_by_face_vector.push_back(node_number[node_id]); - } - } - face_node_number_to_send_by_proc[i_rank] = convert_to_array(node_number_by_face_vector); - } - - std::vector<Array<int>> recv_face_node_number_by_proc(parallel::size()); - for (size_t i_rank=0; i_rank < parallel::size(); ++i_rank) { - recv_face_node_number_by_proc[i_rank] - = Array<int>(sum(recv_face_number_of_node_by_proc[i_rank])); - } - - parallel::exchange(face_node_number_to_send_by_proc, recv_face_node_number_by_proc); - - for (size_t i_rank=0; i_rank < parallel::size(); ++i_rank) { - int l=0; - for (size_t i=0; i<recv_face_number_by_proc[i_rank].size(); ++i) { - std::vector<unsigned int> node_vector; - for (int k=0; k<recv_face_number_of_node_by_proc[i_rank][i]; ++k) { - const auto& searched_node_id = node_number_id_map.find(recv_face_node_number_by_proc[i_rank][l++]); - Assert(searched_node_id != node_number_id_map.end()); - node_vector.push_back(searched_node_id->second); - } - new_descriptor.face_to_node_vector.emplace_back(node_vector); - } - } - - // Getting references - Array<const size_t> number_of_ref_face_list_per_proc - = parallel::allGather(mesh.connectivity().numberOfRefFaceList()); - - const size_t number_of_face_list_sender - = [&] () { - size_t number_of_face_list_sender=0; - for (size_t i_rank = 0; i_rank < parallel::size(); ++i_rank) { - number_of_face_list_sender - += (number_of_ref_face_list_per_proc[i_rank] > 0); - } - return number_of_face_list_sender; - }(); - - if (number_of_face_list_sender > 0) { - if (number_of_face_list_sender > 1) { - perr() << __FILE__ << ':' << __LINE__ << ": " - << rang::fgB::red - <<"need to check that knowing procs know the same ref_face_lists!" - << rang::fg::reset << '\n'; - } - - if (number_of_face_list_sender < parallel::size()) { - const size_t sender_rank - = [&] () { - size_t i_rank = 0; - for (; i_rank < parallel::size(); ++i_rank) { - if (number_of_ref_face_list_per_proc[i_rank] > 0) { - break; - } - } - return i_rank; - }(); - - Assert(number_of_face_list_sender < parallel::size()); - - // sending references tags - Array<RefId::TagNumberType> ref_tag_list{number_of_ref_face_list_per_proc[sender_rank]}; - if (parallel::rank() == sender_rank){ - for (size_t i_ref_face_list=0; i_ref_face_list<mesh.connectivity().numberOfRefFaceList(); - ++i_ref_face_list) { - auto ref_face_list = mesh.connectivity().refFaceList(i_ref_face_list); - ref_tag_list[i_ref_face_list] = ref_face_list.refId().tagNumber(); - } - } - parallel::broadcast(ref_tag_list, sender_rank); - - // sending references name size - Array<size_t> ref_name_size_list{number_of_ref_face_list_per_proc[sender_rank]}; - if (parallel::rank() == sender_rank){ - for (size_t i_ref_face_list=0; i_ref_face_list<mesh.connectivity().numberOfRefFaceList(); - ++i_ref_face_list) { - auto ref_face_list = mesh.connectivity().refFaceList(i_ref_face_list); - ref_name_size_list[i_ref_face_list] = ref_face_list.refId().tagName().size(); - } - } - parallel::broadcast(ref_name_size_list, sender_rank); - - // sending references name size - Array<RefId::TagNameType::value_type> ref_name_cat{sum(ref_name_size_list)}; - if (parallel::rank() == sender_rank){ - size_t i_char=0; - for (size_t i_ref_face_list=0; i_ref_face_list<mesh.connectivity().numberOfRefFaceList(); - ++i_ref_face_list) { - auto ref_face_list = mesh.connectivity().refFaceList(i_ref_face_list); - for (auto c : ref_face_list.refId().tagName()) { - ref_name_cat[i_char++] = c; - } - } - } - parallel::broadcast(ref_name_cat, sender_rank); - - std::vector<RefId> ref_id_list - = [&] () { - std::vector<RefId> ref_id_list; - ref_id_list.reserve(ref_name_size_list.size()); - size_t begining=0; - for (size_t i_ref=0; i_ref < ref_name_size_list.size(); ++i_ref) { - const size_t size = ref_name_size_list[i_ref]; - ref_id_list.emplace_back(ref_tag_list[i_ref], - std::string{&(ref_name_cat[begining]), size}); - begining += size; - } - return ref_id_list; - } (); - - using block_type = int32_t; - constexpr size_t block_size = sizeof(block_type); - const size_t nb_block = ref_id_list.size()/block_size + (ref_id_list.size()%block_size != 0); - for (size_t i_block=0; i_block<nb_block; ++i_block) { - FaceValue<block_type> face_references(mesh.connectivity()); - face_references.fill(0); - - if (mesh.connectivity().numberOfRefFaceList() > 0) { - const size_t max_i_ref = std::min(ref_id_list.size(), block_size*(i_block+1)); - for (size_t i_ref=block_size*i_block, i=0; i_ref<max_i_ref; ++i_ref, ++i) { - block_type ref_bit{1<<i}; - auto ref_face_list = mesh.connectivity().refFaceList(i_ref); - - const auto& face_list = ref_face_list.faceList(); - for (size_t i_face=0; i_face<face_list.size(); ++i_face) { - const FaceId& face_id = face_list[i_face]; - face_references[face_id] |= ref_bit; - } - } - } - - std::vector<Array<const block_type>> send_face_refs_by_proc(parallel::size()); - for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) { - Array<block_type> send_face_refs(nb_face_to_send_by_proc[i_rank]); - const Array<const FaceId> send_face_id = send_face_id_by_proc[i_rank]; - parallel_for(send_face_id.size(), PASTIS_LAMBDA(const size_t& l) { - const FaceId& face_id = send_face_id[l]; - send_face_refs[l] = face_references[face_id]; - }); - send_face_refs_by_proc[i_rank] = send_face_refs; - } - - std::vector<Array<block_type>> recv_face_refs_by_proc(parallel::size()); - for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) { - recv_face_refs_by_proc[i_rank] = Array<block_type>(nb_face_to_recv_by_proc[i_rank]); - } - parallel::exchange(send_face_refs_by_proc, recv_face_refs_by_proc); - - std::vector<block_type> face_refs(new_descriptor.face_number_vector.size()); - for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) { - for (size_t r=0; r<recv_face_owner_by_proc[i_rank].size(); ++r) { - const FaceId& face_id = recv_face_id_correspondance_by_proc[i_rank][r]; - face_refs[face_id] = recv_face_refs_by_proc[i_rank][r]; - } - } - - const size_t max_i_ref = std::min(ref_id_list.size(), block_size*(i_block+1)); - for (size_t i_ref=block_size*i_block, i=0; i_ref<max_i_ref; ++i_ref, ++i) { - block_type ref_bit{1<<i}; - - std::vector<FaceId> face_id_vector; - - for (uint32_t i_face=0; i_face<face_refs.size(); ++i_face) { - const FaceId face_id{i_face}; - if (face_refs[face_id] & ref_bit) { - face_id_vector.push_back(face_id); - } - } - - Array<const FaceId> face_id_array = convert_to_array(face_id_vector); - - new_descriptor.addRefFaceList(RefFaceList(ref_id_list[i_ref], face_id_array)); - } - - pout() << __FILE__ << ':' << __LINE__ << ": remains to build lists\n"; - } - } - } - - } - - using ConnectivityType = Connectivity<Dimension>; - std::shared_ptr p_connectivity - = ConnectivityType::build(new_descriptor); - - const NodeValue<const Rd>& xr = mesh.xr(); - std::vector<Array<const Rd>> send_node_coord_by_proc(parallel::size()); - for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) { - Array<Rd> send_node_coord(nb_node_to_send_by_proc[i_rank]); - const Array<const NodeId> send_node_id = send_node_id_by_proc[i_rank]; - parallel_for(send_node_id.size(), PASTIS_LAMBDA(const size_t& r) { - const NodeId& node_id = send_node_id[r]; - send_node_coord[r] = xr[node_id]; - }); - send_node_coord_by_proc[i_rank] = send_node_coord; - } - - std::vector<Array<Rd>> recv_node_coord_by_proc(parallel::size()); - for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) { - recv_node_coord_by_proc[i_rank] = Array<Rd>(nb_node_to_recv_by_proc[i_rank]); - } - parallel::exchange(send_node_coord_by_proc, recv_node_coord_by_proc); - - Array<const size_t> number_of_ref_node_list_per_proc - = parallel::allGather(mesh.connectivity().numberOfRefNodeList()); - if (max(number_of_ref_node_list_per_proc) > 0) { - perr() << __FILE__ << ':' << __LINE__ << ": " - << rang::fgB::red << "exchange of node ref list not implemented yet!" - << rang::fg::reset << '\n'; - } - - // Finally build the mesh - NodeValue<Rd> new_xr(*p_connectivity); - for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) { - for (size_t r=0; r<recv_node_coord_by_proc[i_rank].size(); ++r) { - const NodeId& node_id = recv_node_id_correspondance_by_proc[i_rank][r]; - new_xr[node_id] = recv_node_coord_by_proc[i_rank][r]; - } - } - - m_mesh = std::make_shared<MeshType>(p_connectivity, new_xr); -} - - GmshReader::GmshReader(const std::string& filename) : m_filename(filename) { diff --git a/src/mesh/MeshDispatcher.cpp b/src/mesh/MeshDispatcher.cpp new file mode 100644 index 0000000000000000000000000000000000000000..970a3e2813ac239de389ca58c179da7183744011 --- /dev/null +++ b/src/mesh/MeshDispatcher.cpp @@ -0,0 +1,779 @@ +#include <MeshDispatcher.hpp> +#include <Partitioner.hpp> + +#include <unordered_map> + +template <int Dimension> +CellValue<int> +MeshDispatcher<Dimension>::_getCellNewOwner() +{ + CSRGraph mesh_graph = m_mesh.cellToCellGraph(); + Partitioner P; + + CellValue<int> cell_new_owner(m_mesh.connectivity()); + cell_new_owner = P.partition(mesh_graph); + return cell_new_owner; +} + +template <int Dimension> +FaceValue<int> +MeshDispatcher<Dimension>::_getFaceNewOwner() +{ + const auto& face_to_cell_matrix + = m_mesh.connectivity().faceToCellMatrix(); + const auto& cell_number = m_mesh.connectivity().cellNumber(); + + FaceValue<int> face_new_owner(m_mesh.connectivity()); + parallel_for(m_mesh.numberOfFaces(), PASTIS_LAMBDA(const FaceId& l) { + const auto& face_to_cell = face_to_cell_matrix[l]; + CellId Jmin = face_to_cell[0]; + + for (size_t j=1; j<face_to_cell.size(); ++j) { + const CellId J = face_to_cell[j]; + if (cell_number[J] < cell_number[Jmin]) { + Jmin=J; + } + } + face_new_owner[l] = m_cell_new_owner[Jmin]; + }); + + synchronize(face_new_owner); + return face_new_owner; +} + +template <int Dimension> +NodeValue<int> +MeshDispatcher<Dimension>::_getNodeNewOwner() +{ + const auto& node_to_cell_matrix + = m_mesh.connectivity().nodeToCellMatrix(); + const auto& cell_number = m_mesh.connectivity().cellNumber(); + + NodeValue<int> node_new_owner(m_mesh.connectivity()); + parallel_for(m_mesh.numberOfNodes(), PASTIS_LAMBDA(const NodeId& r) { + const auto& node_to_cell = node_to_cell_matrix[r]; + CellId Jmin = node_to_cell[0]; + + for (size_t j=1; j<node_to_cell.size(); ++j) { + const CellId J = node_to_cell[j]; + if (cell_number[J] < cell_number[Jmin]) { + Jmin=J; + } + } + node_new_owner[r] = m_cell_new_owner[Jmin]; + }); + + synchronize(node_new_owner); + return node_new_owner; +} + +template <int Dimension> +const typename MeshDispatcher<Dimension>::CellListToSendByProc +MeshDispatcher<Dimension>::_buildCellListToSend() const +{ + const auto& node_to_cell_matrix + = m_mesh.connectivity().nodeToCellMatrix(); + const auto& cell_to_node_matrix + = m_mesh.connectivity().cellToNodeMatrix(); + + std::vector<std::vector<CellId>> cell_vector_to_send_by_proc(parallel::size()); + Array<bool> send_to_rank(parallel::size()); + for (CellId j=0; j<m_mesh.numberOfCells(); ++j) { + send_to_rank.fill(false); + const auto& cell_to_node = cell_to_node_matrix[j]; + + for (size_t R=0; R<cell_to_node.size(); ++R) { + const NodeId& r = cell_to_node[R]; + const auto& node_to_cell = node_to_cell_matrix[r]; + for (size_t K=0; K<node_to_cell.size(); ++K) { + const CellId& k = node_to_cell[K]; + send_to_rank[m_cell_new_owner[k]] = true; + } + } + + for (size_t k=0; k<send_to_rank.size(); ++k) { + if (send_to_rank[k]) { + cell_vector_to_send_by_proc[k].push_back(j); + } + } + } + + std::vector<Array<const CellId>> cell_list_to_send_by_proc(parallel::size()); + for (size_t i=0; i<parallel::size(); ++i) { + cell_list_to_send_by_proc[i] = convert_to_array(cell_vector_to_send_by_proc[i]); + } + + return cell_list_to_send_by_proc; +} + +template <int Dimension> +Array<int> +MeshDispatcher<Dimension>::_buildNbCellToSend() +{ + Array<int> nb_cell_to_send_by_proc(parallel::size()); + for (size_t i=0; i<parallel::size(); ++i) { + nb_cell_to_send_by_proc[i] = m_cell_list_to_send_by_proc[i].size(); + } + return nb_cell_to_send_by_proc; +} + +template <int Dimension> +MeshDispatcher<Dimension>::MeshDispatcher(const MeshType& mesh) + : m_mesh(mesh), + m_cell_new_owner(_getCellNewOwner()), + m_face_new_owner(_getFaceNewOwner()), + m_node_new_owner(_getNodeNewOwner()), + m_cell_list_to_send_by_proc(_buildCellListToSend()), + m_nb_cell_to_send_by_proc(_buildNbCellToSend()), + m_nb_cell_to_recv_by_proc(parallel::allToAll(m_nb_cell_to_send_by_proc)) +{ + std::vector<Array<int>> recv_cell_number_by_proc + = this->exchange(mesh.connectivity().cellNumber()); + + std::vector<Array<CellType>> recv_cell_type_by_proc + = this->exchange(mesh.connectivity().cellType()); + + std::vector<Array<int>> recv_cell_new_owner_by_proc + = this->exchange(m_cell_new_owner); + + CellValue<int> number_of_node_per_cell(mesh.connectivity()); + + const auto& cell_to_node_matrix = mesh.connectivity().cellToNodeMatrix(); + parallel_for(mesh.numberOfCells(), PASTIS_LAMBDA(const CellId& j){ + number_of_node_per_cell[j] = cell_to_node_matrix[j].size(); + }); + + std::vector<Array<int>> recv_number_of_node_per_cell_by_proc + = this->exchange(number_of_node_per_cell); + + const NodeValue<const int>& node_number = mesh.connectivity().nodeNumber(); + std::vector<Array<int>> cell_node_number_to_send_by_proc(parallel::size()); + for (size_t i_rank=0; i_rank < parallel::size(); ++i_rank) { + std::vector<int> node_number_by_cell_vector; + for (size_t j=0; j<m_cell_list_to_send_by_proc[i_rank].size(); ++j) { + const CellId& cell_id = m_cell_list_to_send_by_proc[i_rank][j]; + const auto& cell_node_list = cell_to_node_matrix[cell_id]; + for (size_t r=0; r<cell_node_list.size(); ++r) { + const NodeId& node_id = cell_node_list[r]; + node_number_by_cell_vector.push_back(node_number[node_id]); + } + } + cell_node_number_to_send_by_proc[i_rank] = convert_to_array(node_number_by_cell_vector); + } + + std::vector<Array<int>> recv_cell_node_number_by_proc(parallel::size()); + for (size_t i_rank=0; i_rank < parallel::size(); ++i_rank) { + recv_cell_node_number_by_proc[i_rank] + = Array<int>(sum(recv_number_of_node_per_cell_by_proc[i_rank])); + } + + parallel::exchange(cell_node_number_to_send_by_proc, recv_cell_node_number_by_proc); + + const std::unordered_map<int, int> node_number_id_map + = [&] () { + std::unordered_map<int, int> node_number_id_map; + for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) { + int cpt=0; + for (size_t i=0; i<recv_cell_node_number_by_proc[i_rank].size(); ++i) { + int node_number = recv_cell_node_number_by_proc[i_rank][i]; + auto [iterator, inserted] = node_number_id_map.insert(std::make_pair(node_number, cpt)); + if (inserted) cpt++; + } + } + return node_number_id_map; + } (); + + + ConnectivityDescriptor new_descriptor; + + new_descriptor.node_number_vector.resize(node_number_id_map.size()); + for (const auto& [number, id] : node_number_id_map) { + new_descriptor.node_number_vector[id] = number; + } + + for (size_t i_rank=0; i_rank < parallel::size(); ++i_rank) { + int l=0; + for (size_t i=0; i<recv_cell_number_by_proc[i_rank].size(); ++i) { + std::vector<unsigned int> node_vector; + for (int k=0; k<recv_number_of_node_per_cell_by_proc[i_rank][i]; ++k) { + const auto& searched_node_id = node_number_id_map.find(recv_cell_node_number_by_proc[i_rank][l++]); + Assert(searched_node_id != node_number_id_map.end()); + node_vector.push_back(searched_node_id->second); + } + new_descriptor.cell_by_node_vector.emplace_back(node_vector); + } + } + + const std::unordered_map<int, int> cell_number_id_map + = [&] () { + std::unordered_map<int, int> cell_number_id_map; + for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) { + int cpt=0; + for (size_t i=0; i<recv_cell_number_by_proc[i_rank].size(); ++i) { + int cell_number = recv_cell_number_by_proc[i_rank][i]; + auto [iterator, inserted] = cell_number_id_map.insert(std::make_pair(cell_number, cpt)); + if (inserted) cpt++; + } + } + return cell_number_id_map; + } (); + + new_descriptor.cell_number_vector.resize(cell_number_id_map.size()); + for (const auto& [number, id] : cell_number_id_map) { + new_descriptor.cell_number_vector[id] = number; + } + + new_descriptor.cell_type_vector.resize(cell_number_id_map.size()); + for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) { + for (size_t i=0; i<recv_cell_number_by_proc[i_rank].size(); ++i) { + int cell_number = recv_cell_number_by_proc[i_rank][i]; + const auto& searched_cell_id = cell_number_id_map.find(cell_number); + Assert(searched_cell_id != cell_number_id_map.end()); + new_descriptor.cell_type_vector[searched_cell_id->second] = recv_cell_type_by_proc[i_rank][i]; + } + } + + std::vector<Array<const NodeId>> send_node_id_by_proc(parallel::size()); + for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) { + Array<bool> tag(mesh.numberOfNodes()); + tag.fill(false); + std::vector<NodeId> node_id_vector; + for (size_t j=0; j<m_cell_list_to_send_by_proc[i_rank].size(); ++j) { + const CellId& cell_id = m_cell_list_to_send_by_proc[i_rank][j]; + const auto& cell_node_list = cell_to_node_matrix[cell_id]; + for (size_t r=0; r<cell_node_list.size(); ++r) { + const NodeId& node_id = cell_node_list[r]; + if (not tag[node_id]) { + node_id_vector.push_back(node_id); + tag[node_id] = true; + } + } + } + send_node_id_by_proc[i_rank] = convert_to_array(node_id_vector); + } + + std::vector<Array<const int>> send_node_number_by_proc(parallel::size()); + for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) { + Array<int> send_node_number(send_node_id_by_proc[i_rank].size()); + const Array<const NodeId> send_node_id = send_node_id_by_proc[i_rank]; + parallel_for(send_node_number.size(), PASTIS_LAMBDA(const size_t& j){ + send_node_number[j] = node_number[send_node_id[j]]; + }); + send_node_number_by_proc[i_rank] = send_node_number; + } + + Array<unsigned int> nb_node_to_send_by_proc(parallel::size()); + for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) { + nb_node_to_send_by_proc[i_rank] = send_node_id_by_proc[i_rank].size(); + } + Array<const unsigned int> nb_node_to_recv_by_proc + = parallel::allToAll(nb_node_to_send_by_proc); + + std::vector<Array<int>> recv_node_number_by_proc(parallel::size()); + for (size_t i_rank=0; i_rank < parallel::size(); ++i_rank) { + recv_node_number_by_proc[i_rank] = Array<int>(nb_node_to_recv_by_proc[i_rank]); + } + parallel::exchange(send_node_number_by_proc, recv_node_number_by_proc); + + std::vector<Array<const NodeId>> recv_node_id_correspondance_by_proc(parallel::size()); + for (size_t i_rank=0; i_rank<nb_node_to_recv_by_proc.size(); ++i_rank) { + Array<NodeId> node_id_correspondace(nb_node_to_recv_by_proc[i_rank]); + for (size_t l=0; l<nb_node_to_recv_by_proc[i_rank]; ++l) { + const int& node_number = recv_node_number_by_proc[i_rank][l]; + const auto& searched_node_id = node_number_id_map.find(node_number); + Assert(searched_node_id != node_number_id_map.end()); + node_id_correspondace[l] = searched_node_id->second; + } + recv_node_id_correspondance_by_proc[i_rank] = node_id_correspondace; + } + + new_descriptor.cell_owner_vector.resize(cell_number_id_map.size()); + for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) { + for (size_t i=0; i<recv_cell_number_by_proc[i_rank].size(); ++i) { + int cell_number = recv_cell_number_by_proc[i_rank][i]; + const auto& searched_cell_id = cell_number_id_map.find(cell_number); + Assert(searched_cell_id != cell_number_id_map.end()); + new_descriptor.cell_owner_vector[searched_cell_id->second] = recv_cell_new_owner_by_proc[i_rank][i]; + } + } + + std::vector<Array<const int>> send_node_owner_by_proc(parallel::size()); + for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) { + Array<int> send_node_owner(nb_node_to_send_by_proc[i_rank]); + const Array<const NodeId> send_node_id = send_node_id_by_proc[i_rank]; + parallel_for(send_node_id.size(), PASTIS_LAMBDA(const size_t& r) { + const NodeId& node_id = send_node_id[r]; + send_node_owner[r] = m_node_new_owner[node_id]; + }); + send_node_owner_by_proc[i_rank] = send_node_owner; + } + + std::vector<Array<int>> recv_node_owner_by_proc(parallel::size()); + for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) { + recv_node_owner_by_proc[i_rank] = Array<int>(nb_node_to_recv_by_proc[i_rank]); + } + parallel::exchange(send_node_owner_by_proc, recv_node_owner_by_proc); + + new_descriptor.node_owner_vector.resize(new_descriptor.node_number_vector.size()); + for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) { + for (size_t r=0; r<recv_node_owner_by_proc[i_rank].size(); ++r) { + const NodeId& node_id = recv_node_id_correspondance_by_proc[i_rank][r]; + new_descriptor.node_owner_vector[node_id] = recv_node_owner_by_proc[i_rank][r]; + } + } + + if constexpr(Dimension>1) { // Faces + CellValue<int> number_of_face_per_cell(mesh.connectivity()); + const auto& cell_to_face_matrix = mesh.connectivity().cellToFaceMatrix(); + parallel_for(mesh.numberOfCells(), PASTIS_LAMBDA(const CellId& j){ + number_of_face_per_cell[j] = cell_to_face_matrix[j].size(); + }); + + std::vector<Array<int>> recv_number_of_face_per_cell_by_proc + = this->exchange(number_of_face_per_cell); + + + const FaceValue<const int>& face_number = mesh.connectivity().faceNumber(); + std::vector<Array<int>> cell_face_number_to_send_by_proc(parallel::size()); + for (size_t i_rank=0; i_rank < parallel::size(); ++i_rank) { + std::vector<int> face_number_by_cell_vector; + for (size_t j=0; j<m_cell_list_to_send_by_proc[i_rank].size(); ++j) { + const CellId& cell_id = m_cell_list_to_send_by_proc[i_rank][j]; + const auto& cell_face_list = cell_to_face_matrix[cell_id]; + for (size_t r=0; r<cell_face_list.size(); ++r) { + const FaceId& face_id = cell_face_list[r]; + face_number_by_cell_vector.push_back(face_number[face_id]); + } + } + cell_face_number_to_send_by_proc[i_rank] = convert_to_array(face_number_by_cell_vector); + } + + + std::vector<Array<int>> recv_cell_face_number_by_proc(parallel::size()); + for (size_t i_rank=0; i_rank < parallel::size(); ++i_rank) { + recv_cell_face_number_by_proc[i_rank] + = Array<int>(sum(recv_number_of_face_per_cell_by_proc[i_rank])); + } + + parallel::exchange(cell_face_number_to_send_by_proc, recv_cell_face_number_by_proc); + + const std::unordered_map<int, int> face_number_id_map + = [&] () { + std::unordered_map<int, int> face_number_id_map; + for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) { + int cpt=0; + for (size_t i=0; i<recv_cell_face_number_by_proc[i_rank].size(); ++i) { + int face_number = recv_cell_face_number_by_proc[i_rank][i]; + auto [iterator, inserted] = face_number_id_map.insert(std::make_pair(face_number, cpt)); + if (inserted) cpt++; + } + } + return face_number_id_map; + } (); + + + new_descriptor.face_number_vector.resize(face_number_id_map.size()); + for (const auto& [number, id] : face_number_id_map) { + new_descriptor.face_number_vector[id] = number; + } + + for (size_t i_rank=0; i_rank < parallel::size(); ++i_rank) { + int l=0; + for (size_t i=0; i<recv_cell_number_by_proc[i_rank].size(); ++i) { + std::vector<unsigned int> face_vector; + for (int k=0; k<recv_number_of_face_per_cell_by_proc[i_rank][i]; ++k) { + const auto& searched_face_id = face_number_id_map.find(recv_cell_face_number_by_proc[i_rank][l++]); + Assert(searched_face_id != face_number_id_map.end()); + face_vector.push_back(searched_face_id->second); + } + new_descriptor.cell_to_face_vector.emplace_back(face_vector); + } + } + + std::vector<Array<bool>> cell_face_is_reversed_to_send_by_proc(parallel::size()); + { + const auto& cell_face_is_reversed = mesh.connectivity().cellFaceIsReversed(); + for (size_t i_rank=0; i_rank < parallel::size(); ++i_rank) { + std::vector<bool> face_is_reversed_by_cell_vector; + for (size_t j=0; j<m_cell_list_to_send_by_proc[i_rank].size(); ++j) { + const CellId& cell_id = m_cell_list_to_send_by_proc[i_rank][j]; + const auto& face_is_reversed = cell_face_is_reversed.itemValues(cell_id); + for (size_t L=0; L<face_is_reversed.size(); ++L) { + face_is_reversed_by_cell_vector.push_back(face_is_reversed[L]); + } + } + cell_face_is_reversed_to_send_by_proc[i_rank] = convert_to_array(face_is_reversed_by_cell_vector); + } + } + + std::vector<Array<bool>> recv_cell_face_is_reversed_by_proc(parallel::size()); + for (size_t i_rank=0; i_rank < parallel::size(); ++i_rank) { + recv_cell_face_is_reversed_by_proc[i_rank] + = Array<bool>(sum(recv_number_of_face_per_cell_by_proc[i_rank])); + } + + parallel::exchange(cell_face_is_reversed_to_send_by_proc, recv_cell_face_is_reversed_by_proc); + + for (size_t i_rank=0; i_rank < parallel::size(); ++i_rank) { + int l=0; + for (size_t i=0; i<recv_cell_number_by_proc[i_rank].size(); ++i) { + std::vector<bool> face_is_reversed_vector; + for (int k=0; k<recv_number_of_face_per_cell_by_proc[i_rank][i]; ++k) { + face_is_reversed_vector.push_back(recv_cell_face_is_reversed_by_proc[i_rank][l++]); + } + new_descriptor.cell_face_is_reversed_vector.emplace_back(face_is_reversed_vector); + } + } + + std::vector<Array<const FaceId>> send_face_id_by_proc(parallel::size()); + for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) { + Array<bool> tag(mesh.numberOfFaces()); + tag.fill(false); + std::vector<FaceId> face_id_vector; + for (size_t j=0; j<m_cell_list_to_send_by_proc[i_rank].size(); ++j) { + const CellId& cell_id = m_cell_list_to_send_by_proc[i_rank][j]; + const auto& cell_face_list = cell_to_face_matrix[cell_id]; + for (size_t l=0; l<cell_face_list.size(); ++l) { + const FaceId& face_id = cell_face_list[l]; + if (not tag[face_id]) { + face_id_vector.push_back(face_id); + tag[face_id] = true; + } + } + } + send_face_id_by_proc[i_rank] = convert_to_array(face_id_vector); + } + + std::vector<Array<const int>> send_face_number_by_proc(parallel::size()); + for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) { + Array<int> send_face_number(send_face_id_by_proc[i_rank].size()); + const Array<const FaceId> send_face_id = send_face_id_by_proc[i_rank]; + parallel_for(send_face_number.size(), PASTIS_LAMBDA(const size_t& j){ + send_face_number[j] = face_number[send_face_id[j]]; + }); + send_face_number_by_proc[i_rank] = send_face_number; + } + + Array<unsigned int> nb_face_to_send_by_proc(parallel::size()); + for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) { + nb_face_to_send_by_proc[i_rank] = send_face_id_by_proc[i_rank].size(); + } + Array<const unsigned int> nb_face_to_recv_by_proc + = parallel::allToAll(nb_face_to_send_by_proc); + + std::vector<Array<int>> recv_face_number_by_proc(parallel::size()); + for (size_t i_rank=0; i_rank < parallel::size(); ++i_rank) { + recv_face_number_by_proc[i_rank] = Array<int>(nb_face_to_recv_by_proc[i_rank]); + } + parallel::exchange(send_face_number_by_proc, recv_face_number_by_proc); + + std::vector<Array<const FaceId>> recv_face_id_correspondance_by_proc(parallel::size()); + for (size_t i_rank=0; i_rank<nb_face_to_recv_by_proc.size(); ++i_rank) { + Array<FaceId> face_id_correspondace(nb_face_to_recv_by_proc[i_rank]); + for (size_t l=0; l<nb_face_to_recv_by_proc[i_rank]; ++l) { + const int& face_number = recv_face_number_by_proc[i_rank][l]; + const auto& searched_face_id = face_number_id_map.find(face_number); + Assert(searched_face_id != face_number_id_map.end()); + face_id_correspondace[l] = searched_face_id->second; + } + recv_face_id_correspondance_by_proc[i_rank] = face_id_correspondace; + } + + + std::vector<Array<const int>> send_face_owner_by_proc(parallel::size()); + for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) { + Array<int> send_face_owner(nb_face_to_send_by_proc[i_rank]); + const Array<const FaceId> send_face_id = send_face_id_by_proc[i_rank]; + parallel_for(send_face_id.size(), PASTIS_LAMBDA(const size_t& l) { + const FaceId& face_id = send_face_id[l]; + send_face_owner[l] = m_face_new_owner[face_id]; + }); + send_face_owner_by_proc[i_rank] = send_face_owner; + } + + std::vector<Array<int>> recv_face_owner_by_proc(parallel::size()); + for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) { + recv_face_owner_by_proc[i_rank] = Array<int>(nb_face_to_recv_by_proc[i_rank]); + } + parallel::exchange(send_face_owner_by_proc, recv_face_owner_by_proc); + + new_descriptor.face_owner_vector.resize(new_descriptor.face_number_vector.size()); + for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) { + for (size_t r=0; r<recv_face_owner_by_proc[i_rank].size(); ++r) { + const FaceId& face_id = recv_face_id_correspondance_by_proc[i_rank][r]; + new_descriptor.face_owner_vector[face_id] = recv_face_owner_by_proc[i_rank][r]; + } + } + + new_descriptor.face_owner_vector.resize(new_descriptor.face_number_vector.size()); + for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) { + for (size_t r=0; r<recv_face_owner_by_proc[i_rank].size(); ++r) { + const FaceId& face_id = recv_face_id_correspondance_by_proc[i_rank][r]; + new_descriptor.face_owner_vector[face_id] = recv_face_owner_by_proc[i_rank][r]; + } + } + + FaceValue<int> number_of_node_per_face(mesh.connectivity()); + const auto& face_to_node_matrix = mesh.connectivity().faceToNodeMatrix(); + parallel_for(mesh.numberOfFaces(), PASTIS_LAMBDA(const FaceId& j){ + number_of_node_per_face[j] = face_to_node_matrix[j].size(); + }); + + std::vector<Array<const int>> send_face_number_of_node_by_proc(parallel::size()); + for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) { + Array<int> send_face_number_of_node(nb_face_to_send_by_proc[i_rank]); + const Array<const FaceId> send_face_id = send_face_id_by_proc[i_rank]; + parallel_for(send_face_id.size(), PASTIS_LAMBDA(const size_t& l) { + const FaceId& face_id = send_face_id[l]; + send_face_number_of_node[l] = number_of_node_per_face[face_id]; + }); + send_face_number_of_node_by_proc[i_rank] = send_face_number_of_node; + } + + std::vector<Array<int>> recv_face_number_of_node_by_proc(parallel::size()); + for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) { + recv_face_number_of_node_by_proc[i_rank] = Array<int>(nb_face_to_recv_by_proc[i_rank]); + } + parallel::exchange(send_face_number_of_node_by_proc, recv_face_number_of_node_by_proc); + + std::vector<Array<int>> face_node_number_to_send_by_proc(parallel::size()); + for (size_t i_rank=0; i_rank < parallel::size(); ++i_rank) { + std::vector<int> node_number_by_face_vector; + for (size_t l=0; l<send_face_id_by_proc[i_rank].size(); ++l) { + const FaceId& face_id = send_face_id_by_proc[i_rank][l]; + const auto& face_node_list = face_to_node_matrix[face_id]; + for (size_t r=0; r<face_node_list.size(); ++r) { + const NodeId& node_id = face_node_list[r]; + node_number_by_face_vector.push_back(node_number[node_id]); + } + } + face_node_number_to_send_by_proc[i_rank] = convert_to_array(node_number_by_face_vector); + } + + std::vector<Array<int>> recv_face_node_number_by_proc(parallel::size()); + for (size_t i_rank=0; i_rank < parallel::size(); ++i_rank) { + recv_face_node_number_by_proc[i_rank] + = Array<int>(sum(recv_face_number_of_node_by_proc[i_rank])); + } + + parallel::exchange(face_node_number_to_send_by_proc, recv_face_node_number_by_proc); + + for (size_t i_rank=0; i_rank < parallel::size(); ++i_rank) { + int l=0; + for (size_t i=0; i<recv_face_number_by_proc[i_rank].size(); ++i) { + std::vector<unsigned int> node_vector; + for (int k=0; k<recv_face_number_of_node_by_proc[i_rank][i]; ++k) { + const auto& searched_node_id = node_number_id_map.find(recv_face_node_number_by_proc[i_rank][l++]); + Assert(searched_node_id != node_number_id_map.end()); + node_vector.push_back(searched_node_id->second); + } + new_descriptor.face_to_node_vector.emplace_back(node_vector); + } + } + + // Getting references + Array<const size_t> number_of_ref_face_list_per_proc + = parallel::allGather(mesh.connectivity().numberOfRefFaceList()); + + const size_t number_of_face_list_sender + = [&] () { + size_t number_of_face_list_sender=0; + for (size_t i_rank = 0; i_rank < parallel::size(); ++i_rank) { + number_of_face_list_sender + += (number_of_ref_face_list_per_proc[i_rank] > 0); + } + return number_of_face_list_sender; + }(); + + if (number_of_face_list_sender > 0) { + if (number_of_face_list_sender > 1) { + perr() << __FILE__ << ':' << __LINE__ << ": " + << rang::fgB::red + <<"need to check that knowing procs know the same ref_face_lists!" + << rang::fg::reset << '\n'; + } + + if (number_of_face_list_sender < parallel::size()) { + const size_t sender_rank + = [&] () { + size_t i_rank = 0; + for (; i_rank < parallel::size(); ++i_rank) { + if (number_of_ref_face_list_per_proc[i_rank] > 0) { + break; + } + } + return i_rank; + }(); + + Assert(number_of_face_list_sender < parallel::size()); + + // sending references tags + Array<RefId::TagNumberType> ref_tag_list{number_of_ref_face_list_per_proc[sender_rank]}; + if (parallel::rank() == sender_rank){ + for (size_t i_ref_face_list=0; i_ref_face_list<mesh.connectivity().numberOfRefFaceList(); + ++i_ref_face_list) { + auto ref_face_list = mesh.connectivity().refFaceList(i_ref_face_list); + ref_tag_list[i_ref_face_list] = ref_face_list.refId().tagNumber(); + } + } + parallel::broadcast(ref_tag_list, sender_rank); + + // sending references name size + Array<size_t> ref_name_size_list{number_of_ref_face_list_per_proc[sender_rank]}; + if (parallel::rank() == sender_rank){ + for (size_t i_ref_face_list=0; i_ref_face_list<mesh.connectivity().numberOfRefFaceList(); + ++i_ref_face_list) { + auto ref_face_list = mesh.connectivity().refFaceList(i_ref_face_list); + ref_name_size_list[i_ref_face_list] = ref_face_list.refId().tagName().size(); + } + } + parallel::broadcast(ref_name_size_list, sender_rank); + + // sending references name size + Array<RefId::TagNameType::value_type> ref_name_cat{sum(ref_name_size_list)}; + if (parallel::rank() == sender_rank){ + size_t i_char=0; + for (size_t i_ref_face_list=0; i_ref_face_list<mesh.connectivity().numberOfRefFaceList(); + ++i_ref_face_list) { + auto ref_face_list = mesh.connectivity().refFaceList(i_ref_face_list); + for (auto c : ref_face_list.refId().tagName()) { + ref_name_cat[i_char++] = c; + } + } + } + parallel::broadcast(ref_name_cat, sender_rank); + + std::vector<RefId> ref_id_list + = [&] () { + std::vector<RefId> ref_id_list; + ref_id_list.reserve(ref_name_size_list.size()); + size_t begining=0; + for (size_t i_ref=0; i_ref < ref_name_size_list.size(); ++i_ref) { + const size_t size = ref_name_size_list[i_ref]; + ref_id_list.emplace_back(ref_tag_list[i_ref], + std::string{&(ref_name_cat[begining]), size}); + begining += size; + } + return ref_id_list; + } (); + + using block_type = int32_t; + constexpr size_t block_size = sizeof(block_type); + const size_t nb_block = ref_id_list.size()/block_size + (ref_id_list.size()%block_size != 0); + for (size_t i_block=0; i_block<nb_block; ++i_block) { + FaceValue<block_type> face_references(mesh.connectivity()); + face_references.fill(0); + + if (mesh.connectivity().numberOfRefFaceList() > 0) { + const size_t max_i_ref = std::min(ref_id_list.size(), block_size*(i_block+1)); + for (size_t i_ref=block_size*i_block, i=0; i_ref<max_i_ref; ++i_ref, ++i) { + block_type ref_bit{1<<i}; + auto ref_face_list = mesh.connectivity().refFaceList(i_ref); + + const auto& face_list = ref_face_list.faceList(); + for (size_t i_face=0; i_face<face_list.size(); ++i_face) { + const FaceId& face_id = face_list[i_face]; + face_references[face_id] |= ref_bit; + } + } + } + + std::vector<Array<const block_type>> send_face_refs_by_proc(parallel::size()); + for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) { + Array<block_type> send_face_refs(nb_face_to_send_by_proc[i_rank]); + const Array<const FaceId> send_face_id = send_face_id_by_proc[i_rank]; + parallel_for(send_face_id.size(), PASTIS_LAMBDA(const size_t& l) { + const FaceId& face_id = send_face_id[l]; + send_face_refs[l] = face_references[face_id]; + }); + send_face_refs_by_proc[i_rank] = send_face_refs; + } + + std::vector<Array<block_type>> recv_face_refs_by_proc(parallel::size()); + for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) { + recv_face_refs_by_proc[i_rank] = Array<block_type>(nb_face_to_recv_by_proc[i_rank]); + } + parallel::exchange(send_face_refs_by_proc, recv_face_refs_by_proc); + + std::vector<block_type> face_refs(new_descriptor.face_number_vector.size()); + for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) { + for (size_t r=0; r<recv_face_owner_by_proc[i_rank].size(); ++r) { + const FaceId& face_id = recv_face_id_correspondance_by_proc[i_rank][r]; + face_refs[face_id] = recv_face_refs_by_proc[i_rank][r]; + } + } + + const size_t max_i_ref = std::min(ref_id_list.size(), block_size*(i_block+1)); + for (size_t i_ref=block_size*i_block, i=0; i_ref<max_i_ref; ++i_ref, ++i) { + block_type ref_bit{1<<i}; + + std::vector<FaceId> face_id_vector; + + for (uint32_t i_face=0; i_face<face_refs.size(); ++i_face) { + const FaceId face_id{i_face}; + if (face_refs[face_id] & ref_bit) { + face_id_vector.push_back(face_id); + } + } + + Array<const FaceId> face_id_array = convert_to_array(face_id_vector); + + new_descriptor.addRefFaceList(RefFaceList(ref_id_list[i_ref], face_id_array)); + } + + pout() << __FILE__ << ':' << __LINE__ << ": remains to build lists\n"; + } + } + } + + } + + using ConnectivityType = Connectivity<Dimension>; + std::shared_ptr p_connectivity + = ConnectivityType::build(new_descriptor); + + using Rd = TinyVector<Dimension, double>; + + const NodeValue<const Rd>& xr = m_mesh.xr(); + std::vector<Array<const Rd>> send_node_coord_by_proc(parallel::size()); + for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) { + Array<Rd> send_node_coord(nb_node_to_send_by_proc[i_rank]); + const Array<const NodeId> send_node_id = send_node_id_by_proc[i_rank]; + parallel_for(send_node_id.size(), PASTIS_LAMBDA(const size_t& r) { + const NodeId& node_id = send_node_id[r]; + send_node_coord[r] = xr[node_id]; + }); + send_node_coord_by_proc[i_rank] = send_node_coord; + } + + std::vector<Array<Rd>> recv_node_coord_by_proc(parallel::size()); + for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) { + recv_node_coord_by_proc[i_rank] = Array<Rd>(nb_node_to_recv_by_proc[i_rank]); + } + parallel::exchange(send_node_coord_by_proc, recv_node_coord_by_proc); + + Array<const size_t> number_of_ref_node_list_per_proc + = parallel::allGather(m_mesh.connectivity().numberOfRefNodeList()); + if (max(number_of_ref_node_list_per_proc) > 0) { + perr() << __FILE__ << ':' << __LINE__ << ": " + << rang::fgB::red << "exchange of node ref list not implemented yet!" + << rang::fg::reset << '\n'; + } + + // Finally build the mesh + NodeValue<Rd> new_xr(*p_connectivity); + for (size_t i_rank=0; i_rank<parallel::size(); ++i_rank) { + for (size_t r=0; r<recv_node_coord_by_proc[i_rank].size(); ++r) { + const NodeId& node_id = recv_node_id_correspondance_by_proc[i_rank][r]; + new_xr[node_id] = recv_node_coord_by_proc[i_rank][r]; + } + } + + m_dispatched_mesh = std::make_shared<MeshType>(p_connectivity, new_xr); +} + +template MeshDispatcher<1>::MeshDispatcher(const MeshType& mesh); + +template MeshDispatcher<2>::MeshDispatcher(const MeshType& mesh); + +template MeshDispatcher<3>::MeshDispatcher(const MeshType& mesh); diff --git a/src/mesh/MeshDispatcher.hpp b/src/mesh/MeshDispatcher.hpp new file mode 100644 index 0000000000000000000000000000000000000000..c4d27f086689b82854393bf6f53c4aaeb2337c39 --- /dev/null +++ b/src/mesh/MeshDispatcher.hpp @@ -0,0 +1,98 @@ +#ifndef MESH_DISPATCHER_HPP +#define MESH_DISPATCHER_HPP + +#include <Mesh.hpp> +#include <ItemValue.hpp> +#include <ItemValueUtils.hpp> + +template <int Dimension> +class MeshDispatcher +{ + public: + using MeshType = Mesh<Connectivity<Dimension>>; + + private: + const MeshType& m_mesh; + + std::shared_ptr<MeshType> m_dispatched_mesh; + + CellValue<const int> m_cell_new_owner; + FaceValue<const int> m_face_new_owner; + NodeValue<const int> m_node_new_owner; + + using CellListToSendByProc = std::vector<Array<const CellId>>; + const CellListToSendByProc m_cell_list_to_send_by_proc; + + Array<int> m_nb_cell_to_send_by_proc; + Array<int> m_nb_cell_to_recv_by_proc; + + CellValue<int> _getCellNewOwner(); + FaceValue<int> _getFaceNewOwner(); + NodeValue<int> _getNodeNewOwner(); + + const CellListToSendByProc _buildCellListToSend() const; + + Array<int> _buildNbCellToSend(); + public: + + const std::shared_ptr<MeshType> dispatchedMesh() const + { + return m_dispatched_mesh; + } + + PASTIS_INLINE + const CellValue<const int>& cellNewOwner() const + { + return m_cell_new_owner; + } + + PASTIS_INLINE + const FaceValue<const int>& faceNewOwner() const + { + return m_face_new_owner; + } + + PASTIS_INLINE + const NodeValue<const int>& nodeNewOwner() const + { + return m_node_new_owner; + } + + template<typename DataType, typename ConnectivityPtr> + std::vector<Array<std::remove_const_t<DataType>>> + exchange(ItemValue<DataType, ItemType::cell, ConnectivityPtr> cell_value) const + { + using MutableDataType = std::remove_const_t<DataType>; + std::vector<Array<DataType>> cell_value_to_send_by_proc(parallel::size()); + for (size_t i=0; i<parallel::size(); ++i) { + const Array<const CellId>& cell_list = m_cell_list_to_send_by_proc[i]; + Array<MutableDataType> cell_value_list(cell_list.size()); + parallel_for (cell_list.size(), PASTIS_LAMBDA(const CellId& cell_id) { + cell_value_list[cell_id] = cell_value[cell_list[cell_id]]; + }); + cell_value_to_send_by_proc[i] = cell_value_list; + } + + std::vector<Array<MutableDataType>> recv_cell_value_by_proc(parallel::size()); + for (size_t i=0; i<parallel::size(); ++i) { + recv_cell_value_by_proc[i] = Array<MutableDataType>(m_nb_cell_to_recv_by_proc[i]); + } + + parallel::exchange(cell_value_to_send_by_proc, recv_cell_value_by_proc); + return recv_cell_value_by_proc; + } + + [[deprecated]] + const CellListToSendByProc& cell_list_to_send_by_proc() const + { + return m_cell_list_to_send_by_proc; + } + + MeshDispatcher(const MeshType& mesh); + + MeshDispatcher(const MeshDispatcher&) = delete; + ~MeshDispatcher() = default; +}; + + +#endif // MESH_DISPATCHER_HPP