From 7e6df09e725699dca66c269d64a7e03bec095a90 Mon Sep 17 00:00:00 2001 From: Denis Davydov Date: Fri, 3 Nov 2017 17:21:00 +0100 Subject: [PATCH] move CellDataTransferBuffer out of internal namespace --- include/deal.II/distributed/grid_tools.h | 236 +++++++++++++---------- 1 file changed, 129 insertions(+), 107 deletions(-) diff --git a/include/deal.II/distributed/grid_tools.h b/include/deal.II/distributed/grid_tools.h index 822d34eb8d..8a38dc543c 100644 --- a/include/deal.II/distributed/grid_tools.h +++ b/include/deal.II/distributed/grid_tools.h @@ -134,6 +134,53 @@ namespace parallel exchange_cell_data_to_ghosts (const MeshType &mesh, const std::function (const typename MeshType::active_cell_iterator &)> &pack, const std::function &unpack); + + /** + * A structure that allows the transfer of data of type T from one processor + * to another. It corresponds to a packed buffer that stores a list of + * cells and an array of type T. + * + * The vector @p data is the same size as @p cell_ids. + */ + template + struct CellDataTransferBuffer + { + std::vector cell_ids; + std::vector data; + + /** + * Write the data of this object to a stream for the purpose of + * serialization. + */ + template + void save (Archive &ar, + const unsigned int /*version*/) const; + + /** + * Read the data of this object from a stream for the purpose of + * serialization. Throw away the previous content. + */ + template + void load (Archive &ar, + const unsigned int version); + + BOOST_SERIALIZATION_SPLIT_MEMBER() + + /** + * Pack the data that corresponds to this object into a buffer in + * the form of a vector of chars and return it. + */ + std::vector pack_data () const; + + /** + * Given a buffer in the form of an array of chars, unpack it and + * restore the current object to the state that it was when + * it was packed into said buffer by the pack_data() function. + */ + void unpack_data (const std::vector &buffer); + + }; + } @@ -145,133 +192,108 @@ namespace parallel { namespace GridTools { - namespace internal + + template + template + void + CellDataTransferBuffer::save (Archive &ar, + const unsigned int /*version*/) const { - /** - * A structure that allows the transfer of data of type T from one processor - * to another. It corresponds to a packed buffer that stores a list of - * cells and an array of type T. - * - * The vector @p data is the same size as @p cell_ids. - */ - template - struct CellDataTransferBuffer - { - std::vector cell_ids; - std::vector data; - - /** - * Write the data of this object to a stream for the purpose of - * serialization. - */ - template - void save (Archive &ar, - const unsigned int /*version*/) const + Assert(cell_ids.size() == data.size(), ExcInternalError()); + // archive the cellids in an efficient binary format + const size_t n_cells = cell_ids.size(); + ar &n_cells; + for (auto &it : cell_ids) { - Assert(cell_ids.size() == data.size(), ExcInternalError()); - // archive the cellids in an efficient binary format - const size_t n_cells = cell_ids.size(); - ar &n_cells; - for (auto &it : cell_ids) - { - CellId::binary_type binary_cell_id = it.template to_binary(); - ar &binary_cell_id; - } - - ar &data; + CellId::binary_type binary_cell_id = it.template to_binary(); + ar &binary_cell_id; } - /** - * Read the data of this object from a stream for the purpose of - * serialization. Throw away the previous content. - */ - template - void load (Archive &ar, - const unsigned int /*version*/) + ar &data; + } + + + + template + template + void + CellDataTransferBuffer::load (Archive &ar, + const unsigned int /*version*/) + { + size_t n_cells; + ar &n_cells; + cell_ids.clear(); + cell_ids.reserve(n_cells); + for (unsigned int c=0; c - pack_data () const - { - // set up a buffer and then use it as the target of a compressing - // stream into which we serialize the current object - std::vector buffer; - { + template + std::vector + CellDataTransferBuffer::pack_data () const + { + // set up a buffer and then use it as the target of a compressing + // stream into which we serialize the current object + std::vector buffer; + { #ifdef DEAL_II_WITH_ZLIB - boost::iostreams::filtering_ostream out; - out.push(boost::iostreams::gzip_compressor - (boost::iostreams::gzip_params - (boost::iostreams::gzip::best_compression))); - out.push(boost::iostreams::back_inserter(buffer)); - - boost::archive::binary_oarchive archive(out); - archive << *this; - out.flush(); + boost::iostreams::filtering_ostream out; + out.push(boost::iostreams::gzip_compressor + (boost::iostreams::gzip_params + (boost::iostreams::gzip::best_compression))); + out.push(boost::iostreams::back_inserter(buffer)); + + boost::archive::binary_oarchive archive(out); + archive << *this; + out.flush(); #else - std::ostringstream out; - boost::archive::binary_oarchive archive(out); - archive << *this; - const std::string &s = out.str(); - buffer.reserve(s.size()); - buffer.assign(s.begin(), s.end()); + std::ostringstream out; + boost::archive::binary_oarchive archive(out); + archive << *this; + const std::string &s = out.str(); + buffer.reserve(s.size()); + buffer.assign(s.begin(), s.end()); #endif - } + } - return buffer; - } + return buffer; + } - /** - * Given a buffer in the form of an array of chars, unpack it and - * restore the current object to the state that it was when - * it was packed into said buffer by the pack_data() function. - */ - void unpack_data (const std::vector &buffer) - { - std::string decompressed_buffer; + template + void + CellDataTransferBuffer::unpack_data (const std::vector &buffer) + { + std::string decompressed_buffer; - // first decompress the buffer - { + // first decompress the buffer + { #ifdef DEAL_II_WITH_ZLIB - boost::iostreams::filtering_ostream decompressing_stream; - decompressing_stream.push(boost::iostreams::gzip_decompressor()); - decompressing_stream.push(boost::iostreams::back_inserter(decompressed_buffer)); - decompressing_stream.write (buffer.data(), buffer.size()); + boost::iostreams::filtering_ostream decompressing_stream; + decompressing_stream.push(boost::iostreams::gzip_decompressor()); + decompressing_stream.push(boost::iostreams::back_inserter(decompressed_buffer)); + decompressing_stream.write (buffer.data(), buffer.size()); #else - decompressed_buffer.assign (buffer.begin(), buffer.end()); + decompressed_buffer.assign (buffer.begin(), buffer.end()); #endif - } + } - // then restore the object from the buffer - std::istringstream in(decompressed_buffer); - boost::archive::binary_iarchive archive(in); - - archive >> *this; - } - }; + // then restore the object from the buffer + std::istringstream in(decompressed_buffer); + boost::archive::binary_iarchive archive(in); + archive >> *this; } + template void exchange_cell_data_to_ghosts (const MeshType &mesh, @@ -292,7 +314,7 @@ namespace parallel ExcMessage("The function exchange_cell_data_to_ghosts() only works with parallel triangulations.")); // map neighbor_id -> data_buffer where we accumulate the data to send - typedef std::map > + typedef std::map > DestinationToBufferMap; DestinationToBufferMap destination_to_data_buffer_map; @@ -340,7 +362,7 @@ namespace parallel // or create an empty one otherwise typename DestinationToBufferMap::iterator p = destination_to_data_buffer_map.insert (std::make_pair(subdomain, - internal::CellDataTransferBuffer())) + CellDataTransferBuffer())) .first; p->second.cell_ids.emplace_back(cellid); @@ -362,7 +384,7 @@ namespace parallel it!=ghost_owners.end(); ++it, ++idx) { - internal::CellDataTransferBuffer &data = destination_to_data_buffer_map[*it]; + CellDataTransferBuffer &data = destination_to_data_buffer_map[*it]; // pack all the data into the buffer for this recipient and send it. // keep data around till we can make sure that the packet has been @@ -392,7 +414,7 @@ namespace parallel tria->get_communicator(), &status); AssertThrowMPI(ierr); - internal::CellDataTransferBuffer cellinfo; + CellDataTransferBuffer cellinfo; cellinfo.unpack_data(receive); DataType *data = cellinfo.data.data(); -- 2.39.5