]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Ensure that particles always have a PropertyPool associated with them.
authorWolfgang Bangerth <bangerth@colostate.edu>
Sat, 19 Dec 2020 04:29:06 +0000 (21:29 -0700)
committerWolfgang Bangerth <bangerth@colostate.edu>
Thu, 14 Jan 2021 17:46:30 +0000 (10:46 -0700)
Specifically, do this by providing a global pool for particles not associated
with a dedicated PropertyPool object.

include/deal.II/particles/particle.h
source/particles/particle.cc
source/particles/property_pool.cc

index b6db41ab4161d9c5720a3b93733b5ffb20184fe3..05780de7113c4860e956dd45b3488f938742e8b3 100644 (file)
@@ -151,10 +151,12 @@ namespace Particles
     Particle();
 
     /**
-     * Constructor for Particle, creates a particle with the specified
-     * ID at the specified location. Note that there is no
-     * check for duplicate particle IDs so the user must
-     * make sure the IDs are unique over all processes.
+     * Constructor for Particle. This function creates a particle with the
+     * specified ID at the specified location. Note that there is no check for
+     * duplicate particle IDs so the user must make sure the IDs are unique over
+     * all processes. Data is stored in a global PropertyPool object
+     * (corresponding to the global "heap") but can later be transfered to
+     * another property pool by calling set_property_pool().
      *
      * @param[in] location Initial location of particle.
      * @param[in] reference_location Initial location of the particle
@@ -166,18 +168,20 @@ namespace Particles
              const types::particle_index id);
 
     /**
-     * Copy-Constructor for Particle, creates a particle with exactly the
-     * state of the input argument. Note that since each particle has a
-     * handle for a certain piece of the property memory, and is responsible
-     * for registering and freeing this memory in the property pool this
-     * constructor registers a new chunk, and copies the properties.
+     * Copy-constructor for Particle. This function creates a particle with
+     * exactly the state of the input argument. The copied data is stored in a
+     * global PropertyPool object (corresponding to the global "heap") but can
+     * later be transfered to another property pool by calling
+     * set_property_pool().
      */
     Particle(const Particle<dim, spacedim> &particle);
 
     /**
-     * Constructor for Particle, creates a particle from a data vector.
-     * This constructor is usually called after serializing a particle by
-     * calling the write_data() function.
+     * Constructor for Particle. This function creates a particle from a data
+     * vector. Data is stored in a global PropertyPool object (corresponding to
+     * the global "heap") but can later be transfered to another property pool
+     * by calling set_property_pool(). This constructor is usually called after
+     * serializing a particle by calling the write_data() function.
      *
      * @param[in,out] begin_data A pointer to a memory location from which
      * to read the information that completely describes a particle. This
@@ -186,10 +190,15 @@ namespace Particles
      * the particle information.
      *
      * @param[in,out] property_pool An optional pointer to a property pool
-     * that is used to manage the property data used by this particle. Note that
-     * if a non-null pointer is handed over this constructor assumes @p begin_data
+     * that is used to manage the property data used by this particle. If this
+     * argument is not provided, then a global property pool is used; on the
+     * other hand,
+     * if a non-null pointer is provided, this constructor assumes @p begin_data
      * contains serialized data of the same length and type that is allocated
-     * by @p property_pool.
+     * by @p property_pool. If the data pointer provided here corresponds
+     * to data for a particle that has properties, then this function will only
+     * succeed if a property pool is provided as second argument that is able to
+     * store the correct number of properties per particle.
      */
     Particle(const void *&                      begin_data,
              PropertyPool<dim, spacedim> *const property_pool = nullptr);
@@ -464,6 +473,12 @@ namespace Particles
 #endif
 
   private:
+    /**
+     * A global property pool used when a particle is not associated with
+     * a property pool that belongs to, for example, a ParticleHandler.
+     */
+    static PropertyPool<dim, spacedim> global_property_pool;
+
     /**
      * Current particle location.
      */
@@ -610,27 +625,30 @@ namespace Particles
   {
     // First, we do want to save any properties that may
     // have previously been set, and copy them over to the memory allocated
-    // on the new pool
-    typename PropertyPool<dim, spacedim>::Handle new_handle =
-      PropertyPool<dim, spacedim>::invalid_handle;
-    if (property_pool != nullptr &&
-        property_pool_handle != PropertyPool<dim, spacedim>::invalid_handle)
+    // on the new pool.
+    //
+    // It is possible that a particle currently has no properties -- for
+    // example if it has been created without an associated property
+    // pool (i.e., uses the default global pool which does not store any
+    // properties) but that the new pool has properties. In that case,
+    // there is simply nothing to transfer -- but the register_particle()
+    // call here will make sure that the newly allocated properties are
+    // zero-initialized.
+    const typename PropertyPool<dim, spacedim>::Handle new_handle =
+      new_property_pool.register_particle();
+
+    if (/* old pool */ has_properties())
       {
-        new_handle = new_property_pool.register_particle();
-
-        ArrayView<double> old_properties = this->get_properties();
-        ArrayView<double> new_properties =
+        ArrayView<const double> old_properties = this->get_properties();
+        ArrayView<double>       new_properties =
           new_property_pool.get_properties(new_handle);
         std::copy(old_properties.cbegin(),
                   old_properties.cend(),
                   new_properties.begin());
       }
 
-    // If the particle currently has a reference to properties, then
-    // release those.
-    if (property_pool != nullptr &&
-        property_pool_handle != PropertyPool<dim, spacedim>::invalid_handle)
-      property_pool->deregister_particle(property_pool_handle);
+    // Now release the old memory handle
+    property_pool->deregister_particle(property_pool_handle);
 
 
     // Then set the pointer to the property pool we want to use. Also set the
@@ -645,9 +663,10 @@ namespace Particles
   inline const ArrayView<const double>
   Particle<dim, spacedim>::get_properties() const
   {
-    Assert(has_properties(), ExcInternalError());
-
-    return property_pool->get_properties(property_pool_handle);
+    if (has_properties() == false)
+      return {};
+    else
+      return property_pool->get_properties(property_pool_handle);
   }
 
 
@@ -656,9 +675,15 @@ namespace Particles
   inline bool
   Particle<dim, spacedim>::has_properties() const
   {
-    return (property_pool != nullptr) &&
-           (property_pool_handle !=
-            PropertyPool<dim, spacedim>::invalid_handle);
+    // Particles always have a property pool asssociated with them,
+    // but we can access properties only if there is a valid handle.
+    // The only way a particle can have no valid handle if it has
+    // been moved-from -- but that leaves an object in an invalid
+    // state, and so we can just assert that that can't be the case.
+    Assert((property_pool_handle !=
+            PropertyPool<dim, spacedim>::invalid_handle),
+           ExcInternalError());
+    return (property_pool->n_properties_per_slot() > 0);
   }
 
 } // namespace Particles
index 09038711d682272d5f06e6228473213f497a6ad5..602a513f0dafab050e31f1296ddb9bbfc7bfc3d8 100644 (file)
@@ -21,13 +21,18 @@ DEAL_II_NAMESPACE_OPEN
 
 namespace Particles
 {
+  template <int dim, int spacedim>
+  PropertyPool<dim, spacedim> Particle<dim, spacedim>::global_property_pool = {
+    /* no properties: */ 0};
+
+
   template <int dim, int spacedim>
   Particle<dim, spacedim>::Particle()
     : location(numbers::signaling_nan<Point<spacedim>>())
     , reference_location(numbers::signaling_nan<Point<dim>>())
     , id(0)
-    , property_pool(nullptr)
-    , property_pool_handle(PropertyPool<dim, spacedim>::invalid_handle)
+    , property_pool(&global_property_pool)
+    , property_pool_handle(property_pool->register_particle())
   {}
 
 
@@ -39,8 +44,8 @@ namespace Particles
     : location(location)
     , reference_location(reference_location)
     , id(id)
-    , property_pool(nullptr)
-    , property_pool_handle(PropertyPool<dim, spacedim>::invalid_handle)
+    , property_pool(&global_property_pool)
+    , property_pool_handle(property_pool->register_particle())
   {}
 
 
@@ -51,15 +56,14 @@ namespace Particles
     , reference_location(particle.get_reference_location())
     , id(particle.get_id())
     , property_pool(particle.property_pool)
-    , property_pool_handle(PropertyPool<dim, spacedim>::invalid_handle)
+    , property_pool_handle(property_pool->register_particle())
   {
     if (particle.has_properties())
       {
-        property_pool_handle = property_pool->register_particle();
-        const ArrayView<double> my_properties =
-          property_pool->get_properties(property_pool_handle);
         const ArrayView<const double> their_properties =
           particle.get_properties();
+        const ArrayView<double> my_properties =
+          property_pool->get_properties(property_pool_handle);
 
         std::copy(their_properties.begin(),
                   their_properties.end(),
@@ -73,6 +77,12 @@ namespace Particles
   Particle<dim, spacedim>::Particle(
     const void *&                      data,
     PropertyPool<dim, spacedim> *const new_property_pool)
+    : location(numbers::signaling_nan<Point<spacedim>>())
+    , reference_location(numbers::signaling_nan<Point<dim>>())
+    , id(0)
+    , property_pool(new_property_pool != nullptr ? new_property_pool :
+                                                   &global_property_pool)
+    , property_pool_handle(property_pool->register_particle())
   {
     const types::particle_index *id_data =
       static_cast<const types::particle_index *>(data);
@@ -89,18 +99,11 @@ namespace Particles
       reference_location(i) = *pdata++;
     set_reference_location(reference_location);
 
-    property_pool = new_property_pool;
-    if (property_pool != nullptr)
-      property_pool_handle = property_pool->register_particle();
-    else
-      property_pool_handle = PropertyPool<dim, spacedim>::invalid_handle;
-
     // See if there are properties to load
     if (has_properties())
       {
-        const ArrayView<double> particle_properties =
-          property_pool->get_properties(property_pool_handle);
-        const unsigned int size = particle_properties.size();
+        const ArrayView<double> particle_properties = this->get_properties();
+        const unsigned int      size = particle_properties.size();
         for (unsigned int i = 0; i < size; ++i)
           particle_properties[i] = *pdata++;
       }
@@ -118,6 +121,9 @@ namespace Particles
     , property_pool(std::move(particle.property_pool))
     , property_pool_handle(std::move(particle.property_pool_handle))
   {
+    // We stole the rhs's properties, so we need to invalidate
+    // the handle the rhs holds lest it releases the memory that
+    // we still reference here.
     particle.property_pool_handle = PropertyPool<dim, spacedim>::invalid_handle;
   }
 
@@ -134,9 +140,12 @@ namespace Particles
         id                 = particle.id;
         property_pool      = particle.property_pool;
 
+        Assert(this->property_pool->n_properties_per_slot() ==
+                 particle.property_pool->n_properties_per_slot(),
+               ExcInternalError());
+
         if (particle.has_properties())
           {
-            property_pool_handle = property_pool->register_particle();
             const ArrayView<const double> their_properties =
               particle.get_properties();
             const ArrayView<double> my_properties =
@@ -146,9 +155,8 @@ namespace Particles
                       their_properties.end(),
                       my_properties.begin());
           }
-        else
-          property_pool_handle = PropertyPool<dim, spacedim>::invalid_handle;
       }
+
     return *this;
   }
 
@@ -181,8 +189,10 @@ namespace Particles
   template <int dim, int spacedim>
   Particle<dim, spacedim>::~Particle()
   {
-    if (property_pool != nullptr &&
-        property_pool_handle != PropertyPool<dim, spacedim>::invalid_handle)
+    // If we still hold a handle, release the memory. The only way a
+    // particle can end up not holding a valid handle is if it has been
+    // moved from.
+    if (property_pool_handle != PropertyPool<dim, spacedim>::invalid_handle)
       property_pool->deregister_particle(property_pool_handle);
   }
 
@@ -192,8 +202,12 @@ namespace Particles
   void
   Particle<dim, spacedim>::free_properties()
   {
-    if (property_pool != nullptr &&
-        property_pool_handle != PropertyPool<dim, spacedim>::invalid_handle)
+    // If we still hold a handle, release the memory. The only way a
+    // particle can end up not holding a valid handle is if it has been
+    // moved from.
+    //
+    // deregister_particle() automatically invalidates its argument.
+    if (property_pool_handle != PropertyPool<dim, spacedim>::invalid_handle)
       property_pool->deregister_particle(property_pool_handle);
   }
 
@@ -287,12 +301,6 @@ namespace Particles
   Particle<dim, spacedim>::set_properties(
     const ArrayView<const double> &new_properties)
   {
-    Assert(property_pool != nullptr, ExcInternalError());
-
-    // If we haven't allocated memory yet, do so now
-    if (property_pool_handle == PropertyPool<dim, spacedim>::invalid_handle)
-      property_pool_handle = property_pool->register_particle();
-
     const ArrayView<double> property_values =
       property_pool->get_properties(property_pool_handle);
 
@@ -317,19 +325,6 @@ namespace Particles
   const ArrayView<double>
   Particle<dim, spacedim>::get_properties()
   {
-    Assert(property_pool != nullptr, ExcInternalError());
-
-    // If this particle has no properties yet, allocate and initialize them.
-    if (property_pool_handle == PropertyPool<dim, spacedim>::invalid_handle)
-      {
-        property_pool_handle = property_pool->register_particle();
-
-        ArrayView<double> my_properties =
-          property_pool->get_properties(property_pool_handle);
-
-        std::fill(my_properties.begin(), my_properties.end(), 0.0);
-      }
-
     return property_pool->get_properties(property_pool_handle);
   }
 } // namespace Particles
index 68e2f2f8ff960eb25681eae3a54a28e9ffef793c..8df11ac503f20d9cf946f2bd4bb6c709c93cb8ea 100644 (file)
@@ -83,28 +83,26 @@ namespace Particles
   PropertyPool<dim, spacedim>::register_particle()
   {
     Handle handle = invalid_handle;
-    if (n_properties > 0)
+    if (currently_available_handles.size() > 0)
       {
-        if (currently_available_handles.size() > 0)
-          {
-            handle = currently_available_handles.back();
-            currently_available_handles.pop_back();
-          }
-        else
-          {
-            handle = locations.size();
-
-            // Append new slots to the end of the arrays. Initialize with
-            // invalid points.
-            locations.resize(locations.size() + 1,
-                             numbers::signaling_nan<Point<spacedim>>());
-            reference_locations.resize(reference_locations.size() + 1,
-                                       numbers::signaling_nan<Point<dim>>());
-
-            // In contrast, initialize properties by zero.
-            properties.resize(properties.size() + n_properties, 0.0);
-          }
+        handle = currently_available_handles.back();
+        currently_available_handles.pop_back();
       }
+    else
+      {
+        handle = locations.size();
+
+        locations.resize(locations.size() + 1);
+        reference_locations.resize(reference_locations.size() + 1);
+        properties.resize(properties.size() + n_properties);
+      }
+
+    // Then initialize whatever slot we have taken with invalid locations,
+    // but initialize properties with zero.
+    set_location(handle, numbers::signaling_nan<Point<spacedim>>());
+    set_reference_location(handle, numbers::signaling_nan<Point<dim>>());
+    for (double &x : get_properties(handle))
+      x = 0;
 
     return handle;
   }

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.