From 90e37ef2033acdabc8170ade918b71cc92ae86ba Mon Sep 17 00:00:00 2001 From: Marc Fehling Date: Fri, 6 Jul 2018 12:15:07 -0600 Subject: [PATCH] Unpack functions take buffers as call-by-reference. Pack functions return the packed data. --- include/deal.II/base/quadrature_point_data.h | 22 ++++----- .../deal.II/distributed/solution_transfer.h | 7 ++- include/deal.II/distributed/tria.h | 38 +++++++-------- include/deal.II/particles/particle_handler.h | 7 ++- source/distributed/solution_transfer.cc | 15 +++--- source/distributed/tria.cc | 21 ++++----- source/particles/particle_handler.cc | 46 +++++++++---------- tests/mpi/attach_data_01.cc | 13 ++---- tests/mpi/repartition_02.cc | 13 ++---- 9 files changed, 83 insertions(+), 99 deletions(-) diff --git a/include/deal.II/base/quadrature_point_data.h b/include/deal.II/base/quadrature_point_data.h index 3e07cd955b..de84d8d957 100644 --- a/include/deal.II/base/quadrature_point_data.h +++ b/include/deal.II/base/quadrature_point_data.h @@ -409,13 +409,12 @@ namespace parallel * objects that can later be retrieved after refinement, coarsening and * repartitioning. */ - void + std::vector pack_function( const typename parallel::distributed::Triangulation::cell_iterator &cell, const typename parallel::distributed::Triangulation::CellStatus - status, - std::vector &data); + status); /** * A callback function used to unpack the data on the current mesh that @@ -429,7 +428,7 @@ namespace parallel const typename parallel::distributed::Triangulation::CellStatus status, const boost::iterator_range::const_iterator> - data_range); + &data_range); /** * FiniteElement used to project data from and to quadrature points. @@ -798,8 +797,7 @@ namespace parallel &ContinuousQuadratureDataTransfer::pack_function, this, std::placeholders::_1, - std::placeholders::_2, - std::placeholders::_3)); + std::placeholders::_2)); } @@ -825,13 +823,12 @@ namespace parallel template - void + std::vector ContinuousQuadratureDataTransfer::pack_function( const typename parallel::distributed::Triangulation::cell_iterator &cell, const typename parallel::distributed::Triangulation< - dim>::CellStatus /*status*/, - std::vector &data) + dim>::CellStatus /*status*/) { pack_cell_data(cell, data_storage, matrix_quadrature); @@ -840,7 +837,7 @@ namespace parallel // to get consistent data sizes on each cell for the fixed size transfer, // we won't allow compression - Utilities::pack(matrix_dofs, data, /*allow_compression=*/false); + return Utilities::pack(matrix_dofs, /*allow_compression=*/false); } @@ -851,8 +848,9 @@ namespace parallel const typename parallel::distributed::Triangulation::cell_iterator &cell, const typename parallel::distributed::Triangulation::CellStatus - status, - const boost::iterator_range::const_iterator> data_range) + status, + const boost::iterator_range::const_iterator> + &data_range) { Assert((status != parallel::distributed::Triangulation::CELL_COARSEN), diff --git a/include/deal.II/distributed/solution_transfer.h b/include/deal.II/distributed/solution_transfer.h index 72891ce1cb..756ec39b63 100644 --- a/include/deal.II/distributed/solution_transfer.h +++ b/include/deal.II/distributed/solution_transfer.h @@ -255,13 +255,12 @@ namespace parallel * objects that can later be retrieved after refinement, coarsening and * repartitioning. */ - void + std::vector pack_callback( const typename Triangulation:: cell_iterator &cell, const typename Triangulation:: - CellStatus status, - std::vector &data); + CellStatus status); /** * A callback function used to unpack the data on the current mesh that @@ -275,7 +274,7 @@ namespace parallel const typename Triangulation:: CellStatus status, const boost::iterator_range::const_iterator> - data_range, + & data_range, std::vector &all_out); diff --git a/include/deal.II/distributed/tria.h b/include/deal.II/distributed/tria.h index e3ffca184f..4c0cd6986b 100644 --- a/include/deal.II/distributed/tria.h +++ b/include/deal.II/distributed/tria.h @@ -645,7 +645,7 @@ namespace parallel * of classes that do this is parallel::distributed::SolutionTransfer * where each parallel::distributed::SolutionTransfer object that works * on the current Triangulation object then needs to register its intent. - * Each of these parties registers a call-back function (the second + * Each of these parties registers a callback function (the second * argument here, @p pack_callback) that will be called whenever the * triangulation's execute_coarsening_and_refinement() or save() * functions are called. @@ -715,8 +715,9 @@ namespace parallel * using save() and load(), then the cell status argument with which * the callback is called will always be `CELL_PERSIST`. * - * The third argument given to the callback is a reference to the - * buffer on which the packed data will be appended at the back. + * The callback function is expected to return a memory chunk of the + * format `std::vector`, representing the packed data on a + * certain cell. * * @note The purpose of this function is to register intent to * attach data for a single, subsequent call to @@ -729,9 +730,9 @@ namespace parallel */ unsigned int register_data_attach( - const std::function &)> &pack_callback); + const std::function(const cell_iterator &, + const CellStatus)> + &pack_callback); /** * This function is the opposite of register_data_attach(). It is called @@ -784,10 +785,10 @@ namespace parallel void notify_ready_to_unpack( const unsigned int handle, - const std::function< - void(const cell_iterator &, - const CellStatus, - const boost::iterator_range::const_iterator>)> + const std::function::const_iterator> &)> &unpack_callback); /** @@ -888,10 +889,9 @@ namespace parallel */ unsigned int n_attached_deserialize; - using pack_callback_t = std::function< - void(typename Triangulation::cell_iterator, - CellStatus, - std::vector &)>; + using pack_callback_t = std::function( + typename Triangulation::cell_iterator, + CellStatus)>; /** * These callback functions will be stored in the order on how they have @@ -1010,7 +1010,7 @@ namespace parallel const typename dealii::Triangulation::cell_iterator &, const typename dealii::Triangulation::CellStatus &, - const boost::iterator_range::const_iterator>)> + const boost::iterator_range::const_iterator> &)> &unpack_callback) const; /** @@ -1289,10 +1289,10 @@ namespace parallel */ unsigned int register_data_attach( - const std::function( const typename dealii::Triangulation<1, spacedim>::cell_iterator &, - const typename dealii::Triangulation<1, spacedim>::CellStatus, - std::vector &)> &pack_callback); + const typename dealii::Triangulation<1, spacedim>::CellStatus)> + &pack_callback); /** * This function is not implemented, but needs to be present for the @@ -1304,7 +1304,7 @@ namespace parallel const std::function::cell_iterator &, const typename dealii::Triangulation<1, spacedim>::CellStatus, - const boost::iterator_range::const_iterator>)> + const boost::iterator_range::const_iterator> &)> &unpack_callback); /** diff --git a/include/deal.II/particles/particle_handler.h b/include/deal.II/particles/particle_handler.h index 0b6f8963fc..ade19553b5 100644 --- a/include/deal.II/particles/particle_handler.h +++ b/include/deal.II/particles/particle_handler.h @@ -507,11 +507,10 @@ namespace Particles * before a refinement step. All particles have to be attached to their * cell to be sent around to the new processes. */ - void + std::vector store_particles( const typename Triangulation::cell_iterator &cell, - const typename Triangulation::CellStatus status, - std::vector & data) const; + const typename Triangulation::CellStatus status) const; /** * Called by listener functions after a refinement step. The local map @@ -522,7 +521,7 @@ namespace Particles const typename Triangulation::cell_iterator &cell, const typename Triangulation::CellStatus status, const boost::iterator_range::const_iterator> - data_range); + &data_range); }; /* ---------------------- inline and template functions ------------------ */ diff --git a/source/distributed/solution_transfer.cc b/source/distributed/solution_transfer.cc index add27f1cdd..349ffbae7a 100644 --- a/source/distributed/solution_transfer.cc +++ b/source/distributed/solution_transfer.cc @@ -89,8 +89,7 @@ namespace parallel &SolutionTransfer::pack_callback, this, std::placeholders::_1, - std::placeholders::_2, - std::placeholders::_3)); + std::placeholders::_2)); } @@ -202,13 +201,12 @@ namespace parallel template - void + std::vector SolutionTransfer::pack_callback( const typename Triangulation:: cell_iterator &cell_, const typename Triangulation:: - CellStatus /*status*/, - std::vector &data) + CellStatus /*status*/) { typename DoFHandlerType::cell_iterator cell(*cell_, dof_handler); @@ -228,7 +226,7 @@ namespace parallel // to get consistent data sizes on each cell for the fixed size transfer, // we won't allow compression - Utilities::pack(dofvalues, data, /*allow_compression=*/false); + return Utilities::pack(dofvalues, /*allow_compression=*/false); } @@ -240,8 +238,9 @@ namespace parallel cell_iterator &cell_, const typename Triangulation:: CellStatus /*status*/, - const boost::iterator_range::const_iterator> data_range, - std::vector & all_out) + const boost::iterator_range::const_iterator> + & data_range, + std::vector &all_out) { typename DoFHandlerType::cell_iterator cell(*cell_, dof_handler); diff --git a/source/distributed/tria.cc b/source/distributed/tria.cc index c4b959893b..7fa7c55198 100644 --- a/source/distributed/tria.cc +++ b/source/distributed/tria.cc @@ -1243,7 +1243,7 @@ namespace parallel callback_it != pack_callbacks.cend(); ++callback_it, ++data_it) { - (*callback_it)(dealii_cell, cell_status, *data_it); + *data_it = (*callback_it)(dealii_cell, cell_status); } } } // loop over quad_cell_relations @@ -1420,7 +1420,7 @@ namespace parallel const std::function::cell_iterator &, const typename dealii::Triangulation::CellStatus &, - const boost::iterator_range::const_iterator>)> + const boost::iterator_range::const_iterator> &)> &unpack_callback) const { Assert(cumulative_sizes_fixed.size() > 0, @@ -3848,9 +3848,8 @@ namespace parallel template unsigned int Triangulation::register_data_attach( - const std::function &)> &pack_callback) + const std::function(const cell_iterator &, + const CellStatus)> &pack_callback) { // add new callback_set cell_attached_data.pack_callbacks.push_back(pack_callback); @@ -3869,7 +3868,7 @@ namespace parallel const std::function< void(const cell_iterator &, const CellStatus, - const boost::iterator_range::const_iterator>)> + const boost::iterator_range::const_iterator> &)> &unpack_callback) { # ifdef DEBUG @@ -4458,10 +4457,10 @@ namespace parallel template unsigned int Triangulation<1, spacedim>::register_data_attach( - const std::function< - void(const typename dealii::Triangulation<1, spacedim>::cell_iterator &, - const typename dealii::Triangulation<1, spacedim>::CellStatus, - std::vector &)> & /*pack_callback*/) + const std::function( + const typename dealii::Triangulation<1, spacedim>::cell_iterator &, + const typename dealii::Triangulation<1, spacedim>::CellStatus)> + & /*pack_callback*/) { Assert(false, ExcNotImplemented()); return 0; @@ -4476,7 +4475,7 @@ namespace parallel const std::function< void(const typename dealii::Triangulation<1, spacedim>::cell_iterator &, const typename dealii::Triangulation<1, spacedim>::CellStatus, - const boost::iterator_range::const_iterator>)> + const boost::iterator_range::const_iterator> &)> & /*unpack_callback*/) { Assert(false, ExcNotImplemented()); diff --git a/source/particles/particle_handler.cc b/source/particles/particle_handler.cc index 0941499e2f..2b4131ef04 100644 --- a/source/particles/particle_handler.cc +++ b/source/particles/particle_handler.cc @@ -1037,16 +1037,14 @@ namespace Particles if (global_max_particles_per_cell > 0) { - const std::function< - void(const typename Triangulation::cell_iterator &, - const typename Triangulation::CellStatus, - std::vector &)> + const std::function( + const typename Triangulation::cell_iterator &, + const typename Triangulation::CellStatus)> callback_function = std::bind(&ParticleHandler::store_particles, std::cref(*this), std::placeholders::_1, - std::placeholders::_2, - std::placeholders::_3); + std::placeholders::_2); handle = non_const_triangulation->register_data_attach(callback_function); @@ -1075,16 +1073,14 @@ namespace Particles // data correctly. Only do this if something was actually stored. if (serialization && (global_max_particles_per_cell > 0)) { - const std::function< - void(const typename Triangulation::cell_iterator &, - const typename Triangulation::CellStatus, - std::vector &)> + const std::function( + const typename Triangulation::cell_iterator &, + const typename Triangulation::CellStatus)> callback_function = std::bind(&ParticleHandler::store_particles, std::cref(*this), std::placeholders::_1, - std::placeholders::_2, - std::placeholders::_3); + std::placeholders::_2); handle = non_const_triangulation->register_data_attach(callback_function); @@ -1093,10 +1089,10 @@ namespace Particles // Check if something was stored and load it if (handle != numbers::invalid_unsigned_int) { - const std::function< - void(const typename Triangulation::cell_iterator &, - const typename Triangulation::CellStatus, - const boost::iterator_range::const_iterator>)> + const std::function::cell_iterator &, + const typename Triangulation::CellStatus, + const boost::iterator_range::const_iterator> &)> callback_function = std::bind(&ParticleHandler::load_particles, std::ref(*this), @@ -1117,11 +1113,10 @@ namespace Particles template - void + std::vector ParticleHandler::store_particles( const typename Triangulation::cell_iterator &cell, - const typename Triangulation::CellStatus status, - std::vector &data_vector) const + const typename Triangulation::CellStatus status) const { // ------------------- // HOTFIX: resize memory and static_cast std::vector to void* @@ -1149,9 +1144,8 @@ namespace Particles (size_per_particle * global_max_particles_per_cell) * (serialization ? 1 : std::pow(2, dim)); - const size_t previous_size = data_vector.size(); - data_vector.resize(previous_size + transfer_size_per_cell); - void *data = static_cast(data_vector.data() + previous_size); + std::vector data_vector(transfer_size_per_cell); + void * data = static_cast(data_vector.data()); // -------------------- unsigned int n_particles(0); @@ -1217,14 +1211,16 @@ namespace Particles } else Assert(false, ExcInternalError()); + + return data_vector; } template void ParticleHandler::load_particles( - const typename Triangulation::cell_iterator & cell, - const typename Triangulation::CellStatus status, - const boost::iterator_range::const_iterator> data_range) + const typename Triangulation::cell_iterator & cell, + const typename Triangulation::CellStatus status, + const boost::iterator_range::const_iterator> &data_range) { // ------------------- // HOTFIX: static_cast std::vector to void* diff --git a/tests/mpi/attach_data_01.cc b/tests/mpi/attach_data_01.cc index 06fdafdca4..9c556459bd 100644 --- a/tests/mpi/attach_data_01.cc +++ b/tests/mpi/attach_data_01.cc @@ -33,13 +33,12 @@ template -void +std::vector pack_function( const typename parallel::distributed::Triangulation::cell_iterator &cell, const typename parallel::distributed::Triangulation::CellStatus - status, - std::vector &data) + status) { static int some_number = 0; deallog << "packing cell " << cell->id() << " with data=" << some_number @@ -63,9 +62,7 @@ pack_function( Assert(!cell->has_children(), ExcInternalError()); } - Utilities::pack(some_number, data, /*allow_compression=*/false); - - ++some_number; + return Utilities::pack(some_number++, /*allow_compression=*/false); } template @@ -74,8 +71,8 @@ unpack_function( const typename parallel::distributed::Triangulation::cell_iterator &cell, const typename parallel::distributed::Triangulation::CellStatus - status, - const boost::iterator_range::const_iterator> data_range) + status, + const boost::iterator_range::const_iterator> &data_range) { const int intdata = Utilities::unpack(data_range.begin(), data_range.end(), diff --git a/tests/mpi/repartition_02.cc b/tests/mpi/repartition_02.cc index e01a3fd284..ce953cc88e 100644 --- a/tests/mpi/repartition_02.cc +++ b/tests/mpi/repartition_02.cc @@ -33,13 +33,12 @@ template -void +std::vector pack_function( const typename parallel::distributed::Triangulation::cell_iterator &cell, const typename parallel::distributed::Triangulation::CellStatus - status, - std::vector &data) + status) { static int some_number = cell->index(); deallog << "packing cell " << cell->id() << " with data=" << some_number @@ -63,9 +62,7 @@ pack_function( Assert(!cell->has_children(), ExcInternalError()); } - Utilities::pack(some_number, data, /*allow_compression=*/false); - - ++some_number; + return Utilities::pack(some_number++, /*allow_compression=*/false); } template @@ -74,8 +71,8 @@ unpack_function( const typename parallel::distributed::Triangulation::cell_iterator &cell, const typename parallel::distributed::Triangulation::CellStatus - status, - const boost::iterator_range::const_iterator> data_range) + status, + const boost::iterator_range::const_iterator> &data_range) { const int number = Utilities::unpack(data_range.begin(), data_range.end(), -- 2.39.5