From: Martin Kronbichler Date: Tue, 24 May 2022 11:18:42 +0000 (+0200) Subject: Avoid usage of CellDataTransferBuffer in GridTools::exchange_cell_data_to_ghosts X-Git-Tag: v9.4.0-rc1~152^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=76fcf5d26a00de6b80725c3a249207a806ad5458;p=dealii.git Avoid usage of CellDataTransferBuffer in GridTools::exchange_cell_data_to_ghosts --- diff --git a/include/deal.II/grid/grid_tools.h b/include/deal.II/grid/grid_tools.h index f3f85b5ee9..5d7f894d83 100644 --- a/include/deal.II/grid/grid_tools.h +++ b/include/deal.II/grid/grid_tools.h @@ -4407,7 +4407,7 @@ namespace GridTools const int mpi_tag_reply = Utilities::MPI::internal::Tags::exchange_cell_data_reply; - // send our requests: + // send our requests std::vector requests(ghost_owners.size()); { unsigned int idx = 0; @@ -4426,16 +4426,7 @@ namespace GridTools } } - using DestinationToBufferMap = - std::map>; - DestinationToBufferMap destination_to_data_buffer_map; - - // receive requests and reply with the ghost indices - std::vector> cell_data_to_send( - ghost_owners.size()); - std::vector> - send_dof_numbers_and_indices(ghost_owners.size()); + // receive requests and reply with the results std::vector reply_requests(ghost_owners.size()); std::vector> sendbuffers(ghost_owners.size()); @@ -4451,14 +4442,18 @@ namespace GridTools int len; ierr = MPI_Get_count(&status, MPI_BYTE, &len); AssertThrowMPI(ierr); - Assert(len % sizeof(cell_data_to_send[idx][0]) == 0, + Assert(len % sizeof(typename CellId::binary_type) == 0, ExcInternalError()); const unsigned int n_cells = len / sizeof(typename CellId::binary_type); - cell_data_to_send[idx].resize(n_cells); + std::vector cells_with_requests( + n_cells); + std::vector data_to_send; + data_to_send.reserve(n_cells); + std::vector cell_carries_data(n_cells, false); - ierr = MPI_Recv(cell_data_to_send[idx].data(), + ierr = MPI_Recv(cells_with_requests.data(), len, MPI_BYTE, status.MPI_SOURCE, @@ -4471,34 +4466,46 @@ namespace GridTools for (unsigned int c = 0; c < static_cast(n_cells); ++c) { const auto cell = - tria->create_cell_iterator(CellId(cell_data_to_send[idx][c])); + tria->create_cell_iterator(CellId(cells_with_requests[c])); - MeshCellIteratorType mesh_it(tria, + MeshCellIteratorType mesh_it(tria, cell->level(), cell->index(), &mesh); - const std_cxx17::optional data = pack(mesh_it); + const std_cxx17::optional data = pack(mesh_it); if (data) { - typename DestinationToBufferMap::iterator p = - destination_to_data_buffer_map - .insert(std::make_pair( - idx, - GridTools::CellDataTransferBuffer())) - .first; - - p->second.cell_ids.emplace_back(cell->id()); - p->second.data.emplace_back(*data); + data_to_send.emplace_back(*data); + cell_carries_data[c] = true; } } - // send reply - GridTools::CellDataTransferBuffer &data = - destination_to_data_buffer_map[idx]; - - sendbuffers[idx] = - Utilities::pack(data, /*enable_compression*/ false); + // collect data for sending the reply in a buffer + + // (a) make room for storing the local offsets in case we receive + // other data + sendbuffers[idx].resize(sizeof(std::size_t)); + + // (b) append the actual data and store how much memory it + // corresponds to, which we then insert into the leading position of + // the sendbuffer + std::size_t size_of_send = + Utilities::pack(data_to_send, + sendbuffers[idx], + /*enable_compression*/ false); + std::memcpy(sendbuffers[idx].data(), + &size_of_send, + sizeof(std::size_t)); + + // (c) append information of certain cells that got left out in case + // we need it + if (data_to_send.size() < n_cells) + Utilities::pack(cell_carries_data, + sendbuffers[idx], + /*enable_compression*/ false); + + // send data ierr = MPI_Isend(sendbuffers[idx].data(), sendbuffers[idx].size(), MPI_BYTE, @@ -4511,7 +4518,7 @@ namespace GridTools // finally receive the replies std::vector receive; - for (unsigned int idx = 0; idx < ghost_owners.size(); ++idx) + for (unsigned int id = 0; id < neighbor_cell_list.size(); ++id) { MPI_Status status; int ierr = MPI_Probe(MPI_ANY_SOURCE, @@ -4526,8 +4533,7 @@ namespace GridTools receive.resize(len); - char *ptr = receive.data(); - ierr = MPI_Recv(ptr, + ierr = MPI_Recv(receive.data(), len, MPI_BYTE, status.MPI_SOURCE, @@ -4536,23 +4542,52 @@ namespace GridTools &status); AssertThrowMPI(ierr); - auto cellinfo = - Utilities::unpack>( - receive, /*enable_compression*/ false); - - DataType *data = cellinfo.data.data(); - for (unsigned int c = 0; c < cellinfo.cell_ids.size(); ++c, ++data) + // (a) first determine the length of the data section in the + // received buffer + auto data_iterator = receive.begin(); + std::size_t size_of_received_data = + Utilities::unpack(data_iterator, + data_iterator + sizeof(std::size_t)); + data_iterator += sizeof(std::size_t); + + // (b) unpack the data section in the indicated region + auto received_data = Utilities::unpack>( + data_iterator, + data_iterator + size_of_received_data, + /*enable_compression*/ false); + data_iterator += size_of_received_data; + + // (c) check if the received data contained fewer entries than the + // number of cells we identified in the beginning, in which case we + // need to extract the boolean vector with the relevant information + const std::vector &this_cell_list = + neighbor_cell_list[status.MPI_SOURCE]; + AssertIndexRange(received_data.size(), this_cell_list.size() + 1); + std::vector cells_with_data; + if (received_data.size() < this_cell_list.size()) { - const typename Triangulation::cell_iterator - tria_cell = tria->create_cell_iterator(cellinfo.cell_ids[c]); + cells_with_data = Utilities::unpack>( + data_iterator, receive.end(), /*enable_compression*/ false); + AssertDimension(cells_with_data.size(), this_cell_list.size()); + } - MeshCellIteratorType cell(tria, - tria_cell->level(), - tria_cell->index(), - &mesh); + // (d) go through the received data and call the user-provided + // unpack function + auto received_data_iterator = received_data.begin(); + for (unsigned int c = 0; c < this_cell_list.size(); ++c) + if (cells_with_data.empty() || cells_with_data[c]) + { + const typename Triangulation::cell_iterator + tria_cell = tria->create_cell_iterator(this_cell_list[c]); - unpack(cell, *data); - } + MeshCellIteratorType cell(tria, + tria_cell->level(), + tria_cell->index(), + &mesh); + + unpack(cell, *received_data_iterator); + ++received_data_iterator; + } } // make sure that all communication is finished