]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Optimize particle property pool 11314/head
authorRene Gassmoeller <rene.gassmoeller@mailbox.org>
Thu, 3 Dec 2020 21:33:12 +0000 (16:33 -0500)
committerRene Gassmoeller <rene.gassmoeller@mailbox.org>
Sat, 12 Dec 2020 00:29:51 +0000 (19:29 -0500)
include/deal.II/particles/particle.h
include/deal.II/particles/property_pool.h
source/particles/particle.cc
source/particles/particle_handler.cc
source/particles/property_pool.cc
tests/particles/particle_handler_18.cc
tests/serialization/particle_01.cc [new file with mode: 0644]
tests/serialization/particle_01.with_p4est=true.mpirun=1.output [new file with mode: 0644]

index e43dc03b3197650a1ff4fb0da216d2488a206dc1..dcf43716eaa4d0de113dc3ff301112b755dd6855 100644 (file)
@@ -429,7 +429,10 @@ namespace Particles
 
     /**
      * Read the data of this object from a stream for the purpose of
-     * serialization.
+     * serialization. Note that in order to store the properties
+     * correctly, the property pool of this particle has to
+     * be known at the time of reading, i.e. set_property_pool()
+     * has to have been called, before this function is called.
      */
     template <class Archive>
     void
@@ -496,8 +499,18 @@ namespace Particles
 
     if (n_properties > 0)
       {
-        properties = new double[n_properties];
-        ar &boost::serialization::make_array(properties, n_properties);
+        ArrayView<double> properties(get_properties());
+        Assert(
+          properties.size() == n_properties,
+          ExcMessage(
+            "This particle was serialized with " +
+            std::to_string(n_properties) +
+            " properties, but the new property handler provides space for " +
+            std::to_string(properties.size()) +
+            " properties. Deserializing a particle only works for matching property sizes."));
+
+        ar &boost::serialization::make_array(get_properties().data(),
+                                             n_properties);
       }
   }
 
@@ -516,7 +529,8 @@ namespace Particles
     ar &location &reference_location &id &n_properties;
 
     if (n_properties > 0)
-      ar &boost::serialization::make_array(properties, n_properties);
+      ar &boost::serialization::make_array(get_properties().data(),
+                                           n_properties);
   }
 
 
@@ -589,7 +603,7 @@ namespace Particles
 
         ArrayView<double> old_properties = this->get_properties();
         ArrayView<double> new_properties =
-          property_pool->get_properties(new_handle);
+          new_property_pool.get_properties(new_handle);
         std::copy(old_properties.cbegin(),
                   old_properties.cend(),
                   new_properties.begin());
index 9e007e46e73fdd80960b6d789b9a0c4ef8ca2794..0ca579208da398b87ce28f6b7423aea672952e3e 100644 (file)
@@ -20,8 +20,6 @@
 
 #include <deal.II/base/array_view.h>
 
-#include <set>
-
 
 DEAL_II_NAMESPACE_OPEN
 
@@ -54,7 +52,7 @@ namespace Particles
      * uniquely identifies the slot of memory that is reserved for this
      * particle.
      */
-    using Handle = double *;
+    using Handle = unsigned int;
 
     /**
      * Define a default (invalid) value for handles.
@@ -74,20 +72,33 @@ namespace Particles
     ~PropertyPool();
 
     /**
-     * Return a new handle that allows accessing the reserved block
-     * of memory. If the number of properties is zero this will return an
-     * invalid handle.
+     * Clear the dynamic memory allocated by this class. This function
+     * ensures that all memory that had previously been allocated using
+     * allocate_properties_array() has also been returned via
+     * deallocate_properties_array().
+     */
+    void
+    clear();
+
+    /**
+     * Return a new handle that allows accessing the reserved particle
+     * properties. If the number of properties is zero this will return an
+     * invalid handle. Handles can be copied, but after deallocating one
+     * of the copied handles using any of the other copies will cause
+     * undefined behavior.
      */
     Handle
     allocate_properties_array();
 
     /**
      * Mark the properties corresponding to the handle @p handle as
-     * deleted. Calling this function more than once for the same
-     * handle causes undefined behavior.
+     * deleted and invalidate the @p handle. Calling this function
+     * more than once for the same handle (e.g. by calling it again
+     * with a copy of the previously invalidated handle) causes
+     * undefined behavior.
      */
     void
-    deallocate_properties_array(const Handle handle);
+    deallocate_properties_array(Handle &handle);
 
     /**
      * Return an ArrayView to the properties that correspond to the given
@@ -115,14 +126,44 @@ namespace Particles
      */
     const unsigned int n_properties;
 
+    /**
+     * The currently allocated properties (whether assigned to
+     * a particle or available for assignment).
+     */
+    std::vector<double> properties;
+
     /**
      * A collection of handles that have been created by
-     * allocate_properties_array() and have not been destroyed by
-     * deallocate_properties_array().
+     * allocate_properties_array() and have been destroyed by
+     * deallocate_properties_array(). Since the memory is still
+     * allocated these handles can be reused for new particles
+     * to avoid memory allocation.
      */
-    std::set<Handle> currently_open_handles;
+    std::vector<Handle> currently_available_handles;
   };
 
+  /* ---------------------- inline and template functions ------------------ */
+
+  inline ArrayView<double>
+  PropertyPool::get_properties(const Handle handle)
+  {
+    const std::vector<double>::size_type data_index =
+      (handle != invalid_handle) ? handle * n_properties : 0;
+
+    // Ideally we would need to assert that 'handle' has not been deallocated
+    // by searching through 'currently_available_handles'. However, this
+    // is expensive and this function is performance critical, so instead
+    // just check against the properties range, and rely on the fact
+    // that handles are invalidated when handed over to
+    // deallocate_properties_array().
+    Assert(data_index <= properties.size() - n_properties,
+           ExcMessage("Invalid property handle. This can happen if the "
+                      "handle was duplicated and then one copy was deallocated "
+                      "before trying to access the properties."));
+
+    return ArrayView<double>(properties.data() + data_index, n_properties);
+  }
+
 
 } // namespace Particles
 
index 0a68480dc1382eb50aafd9bdff3ad9bb52b1ea4b..984b327dad46982a97caf56d0918bf89d0a51cfa 100644 (file)
@@ -51,12 +51,11 @@ namespace Particles
     , reference_location(particle.get_reference_location())
     , id(particle.get_id())
     , property_pool(particle.property_pool)
-    , properties((particle.has_properties()) ?
-                   property_pool->allocate_properties_array() :
-                   PropertyPool::invalid_handle)
+    , properties(PropertyPool::invalid_handle)
   {
     if (particle.has_properties())
       {
+        properties = property_pool->allocate_properties_array();
         const ArrayView<double> my_properties =
           property_pool->get_properties(properties);
         const ArrayView<const double> their_properties =
index fd3fc0387439c9773c62f389d31f465217e85f51..539415d8c408920b7360dfd18b9c8c6727039ef8 100644 (file)
@@ -87,6 +87,7 @@ namespace Particles
   template <int dim, int spacedim>
   ParticleHandler<dim, spacedim>::ParticleHandler()
     : triangulation()
+    , mapping()
     , property_pool(std::make_unique<PropertyPool>(0))
     , particles()
     , ghost_particles()
@@ -159,6 +160,8 @@ namespace Particles
     initialize(*particle_handler.triangulation,
                *particle_handler.mapping,
                n_properties);
+    property_pool->reserve(particle_handler.particles.size() +
+                           particle_handler.ghost_particles.size());
 
     // copy static members
     global_number_of_particles = particle_handler.global_number_of_particles;
@@ -177,7 +180,6 @@ namespace Particles
     for (auto &particle : *this)
       {
         particle.set_property_pool(*property_pool);
-        particle.set_properties(from_particle->get_properties());
         ++from_particle;
       }
 
@@ -186,7 +188,6 @@ namespace Particles
          ++ghost, ++from_ghost)
       {
         ghost->set_property_pool(*property_pool);
-        ghost->set_properties(from_ghost->get_properties());
       }
   }
 
@@ -209,6 +210,11 @@ namespace Particles
   ParticleHandler<dim, spacedim>::clear_particles()
   {
     particles.clear();
+    ghost_particles.clear();
+
+    // the particle properties have already been deleted by their destructor,
+    // but the memory is still allocated. Return the memory as well.
+    property_pool->clear();
   }
 
 
index e429b9ff6c7445f6ac027d9cf6d1d683a3056634..0ddb28dfdab4a87a389f316111bfe0d68eccc0be 100644 (file)
@@ -20,7 +20,8 @@ DEAL_II_NAMESPACE_OPEN
 
 namespace Particles
 {
-  const PropertyPool::Handle PropertyPool::invalid_handle = nullptr;
+  const PropertyPool::Handle PropertyPool::invalid_handle =
+    static_cast<Handle>(-1);
 
 
   PropertyPool::PropertyPool(const unsigned int n_properties_per_slot)
@@ -30,45 +31,80 @@ namespace Particles
 
   PropertyPool::~PropertyPool()
   {
-    Assert(currently_open_handles.size() == 0,
-           ExcMessage("This property pool currently still holds " +
-                      std::to_string(currently_open_handles.size()) +
-                      " open handles to memory that was allocated "
-                      "via allocate_properties_array() but that has "
-                      "not been returned via deallocate_properties_array()."));
+    clear();
   }
 
 
 
-  PropertyPool::Handle
-  PropertyPool::allocate_properties_array()
+  void
+  PropertyPool::clear()
   {
-    PropertyPool::Handle handle = PropertyPool::invalid_handle;
     if (n_properties > 0)
       {
-        handle = new double[n_properties];
-        currently_open_handles.insert(handle);
+        const unsigned int n_open_handles =
+          properties.size() / n_properties - currently_available_handles.size();
+        (void)n_open_handles;
+        AssertThrow(n_open_handles == 0,
+                    ExcMessage("This property pool currently still holds " +
+                               std::to_string(n_open_handles) +
+                               " open handles to memory that was allocated "
+                               "via allocate_properties_array() but that has "
+                               "not been returned via "
+                               "deallocate_properties_array()."));
       }
 
-    return handle;
+    // Clear vectors and ensure deallocation of memory
+    properties.clear();
+    properties.shrink_to_fit();
+
+    currently_available_handles.clear();
+    currently_available_handles.shrink_to_fit();
   }
 
 
 
-  void
-  PropertyPool::deallocate_properties_array(Handle handle)
+  PropertyPool::Handle
+  PropertyPool::allocate_properties_array()
   {
-    Assert(currently_open_handles.count(handle) == 1, ExcInternalError());
-    currently_open_handles.erase(handle);
-    delete[] handle;
+    Handle handle = invalid_handle;
+    if (n_properties > 0)
+      {
+        if (currently_available_handles.size() > 0)
+          {
+            handle = currently_available_handles.back();
+            currently_available_handles.pop_back();
+          }
+        else
+          {
+            handle = properties.size() / n_properties;
+            properties.resize(properties.size() + n_properties, 0.0);
+          }
+      }
+
+    return handle;
   }
 
 
 
-  ArrayView<double>
-  PropertyPool::get_properties(const Handle handle)
+  void
+  PropertyPool::deallocate_properties_array(Handle &handle)
   {
-    return ArrayView<double>(handle, n_properties);
+    Assert(
+      handle != invalid_handle,
+      ExcMessage(
+        "This handle is invalid and cannot be deallocated. This can happen if the "
+        "handle was deallocated already before calling this function."));
+    currently_available_handles.push_back(handle);
+    handle = invalid_handle;
+
+    // If this was the last handle, resize containers to 0.
+    // This guarantees that all future properties
+    // are allocated in a sorted way without requiring reallocation.
+    if (currently_available_handles.size() * n_properties == properties.size())
+      {
+        currently_available_handles.clear();
+        properties.clear();
+      }
   }
 
 
@@ -76,7 +112,7 @@ namespace Particles
   void
   PropertyPool::reserve(const std::size_t size)
   {
-    (void)size;
+    properties.reserve(size * n_properties);
   }
 
 
index 9f8088f0a648c16addccbf74af17cea1b636b335..dd14775c7696781fa3b5c520c3af55804a6b6eef 100644 (file)
@@ -63,6 +63,12 @@ test()
     Particles::ParticleHandler<dim, spacedim> particle_handler_copy;
     particle_handler_copy.copy_from(particle_handler);
 
+    // Make sure the old particle handler and the property pool
+    // are cleared. This catches problems if the new particles try to access
+    // old memory addresses (this was a bug, fixed in
+    // https://github.com/dealii/dealii/pull/11314)
+    particle_handler.clear();
+
     for (const auto &particle : particle_handler_copy)
       deallog << "After copying particle id " << particle.get_id()
               << " has first property " << particle.get_properties()[0]
diff --git a/tests/serialization/particle_01.cc b/tests/serialization/particle_01.cc
new file mode 100644 (file)
index 0000000..0ddcf88
--- /dev/null
@@ -0,0 +1,102 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2010 - 2018 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+
+// check and illustrate the serialization process for individual particles
+// outside a ParticleHandler class
+
+#include <deal.II/particles/particle.h>
+#include <deal.II/particles/property_pool.h>
+
+#include "serialization.h"
+
+
+
+template <int dim, int spacedim>
+void
+test()
+{
+  // save data to archive
+  std::ostringstream oss;
+  {
+    Point<spacedim> location;
+    Point<dim>      reference_location;
+
+    location[spacedim - 1]      = 0.3;
+    reference_location[dim - 1] = 0.5;
+    const unsigned int      id  = 6;
+    Particles::PropertyPool property_pool(2);
+
+    Particles::Particle<dim, spacedim> particle(location,
+                                                reference_location,
+                                                id);
+    particle.set_property_pool(property_pool);
+    particle.get_properties()[0] = 0.1;
+    particle.get_properties()[1] = 0.7;
+
+    deallog << "Before serialization particle id " << particle.get_id()
+            << " has location " << particle.get_location()
+            << ", has reference location " << particle.get_reference_location()
+            << ", and has properties " << particle.get_properties()[0] << " "
+            << particle.get_properties()[1] << std::endl;
+
+    boost::archive::text_oarchive oa(oss, boost::archive::no_header);
+    oa << particle;
+
+    // archive and stream closed when
+    // destructors are called
+  }
+  deallog << "Serialized string: " << oss.str() << std::endl;
+
+  // verify correctness of the serialization.
+  {
+    std::istringstream            iss(oss.str());
+    boost::archive::text_iarchive ia(iss, boost::archive::no_header);
+
+    Particles::PropertyPool property_pool(2);
+
+    Particles::Particle<dim, spacedim> particle;
+    particle.set_property_pool(property_pool);
+
+    ia >> particle;
+
+    deallog << "After serialization particle id " << particle.get_id()
+            << " has location " << particle.get_location()
+            << ", has reference location " << particle.get_reference_location()
+            << ", and has properties " << particle.get_properties()[0] << " "
+            << particle.get_properties()[1] << std::endl;
+  }
+
+  deallog << "OK" << std::endl << std::endl;
+}
+
+
+int
+main(int argc, char *argv[])
+{
+  Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1);
+
+  MPILogInitAll all;
+
+  deallog.push("2d/2d");
+  test<2, 2>();
+  deallog.pop();
+  deallog.push("2d/3d");
+  test<2, 3>();
+  deallog.pop();
+  deallog.push("3d/3d");
+  test<3, 3>();
+  deallog.pop();
+}
diff --git a/tests/serialization/particle_01.with_p4est=true.mpirun=1.output b/tests/serialization/particle_01.with_p4est=true.mpirun=1.output
new file mode 100644 (file)
index 0000000..ab95148
--- /dev/null
@@ -0,0 +1,19 @@
+
+DEAL:0:2d/2d::Before serialization particle id 6 has location 0.00000 0.300000, has reference location 0.00000 0.500000, and has properties 0.100000 0.700000
+DEAL:0:2d/2d::Serialized string: 0 0 0 0 0 0 2 0 0 0.00000000000000000e+00 2.99999999999999989e-01 2 0.00000000000000000e+00 5.00000000000000000e-01 6 2 1.00000000000000006e-01 6.99999999999999956e-01
+
+DEAL:0:2d/2d::After serialization particle id 6 has location 0.00000 0.300000, has reference location 0.00000 0.500000, and has properties 0.100000 0.700000
+DEAL:0:2d/2d::OK
+DEAL:0:2d/2d::
+DEAL:0:2d/3d::Before serialization particle id 6 has location 0.00000 0.00000 0.300000, has reference location 0.00000 0.500000, and has properties 0.100000 0.700000
+DEAL:0:2d/3d::Serialized string: 0 0 0 0 0 0 3 0 0 0.00000000000000000e+00 0.00000000000000000e+00 2.99999999999999989e-01 0 0 0 0 2 0 0 0.00000000000000000e+00 5.00000000000000000e-01 6 2 1.00000000000000006e-01 6.99999999999999956e-01
+
+DEAL:0:2d/3d::After serialization particle id 6 has location 0.00000 0.00000 0.300000, has reference location 0.00000 0.500000, and has properties 0.100000 0.700000
+DEAL:0:2d/3d::OK
+DEAL:0:2d/3d::
+DEAL:0:3d/3d::Before serialization particle id 6 has location 0.00000 0.00000 0.300000, has reference location 0.00000 0.00000 0.500000, and has properties 0.100000 0.700000
+DEAL:0:3d/3d::Serialized string: 0 0 0 0 0 0 3 0 0 0.00000000000000000e+00 0.00000000000000000e+00 2.99999999999999989e-01 3 0.00000000000000000e+00 0.00000000000000000e+00 5.00000000000000000e-01 6 2 1.00000000000000006e-01 6.99999999999999956e-01
+
+DEAL:0:3d/3d::After serialization particle id 6 has location 0.00000 0.00000 0.300000, has reference location 0.00000 0.00000 0.500000, and has properties 0.100000 0.700000
+DEAL:0:3d/3d::OK
+DEAL:0:3d/3d::

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.