]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Simplify particle transfer during refinement 12450/head
authorRene Gassmoeller <rene.gassmoeller@mailbox.org>
Fri, 11 Jun 2021 19:52:02 +0000 (15:52 -0400)
committerRene Gassmoeller <rene.gassmoeller@mailbox.org>
Tue, 22 Jun 2021 18:24:55 +0000 (14:24 -0400)
doc/news/changes/major/20210606Gassmoeller
examples/step-68/step-68.cc
examples/step-70/step-70.cc
include/deal.II/particles/particle_handler.h
source/particles/particle_handler.cc
tests/particles/particle_handler_05.cc
tests/particles/particle_handler_08.cc
tests/particles/particle_handler_14.cc
tests/particles/step-68.cc
tests/serialization/particle_handler_01.cc
tests/serialization/particle_handler_02.cc

index 882398fb094dadcc6edc82d640c99fccc1279357..1f9a96f5a1135b5e300266b38ad8bb227ea7bf5f 100644 (file)
@@ -2,4 +2,4 @@ Improved: The internal particle handler storage structure has been reworked to
 allow for significantly faster particle sorting and iteration as well as for
 future improvements and extensions.
 <br>
-(Rene Gassmoeller, 2021/06/06)
+(Rene Gassmoeller, Bruno Blais, Martin Kronbichler, Peter Munch, 2021/06/06)
index 5c89f7b854f37752180e97efa92b2c153aaccc0f..2b9efafe7036ce2f246cbb4be768b11567c32298 100644 (file)
@@ -392,17 +392,13 @@ namespace Step68
     background_triangulation.refine_global(par.fluid_refinement);
 
     // In order to consider the particles when repartitioning the triangulation
-    // the algorithm needs to know three things:
+    // the algorithm needs to know how much weight to assign to each cell
+    // (how many particles are in there).
     //
-    // 1. How much weight to assign to each cell (how many particles are in
-    // there);
-    // 2. How to pack the particles before shipping data around;
-    // 3. How to unpack the particles after repartitioning.
-    //
-    // We attach the correct functions to the signals inside
-    // parallel::distributed::Triangulation. These signal will be called every
-    // time the repartition() function is called. These connections only need to
-    // be created once, so we might as well have set them up in the constructor
+    // We attach a weight function to the signal inside
+    // parallel::distributed::Triangulation. This signal will be called every
+    // time the repartition() function is called. This connection only needs to
+    // be created once, so we might as well have set it up in the constructor
     // of this class, but for the purpose of this example we want to group the
     // particle related instructions.
     background_triangulation.signals.cell_weight.connect(
@@ -412,12 +408,6 @@ namespace Step68
         const typename parallel::distributed::Triangulation<dim>::CellStatus
           status) -> unsigned int { return this->cell_weight(cell, status); });
 
-    background_triangulation.signals.pre_distributed_repartition.connect(
-      [this]() { this->particle_handler.register_store_callback_function(); });
-
-    background_triangulation.signals.post_distributed_repartition.connect(
-      [&]() { this->particle_handler.register_load_callback_function(false); });
-
     // This initializes the background triangulation where the particles are
     // living and the number of properties of the particles.
     particle_handler.initialize(background_triangulation, mapping, 1 + dim);
@@ -704,7 +694,10 @@ namespace Step68
 
     pcout << "Repartitioning triangulation after particle generation"
           << std::endl;
+
+    particle_handler.prepare_for_coarsening_and_refinement();
     background_triangulation.repartition();
+    particle_handler.unpack_after_coarsening_and_refinement();
 
     // We set the initial property of the particles by doing an
     // explicit Euler iteration with a time-step of 0 both in the case
@@ -730,7 +723,10 @@ namespace Step68
 
         if ((discrete_time.get_step_number() % par.repartition_frequency) == 0)
           {
+            particle_handler.prepare_for_coarsening_and_refinement();
             background_triangulation.repartition();
+            particle_handler.unpack_after_coarsening_and_refinement();
+
             if (interpolated_velocity)
               setup_background_dofs();
           }
index 2ac8003840e5d7f20dc2299af9c75a92f3af818b..a91f6492692f116d228b1548ff1a9cd0d85d3da3 100644 (file)
@@ -1082,25 +1082,6 @@ namespace Step70
     // would be coupled to what is occurring in the fluid domain.
     locally_relevant_tracer_particle_coordinates =
       locally_owned_tracer_particle_coordinates;
-
-    // Finally, we make sure that upon refinement, particles are correctly
-    // transferred. When performing local refinement or coarsening, particles
-    // will land in another cell. We could in principle redistribute all
-    // particles after refining, however this would be overly expensive.
-    //
-    // The Particles::ParticleHandler class has a way to transfer information
-    // from a cell to its children or to its parent upon refinement, without the
-    // need to reconstruct the entire data structure. This is done by
-    // registering two callback functions to the triangulation. These
-    // functions will receive a signal when refinement is about to happen, and
-    // when it has just happened, and will take care of transferring all
-    // information to the newly refined grid with minimal computational cost.
-    fluid_tria.signals.pre_distributed_refinement.connect(
-      [&]() { tracer_particle_handler.register_store_callback_function(); });
-
-    fluid_tria.signals.post_distributed_refinement.connect([&]() {
-      tracer_particle_handler.register_load_callback_function(false);
-    });
   }
 
 
@@ -1185,15 +1166,6 @@ namespace Step70
                                                    global_fluid_bounding_boxes,
                                                    properties);
 
-
-    // As in the previous function, we end by making sure that upon refinement,
-    // particles are correctly transferred:
-    fluid_tria.signals.pre_distributed_refinement.connect(
-      [&]() { solid_particle_handler.register_store_callback_function(); });
-
-    fluid_tria.signals.post_distributed_refinement.connect(
-      [&]() { solid_particle_handler.register_load_callback_function(false); });
-
     pcout << "Solid particles: " << solid_particle_handler.n_global_particles()
           << std::endl;
   }
@@ -1634,7 +1606,18 @@ namespace Step70
 
   // @sect4{Mesh refinement}
 
-  // We deal with mesh refinement in a completely standard way:
+  // We deal with mesh refinement in a completely standard way, except
+  // we now also transfer the particles of the two particle handlers from the
+  // existing to the refined mesh. When performing local refinement or
+  // coarsening, particles will land in another cell. We could in principle
+  // redistribute all particles after refining, however this would be overly
+  // expensive.
+  //
+  // The Particles::ParticleHandler class has a way to transfer information
+  // from a cell to its children or to its parent upon refinement, without the
+  // need to reconstruct the entire data structure. This is done similarly
+  // to the SolutionTransfer class by calling two functions, one to prepare
+  // for refinement, and one to transfer the information after refinement.
   template <int dim, int spacedim>
   void StokesImmersedProblem<dim, spacedim>::refine_and_transfer()
   {
@@ -1680,13 +1663,20 @@ namespace Step70
 
     parallel::distributed::SolutionTransfer<spacedim, LA::MPI::BlockVector>
       transfer(fluid_dh);
+
     fluid_tria.prepare_coarsening_and_refinement();
     transfer.prepare_for_coarsening_and_refinement(locally_relevant_solution);
+    tracer_particle_handler.prepare_for_coarsening_and_refinement();
+    solid_particle_handler.prepare_for_coarsening_and_refinement();
+
     fluid_tria.execute_coarsening_and_refinement();
 
     setup_dofs();
 
     transfer.interpolate(solution);
+    tracer_particle_handler.unpack_after_coarsening_and_refinement();
+    solid_particle_handler.unpack_after_coarsening_and_refinement();
+
     constraints.distribute(solution);
     locally_relevant_solution = solution;
   }
index a674d73d546d37d9a8843e4d667f8fb871025c0f..70ea6e864ec8d7ee50369c982c9f54788f613ee4 100644 (file)
@@ -702,21 +702,89 @@ namespace Particles
     void
     update_ghost_particles();
 
+    /**
+     * This function prepares the particle handler for a coarsening and
+     * refinement cycle, by storing the necessary information to transfer
+     * particles to their new cells. The implementation depends on the
+     * triangulation type that is connected to the particle handler and differs
+     * between shared and distributed triangulations. This function should be
+     * used like the corresponding function with the same name in the
+     * SolutionTransfer() class.
+     */
+    void
+    prepare_for_coarsening_and_refinement();
+
+    /**
+     * This function unpacks the particle data after a coarsening and
+     * refinement cycle, by reading the necessary information to transfer
+     * particles to their new cells. This function should be
+     * used like the SolutionTransfer::interpolate() function after mesh
+     * refinement has finished. Note that this function requires a working
+     * mapping, i.e. if you use a mapping class that requires setup after
+     * a mesh refinement (e.g. MappingQCache(), or MappingEulerian()), the
+     * mapping has to be ready before you can call this function.
+     *
+     * @note It is important to note, that if you do not call this function
+     * after a refinement/coarsening operation (and a repartitioning operation
+     * in a distributed triangulation) and the particle handler contained
+     * particles before, the particle handler will not be usable any more. Not
+     * calling this function after refinement would be similar to trying to
+     * access a solution vector after mesh refinement without first using a
+     * SolutionTransfer class to transfer the vector to the new mesh.
+     */
+    void
+    unpack_after_coarsening_and_refinement();
+
+    /**
+     * This function prepares the particle handler for serialization. This needs
+     * to be done before calling Triangulation::save(). This function should be
+     * used like the corresponding function with the same name in the
+     * SolutionTransfer() class.
+     *
+     * @note It is important to note, that if you do not call this function
+     * before a serialization operation in a distributed triangulation, the
+     * particles will not be serialized, and therefore will disappear after
+     * deserialization.
+     */
+    void
+    prepare_for_serialization();
+
+    /**
+     * Execute the deserialization of the particle data. This needs to be
+     * done after calling Triangulation::load(). The data must have been stored
+     * before the serialization of the triangulation using the
+     * prepare_for_serialization() function. This function should be used like
+     * the corresponding function with the same name in the SolutionTransfer
+     * class.
+     */
+    void
+    deserialize();
+
     /**
      * Callback function that should be called before every refinement
      * and when writing checkpoints. This function is used to
-     * register store_particles() with the triangulation. This function
+     * register pack_callback() with the triangulation. This function
      * is used in step-70.
+     *
+     * @deprecated Please use prepare_for_coarsening_and_refinement() or
+     * prepare_for_serialization() instead. See there for further information
+     * about the purpose of this function.
      */
+    DEAL_II_DEPRECATED_EARLY
     void
     register_store_callback_function();
 
     /**
      * Callback function that should be called after every refinement
      * and after resuming from a checkpoint.  This function is used to
-     * register load_particles() with the triangulation. This function
+     * register unpack_callback() with the triangulation. This function
      * is used in step-70.
+     *
+     * @deprecated Please use unpack_after_coarsening_and_refinement() or
+     * deserialize() instead. See there for further information about the
+     * purpose of this function.
      */
+    DEAL_II_DEPRECATED_EARLY
     void
     register_load_callback_function(const bool serialization);
 
@@ -890,8 +958,8 @@ namespace Particles
       load_callback;
 
     /**
-     * This variable is set by the register_store_callback_function()
-     * function and used by the register_load_callback_function() function
+     * This variable is set by the register_data_attach()
+     * function and used by the notify_ready_to_unpack() function
      * to check where the particle data was registered in the corresponding
      * triangulation object.
      */
@@ -984,22 +1052,64 @@ namespace Particles
      */
     internal::GhostParticlePartitioner<dim, spacedim> ghost_particles_cache;
 
+    /**
+     * Connect the particle handler to the relevant triangulation signals to
+     * appropriately react to changes in the underlying triangulation.
+     */
+    void
+    connect_to_triangulation_signals();
+
+    /**
+     * A list of connections with which this object connects to the
+     * triangulation to get information about when the triangulation changes.
+     */
+    std::vector<boost::signals2::connection> tria_listeners;
+
+    /**
+     * Function that gets called by the triangulation signals if
+     * the structure of the mesh has changed. This function is
+     * responsible for resizing the particle container if no
+     * particles are stored, i.e. if the usual call to
+     * prepare_for_..., and transfer_particles_after_... functions
+     * is not necessary.
+     */
+    void
+    post_mesh_change_action();
+
+    /**
+     * Function that should be called before every refinement
+     * and when writing checkpoints. This function is used to
+     * register pack_callback() with the triangulation.
+     */
+    void
+    register_data_attach();
+
+    /**
+     * Callback function that should be called after every refinement
+     * and after resuming from a checkpoint.  This function is used to
+     * call the notify_ready_to_unpack() function of the triangulation
+     * and hand over the unpack_callback() function, which will
+     * unpack the particle data for each cell.
+     */
+    void
+    notify_ready_to_unpack(const bool serialization);
+
     /**
      * Called by listener functions from Triangulation for every cell
      * before a refinement step. All particles have to be attached to their
-     * cell to be sent around to the new processes.
+     * cell to be sent around to the new cell and owning process.
      */
     std::vector<char>
-    store_particles(
+    pack_callback(
       const typename Triangulation<dim, spacedim>::cell_iterator &cell,
       const typename Triangulation<dim, spacedim>::CellStatus     status) const;
 
     /**
-     * Called by listener functions after a refinement step to receive
-     * particles and insert them into the particle container.
+     * Called by listener functions after a refinement step for each cell
+     * to unpack the particle data and transfer it to the local container.
      */
     void
-    load_particles(
+    unpack_callback(
       const typename Triangulation<dim, spacedim>::cell_iterator &cell,
       const typename Triangulation<dim, spacedim>::CellStatus     status,
       const boost::iterator_range<std::vector<char>::const_iterator>
index 2269b73c74bba32315b5b7f2d3e173a96ef986f1..a5de7223df613e71e1f1e1bcd9c6fcd04f82460c 100644 (file)
@@ -65,6 +65,7 @@ namespace Particles
     , store_callback()
     , load_callback()
     , handle(numbers::invalid_unsigned_int)
+    , tria_listeners()
   {}
 
 
@@ -86,17 +87,12 @@ namespace Particles
     , store_callback()
     , load_callback()
     , handle(numbers::invalid_unsigned_int)
+    , triangulation_cache(
+        std::make_unique<GridTools::Cache<dim, spacedim>>(triangulation,
+                                                          mapping))
+    , tria_listeners()
   {
-    triangulation_cache =
-      std::make_unique<GridTools::Cache<dim, spacedim>>(triangulation, mapping);
-
-    triangulation.signals.create.connect([&]() { this->clear_particles(); });
-    triangulation.signals.post_refinement.connect(
-      [&]() { this->clear_particles(); });
-    triangulation.signals.post_distributed_repartition.connect(
-      [&]() { this->clear_particles(); });
-    triangulation.signals.post_distributed_load.connect(
-      [&]() { this->clear_particles(); });
+    connect_to_triangulation_signals();
   }
 
 
@@ -105,6 +101,9 @@ namespace Particles
   ParticleHandler<dim, spacedim>::~ParticleHandler()
   {
     clear_particles();
+
+    for (const auto &connection : tria_listeners)
+      connection.disconnect();
   }
 
 
@@ -121,19 +120,6 @@ namespace Particles
     triangulation = &new_triangulation;
     mapping       = &new_mapping;
 
-    // Note that we connect all signals 'at_front' in case user code
-    // connected the register_store_callback_function and
-    // register_load_callback_function functions already to the new
-    // triangulation. We want to clear before loading.
-    triangulation->signals.create.connect([&]() { this->clear_particles(); },
-                                          boost::signals2::at_front);
-    triangulation->signals.post_distributed_refinement.connect(
-      [&]() { this->clear_particles(); }, boost::signals2::at_front);
-    triangulation->signals.post_distributed_repartition.connect(
-      [&]() { this->clear_particles(); }, boost::signals2::at_front);
-    triangulation->signals.post_distributed_load.connect(
-      [&]() { this->clear_particles(); }, boost::signals2::at_front);
-
     // Create the memory pool that will store all particle properties
     property_pool = std::make_unique<PropertyPool<dim, spacedim>>(n_properties);
 
@@ -144,6 +130,8 @@ namespace Particles
                                                         new_mapping);
 
     particles.resize(triangulation->n_active_cells());
+
+    connect_to_triangulation_signals();
   }
 
 
@@ -1862,76 +1850,194 @@ namespace Particles
   }
 
 
+  template <int dim, int spacedim>
+  void
+  ParticleHandler<dim, spacedim>::connect_to_triangulation_signals()
+  {
+    // First disconnect existing connections
+    for (const auto &connection : tria_listeners)
+      connection.disconnect();
+
+    tria_listeners.clear();
+
+    tria_listeners.push_back(triangulation->signals.create.connect([&]() {
+      this->initialize(*(this->triangulation),
+                       *(this->mapping),
+                       this->property_pool->n_properties_per_slot());
+    }));
+
+    this->tria_listeners.push_back(
+      this->triangulation->signals.clear.connect([&]() { this->clear(); }));
+
+    // for distributed triangulations, connect to distributed signals
+    if (dynamic_cast<const parallel::DistributedTriangulationBase<dim, spacedim>
+                       *>(&(*triangulation)) != nullptr)
+      {
+        tria_listeners.push_back(
+          triangulation->signals.post_distributed_refinement.connect(
+            [&]() { this->post_mesh_change_action(); }));
+        tria_listeners.push_back(
+          triangulation->signals.post_distributed_repartition.connect(
+            [&]() { this->post_mesh_change_action(); }));
+        tria_listeners.push_back(
+          triangulation->signals.post_distributed_load.connect(
+            [&]() { this->post_mesh_change_action(); }));
+      }
+    else
+      {
+        tria_listeners.push_back(triangulation->signals.post_refinement.connect(
+          [&]() { this->post_mesh_change_action(); }));
+      }
+  }
+
+
+
+  template <int dim, int spacedim>
+  void
+  ParticleHandler<dim, spacedim>::post_mesh_change_action()
+  {
+    Assert(triangulation != nullptr, ExcInternalError());
+
+    const bool distributed_triangulation =
+      dynamic_cast<
+        const parallel::DistributedTriangulationBase<dim, spacedim> *>(
+        &(*triangulation)) != nullptr;
+    (void)distributed_triangulation;
+
+    Assert(
+      distributed_triangulation || local_number_of_particles == 0,
+      ExcMessage(
+        "Mesh refinement in a non-distributed triangulation is not supported "
+        "by the ParticleHandler class. Either insert particles after mesh "
+        "creation, or use a distributed triangulation."));
+
+    // Resize the container if it is possible without
+    // transferring particles
+    if (local_number_of_particles == 0)
+      particles.resize(triangulation->n_active_cells());
+  }
+
+
+
+  template <int dim, int spacedim>
+  void
+  ParticleHandler<dim, spacedim>::prepare_for_coarsening_and_refinement()
+  {
+    register_data_attach();
+  }
+
+
+
+  template <int dim, int spacedim>
+  void
+  ParticleHandler<dim, spacedim>::prepare_for_serialization()
+  {
+    register_data_attach();
+  }
+
+
 
   template <int dim, int spacedim>
   void
   ParticleHandler<dim, spacedim>::register_store_callback_function()
+  {
+    register_data_attach();
+  }
+
+
+
+  template <int dim, int spacedim>
+  void
+  ParticleHandler<dim, spacedim>::register_data_attach()
   {
     parallel::distributed::Triangulation<dim, spacedim>
-      *non_const_triangulation =
+      *distributed_triangulation =
         const_cast<parallel::distributed::Triangulation<dim, spacedim> *>(
           dynamic_cast<const parallel::distributed::Triangulation<dim, spacedim>
                          *>(&(*triangulation)));
-    (void)non_const_triangulation;
+    (void)distributed_triangulation;
 
-    Assert(non_const_triangulation != nullptr, dealii::ExcNotImplemented());
+    Assert(
+      distributed_triangulation != nullptr,
+      ExcMessage(
+        "Mesh refinement in a non-distributed triangulation is not supported "
+        "by the ParticleHandler class. Either insert particles after mesh "
+        "creation and do not refine afterwards, or use a distributed triangulation."));
 
 #ifdef DEAL_II_WITH_P4EST
-    // Only save and load particles if there are any, we might get here for
-    // example if somebody created a ParticleHandler but generated 0
-    // particles.
-    update_cached_numbers();
+    const auto callback_function =
+      [this](
+        const typename Triangulation<dim, spacedim>::cell_iterator
+          &                                                     cell_iterator,
+        const typename Triangulation<dim, spacedim>::CellStatus cell_status) {
+        return this->pack_callback(cell_iterator, cell_status);
+      };
+
+    handle = distributed_triangulation->register_data_attach(
+      callback_function, /*returns_variable_size_data=*/true);
+#endif
+  }
 
-    if (global_max_particles_per_cell > 0)
-      {
-        const auto callback_function =
-          [this](const typename Triangulation<dim, spacedim>::cell_iterator
-                   &cell_iterator,
-                 const typename Triangulation<dim, spacedim>::CellStatus
-                   cell_status) {
-            return this->store_particles(cell_iterator, cell_status);
-          };
 
-        handle = non_const_triangulation->register_data_attach(
-          callback_function, /*returns_variable_size_data=*/true);
-      }
-#endif
+
+  template <int dim, int spacedim>
+  void
+  ParticleHandler<dim, spacedim>::unpack_after_coarsening_and_refinement()
+  {
+    const bool serialization = false;
+    notify_ready_to_unpack(serialization);
   }
 
 
 
+  template <int dim, int spacedim>
+  void
+  ParticleHandler<dim, spacedim>::deserialize()
+  {
+    const bool serialization = true;
+    notify_ready_to_unpack(serialization);
+  }
+
+
   template <int dim, int spacedim>
   void
   ParticleHandler<dim, spacedim>::register_load_callback_function(
     const bool serialization)
+  {
+    notify_ready_to_unpack(serialization);
+  }
+
+
+
+  template <int dim, int spacedim>
+  void
+  ParticleHandler<dim, spacedim>::notify_ready_to_unpack(
+    const bool serialization)
   {
     parallel::distributed::Triangulation<dim, spacedim>
-      *non_const_triangulation =
+      *distributed_triangulation =
         const_cast<parallel::distributed::Triangulation<dim, spacedim> *>(
           dynamic_cast<const parallel::distributed::Triangulation<dim, spacedim>
                          *>(&(*triangulation)));
-    (void)non_const_triangulation;
+    (void)distributed_triangulation;
+
+    Assert(
+      distributed_triangulation != nullptr,
+      ExcMessage(
+        "Mesh refinement in a non-distributed triangulation is not supported "
+        "by the ParticleHandler class. Either insert particles after mesh "
+        "creation and do not refine afterwards, or use a distributed triangulation."));
 
-    Assert(non_const_triangulation != nullptr, dealii::ExcNotImplemented());
+    // First prepare container for insertion
+    clear();
 
 #ifdef DEAL_II_WITH_P4EST
     // If we are resuming from a checkpoint, we first have to register the
-    // store function again, to set the triangulation in the same state as
-    // before the serialization. Only by this it knows how to deserialize the
-    // data correctly. Only do this if something was actually stored.
-    if (serialization && (global_max_particles_per_cell > 0))
-      {
-        const auto callback_function =
-          [this](const typename Triangulation<dim, spacedim>::cell_iterator
-                   &cell_iterator,
-                 const typename Triangulation<dim, spacedim>::CellStatus
-                   cell_status) {
-            return this->store_particles(cell_iterator, cell_status);
-          };
-
-        handle = non_const_triangulation->register_data_attach(
-          callback_function, /*returns_variable_size_data=*/true);
-      }
+    // store function again, to set the triangulation to the same state as
+    // before the serialization. Only afterwards we know how to deserialize the
+    // data correctly.
+    if (serialization)
+      register_data_attach();
 
     // Check if something was stored and load it
     if (handle != numbers::invalid_unsigned_int)
@@ -1943,14 +2049,13 @@ namespace Particles
             const typename Triangulation<dim, spacedim>::CellStatus cell_status,
             const boost::iterator_range<std::vector<char>::const_iterator>
               &range_iterator) {
-            this->load_particles(cell_iterator, cell_status, range_iterator);
+            this->unpack_callback(cell_iterator, cell_status, range_iterator);
           };
 
-        non_const_triangulation->notify_ready_to_unpack(handle,
-                                                        callback_function);
+        distributed_triangulation->notify_ready_to_unpack(handle,
+                                                          callback_function);
 
-        // Reset handle and update global number of particles. The number
-        // can change because of discarded or newly generated particles
+        // Reset handle and update global numbers.
         handle = numbers::invalid_unsigned_int;
         update_cached_numbers();
       }
@@ -1963,7 +2068,7 @@ namespace Particles
 
   template <int dim, int spacedim>
   std::vector<char>
-  ParticleHandler<dim, spacedim>::store_particles(
+  ParticleHandler<dim, spacedim>::pack_callback(
     const typename Triangulation<dim, spacedim>::cell_iterator &cell,
     const typename Triangulation<dim, spacedim>::CellStatus     status) const
   {
@@ -2016,7 +2121,7 @@ namespace Particles
 
   template <int dim, int spacedim>
   void
-  ParticleHandler<dim, spacedim>::load_particles(
+  ParticleHandler<dim, spacedim>::unpack_callback(
     const typename Triangulation<dim, spacedim>::cell_iterator &    cell,
     const typename Triangulation<dim, spacedim>::CellStatus         status,
     const boost::iterator_range<std::vector<char>::const_iterator> &data_range)
index 015448f0c694bf34c0ebba1a1ff7ae9f2719318f..146058b64ab7df20dc15764b78f34d81e8a7424b 100644 (file)
@@ -98,22 +98,6 @@ test()
 
     Particles::ParticleHandler<dim, spacedim> particle_handler(tr, mapping);
 
-    // TODO: Move this into the Particle handler class. Unfortunately, there are
-    // some interactions with the SolutionTransfer class that prevent us from
-    // doing this at the moment. When doing this, check that transferring a
-    // solution and particles during the same refinement is possible (in
-    // particular that the order of serialization/deserialization is preserved).
-    tr.signals.pre_distributed_refinement.connect(std::bind(
-      &Particles::ParticleHandler<dim,
-                                  spacedim>::register_store_callback_function,
-      &particle_handler));
-
-    tr.signals.post_distributed_refinement.connect(std::bind(
-      &Particles::ParticleHandler<dim,
-                                  spacedim>::register_load_callback_function,
-      &particle_handler,
-      false));
-
     create_regular_particle_distribution(particle_handler, tr);
 
     for (const auto &particle : particle_handler)
@@ -122,7 +106,9 @@ test()
               << std::endl;
 
     // Check that all particles are moved to children
+    particle_handler.prepare_for_coarsening_and_refinement();
     tr.refine_global(1);
+    particle_handler.unpack_after_coarsening_and_refinement();
 
     for (const auto &particle : particle_handler)
       deallog << "After refinement particle id " << particle.get_id()
@@ -133,7 +119,9 @@ test()
     for (auto cell = tr.begin_active(); cell != tr.end(); ++cell)
       cell->set_coarsen_flag();
 
+    particle_handler.prepare_for_coarsening_and_refinement();
     tr.execute_coarsening_and_refinement();
+    particle_handler.unpack_after_coarsening_and_refinement();
 
     for (const auto &particle : particle_handler)
       deallog << "After coarsening particle id " << particle.get_id()
index 1be745e584fe46901be19f6798bbeabdfecc02f5..bec558e0938b11e9ed92bb0c2b5503d067455b02 100644 (file)
@@ -107,22 +107,6 @@ test()
                                                                mapping,
                                                                n_properties);
 
-    // TODO: Move this into the Particle handler class. Unfortunately, there are
-    // some interactions with the SolutionTransfer class that prevent us from
-    // doing this at the moment. When doing this, check that transferring a
-    // solution and particles during the same refinement is possible (in
-    // particular that the order of serialization/deserialization is preserved).
-    tr.signals.pre_distributed_refinement.connect(std::bind(
-      &Particles::ParticleHandler<dim,
-                                  spacedim>::register_store_callback_function,
-      &particle_handler));
-
-    tr.signals.post_distributed_refinement.connect(std::bind(
-      &Particles::ParticleHandler<dim,
-                                  spacedim>::register_load_callback_function,
-      &particle_handler,
-      false));
-
     create_regular_particle_distribution(particle_handler, tr);
 
     for (const auto &particle : particle_handler)
@@ -132,7 +116,9 @@ test()
               << std::endl;
 
     // Check that all particles are moved to children
+    particle_handler.prepare_for_coarsening_and_refinement();
     tr.refine_global(1);
+    particle_handler.unpack_after_coarsening_and_refinement();
 
     for (const auto &particle : particle_handler)
       deallog << "After refinement particle id " << particle.get_id()
@@ -144,7 +130,9 @@ test()
     for (auto cell = tr.begin_active(); cell != tr.end(); ++cell)
       cell->set_coarsen_flag();
 
+    particle_handler.prepare_for_coarsening_and_refinement();
     tr.execute_coarsening_and_refinement();
+    particle_handler.unpack_after_coarsening_and_refinement();
 
     for (const auto &particle : particle_handler)
       deallog << "After coarsening particle id " << particle.get_id()
index e2aef329177aa7ae52c07115e26cfb16ea188a05..84357646070846ddbacc24ef57cd4da354abc6f1 100644 (file)
@@ -54,16 +54,9 @@ test()
   particle_handler.insert_particle(particle, cell);
   particle_handler.update_cached_numbers();
 
-  // initiate data transfer
-  tr.signals.pre_distributed_repartition.connect([&particle_handler]() {
-    particle_handler.register_store_callback_function();
-  });
-
-  tr.signals.post_distributed_repartition.connect([&particle_handler]() {
-    particle_handler.register_load_callback_function(false);
-  });
-
+  particle_handler.prepare_for_coarsening_and_refinement();
   tr.repartition();
+  particle_handler.unpack_after_coarsening_and_refinement();
 
   deallog << "OK" << std::endl;
 }
index 868722978b4ebb1b4e805ffe61ab3819d2102632..dc05749e2e3389c088a815bab8214545d05d1e77 100644 (file)
@@ -321,17 +321,13 @@ namespace Step68
     background_triangulation.refine_global(fluid_refinement);
 
     // In order to consider the particles when repartitioning the triangulation
-    // the algorithm needs to know three things:
+    // the algorithm needs to know how much weight to assign to each cell
+    // (how many particles are in there).
     //
-    // 1. How much weight to assign to each cell (how many particles are in
-    // there);
-    // 2. How to pack the particles before shipping data around;
-    // 3. How to unpack the particles after repartitioning.
-    //
-    // We attach the correct functions to the signals inside
-    // parallel::distributed::Triangulation. These signal will be called every
-    // time the repartition() function is called. These connections only need to
-    // be created once, so we might as well have set them up in the constructor
+    // We attach a weight function to the signal inside
+    // parallel::distributed::Triangulation. This signal will be called every
+    // time the repartition() function is called. This connection only needs to
+    // be created once, so we might as well have set it up in the constructor
     // of this class, but for the purpose of this example we want to group the
     // particle related instructions.
     background_triangulation.signals.cell_weight.connect(
@@ -341,12 +337,6 @@ namespace Step68
         const typename parallel::distributed::Triangulation<dim>::CellStatus
           status) -> unsigned int { return this->cell_weight(cell, status); });
 
-    background_triangulation.signals.pre_distributed_repartition.connect(
-      [this]() { this->particle_handler.register_store_callback_function(); });
-
-    background_triangulation.signals.post_distributed_repartition.connect(
-      [&]() { this->particle_handler.register_load_callback_function(false); });
-
     // This initializes the background triangulation where the particles are
     // living and the number of properties of the particles.
     particle_handler.initialize(background_triangulation, mapping, 1 + dim);
@@ -644,7 +634,10 @@ namespace Step68
     if (Utilities::MPI::this_mpi_process(mpi_communicator) == 0)
       deallog << "Repartitioning triangulation after particle generation"
               << std::endl;
+
+    particle_handler.prepare_for_coarsening_and_refinement();
     background_triangulation.repartition();
+    particle_handler.unpack_after_coarsening_and_refinement();
 
     // We set the initial property of the particles by doing an
     // explicit Euler iteration with a time-step of 0 both in the case
@@ -670,7 +663,10 @@ namespace Step68
 
         if ((discrete_time.get_step_number() % repartition_frequency) == 0)
           {
+            particle_handler.prepare_for_coarsening_and_refinement();
             background_triangulation.repartition();
+            particle_handler.unpack_after_coarsening_and_refinement();
+
             if (interpolated_velocity)
               setup_background_dofs();
           }
index 4ae2991c61340e2258f6f243214bf91264bdb5ba..3562875214f89020f6f6623fe4b05a1b42255d97 100644 (file)
@@ -103,21 +103,13 @@ test()
             << " is in cell " << particle->get_surrounding_cell(tr)
             << std::endl;
 
-  // TODO: Move this into the Particle handler class. Unfortunately, there are
-  // some interactions with the SolutionTransfer class that prevent us from
-  // doing this at the moment. When doing this, check that transferring a
-  // solution and particles during the same refinement is possible (in
-  // particular that the order of serialization/deserialization is preserved).
-  tr.signals.pre_distributed_save.connect(std::bind(
-    &Particles::ParticleHandler<dim,
-                                spacedim>::register_store_callback_function,
-    &particle_handler));
-
   // save data to archive
   std::ostringstream oss;
   {
     boost::archive::text_oarchive oa(oss, boost::archive::no_header);
+
     oa << particle_handler;
+    particle_handler.prepare_for_serialization();
     tr.save("checkpoint");
 
     // archive and stream closed when
@@ -140,17 +132,6 @@ test()
     deallog << "In between particle id " << particle->get_id() << " is in cell "
             << particle->get_surrounding_cell(tr) << std::endl;
 
-
-  // TODO: Move this into the Particle handler class. Unfortunately, there are
-  // some interactions with the SolutionTransfer class that prevent us from
-  // doing this at the moment. When doing this, check that transferring a
-  // solution and particles during the same refinement is possible (in
-  // particular that the order of serialization/deserialization is preserved).
-  tr.signals.post_distributed_load.connect(std::bind(
-    &Particles::ParticleHandler<dim, spacedim>::register_load_callback_function,
-    &particle_handler,
-    true));
-
   // verify correctness of the serialization. Note that the deserialization of
   // the particle handler has to happen before the triangulation (otherwise it
   // does not know if something was stored in the user data of the
@@ -161,6 +142,7 @@ test()
 
     ia >> particle_handler;
     tr.load("checkpoint");
+    particle_handler.deserialize();
   }
 
   for (auto particle = particle_handler.begin();
index 31b49642851b96dc8efa68c121bec38307776419..11aae590990daf05f3e462feda76c6c019f1c368 100644 (file)
@@ -107,16 +107,6 @@ test()
             << " is in cell " << particle->get_surrounding_cell(tr)
             << std::endl;
 
-  // TODO: Move this into the Particle handler class. Unfortunately, there are
-  // some interactions with the SolutionTransfer class that prevent us from
-  // doing this at the moment. When doing this, check that transferring a
-  // solution and particles during the same refinement is possible (in
-  // particular that the order of serialization/deserialization is preserved).
-  tr.signals.pre_distributed_save.connect(std::bind(
-    &Particles::ParticleHandler<dim,
-                                spacedim>::register_store_callback_function,
-    &particle_handler));
-
   // save data to archive
   std::ostringstream oss;
   {
@@ -136,6 +126,7 @@ test()
           oa << particle;
       }
 
+    particle_handler.prepare_for_serialization();
     tr.save("checkpoint");
 
     // archive and stream closed when
@@ -158,17 +149,6 @@ test()
     deallog << "In between particle id " << particle->get_id() << " is in cell "
             << particle->get_surrounding_cell(tr) << std::endl;
 
-
-  // TODO: Move this into the Particle handler class. Unfortunately, there are
-  // some interactions with the SolutionTransfer class that prevent us from
-  // doing this at the moment. When doing this, check that transferring a
-  // solution and particles during the same refinement is possible (in
-  // particular that the order of serialization/deserialization is preserved).
-  tr.signals.post_distributed_load.connect(std::bind(
-    &Particles::ParticleHandler<dim, spacedim>::register_load_callback_function,
-    &particle_handler,
-    true));
-
   // verify correctness of the serialization. Note that the deserialization of
   // the particle handler has to happen before the triangulation (otherwise it
   // does not know if something was stored in the user data of the
@@ -180,6 +160,7 @@ test()
     ia >> particle_handler;
 
     tr.load("checkpoint");
+    particle_handler.deserialize();
 
     // now also read the particles by hand into the stream. the same
     // comment applied as where we were writing them

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.