]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Clarify locally owned particles 12485/head
authorRene Gassmoeller <rene.gassmoeller@mailbox.org>
Tue, 22 Jun 2021 14:31:27 +0000 (10:31 -0400)
committerRene Gassmoeller <rene.gassmoeller@mailbox.org>
Wed, 23 Jun 2021 18:22:14 +0000 (14:22 -0400)
include/deal.II/particles/particle_handler.h
include/deal.II/particles/property_pool.h
source/particles/particle_handler.cc
source/particles/property_pool.cc
tests/particles/exchange_ghosts_01.cc [new file with mode: 0644]
tests/particles/exchange_ghosts_01.with_p4est=true.mpirun=2.output [new file with mode: 0644]

index a674d73d546d37d9a8843e4d667f8fb871025c0f..865d94b20558c8114eccb9b6a56f8cbd662640f3 100644 (file)
@@ -830,10 +830,10 @@ namespace Particles
     types::particle_index global_number_of_particles;
 
     /**
-     * This variable stores how many particles are stored locally. It is
+     * This variable stores how many particles are locally owned. It is
      * calculated by update_cached_numbers().
      */
-    types::particle_index local_number_of_particles;
+    types::particle_index number_of_locally_owned_particles;
 
     /**
      * The maximum number of particles per cell in the global domain. This
index 3f74a8167c4b8d7b83448fc8e33179a477f1d3d6..944d2ca6432e177eeed00b8784e319cf6060f0ea 100644 (file)
@@ -211,6 +211,12 @@ namespace Particles
     unsigned int
     n_properties_per_slot() const;
 
+    /**
+     * Return how many slots are currently registered in the pool.
+     */
+    unsigned int
+    n_registered_slots() const;
+
     /**
      * This function makes sure that all internally stored memory blocks
      * are sorted in the same order as one would loop over the @p handles_to_sort
index 2269b73c74bba32315b5b7f2d3e173a96ef986f1..a8cb1e29b3e3eee32c5add65b5460792c1fb906b 100644 (file)
@@ -58,7 +58,7 @@ namespace Particles
     , property_pool(std::make_unique<PropertyPool<dim, spacedim>>(0))
     , particles()
     , global_number_of_particles(0)
-    , local_number_of_particles(0)
+    , number_of_locally_owned_particles(0)
     , global_max_particles_per_cell(0)
     , next_free_particle_index(0)
     , size_callback()
@@ -79,7 +79,7 @@ namespace Particles
     , property_pool(std::make_unique<PropertyPool<dim, spacedim>>(n_properties))
     , particles(triangulation.n_active_cells())
     , global_number_of_particles(0)
-    , local_number_of_particles(0)
+    , number_of_locally_owned_particles(0)
     , global_max_particles_per_cell(0)
     , next_free_particle_index(0)
     , size_callback()
@@ -167,7 +167,9 @@ namespace Particles
 
     // copy static members
     global_number_of_particles = particle_handler.global_number_of_particles;
-    local_number_of_particles  = particle_handler.local_number_of_particles;
+    number_of_locally_owned_particles =
+      particle_handler.number_of_locally_owned_particles;
+
     global_max_particles_per_cell =
       particle_handler.global_max_particles_per_cell;
     next_free_particle_index = particle_handler.next_free_particle_index;
@@ -185,10 +187,10 @@ namespace Particles
   ParticleHandler<dim, spacedim>::clear()
   {
     clear_particles();
-    global_number_of_particles    = 0;
-    local_number_of_particles     = 0;
-    next_free_particle_index      = 0;
-    global_max_particles_per_cell = 0;
+    global_number_of_particles        = 0;
+    number_of_locally_owned_particles = 0;
+    next_free_particle_index          = 0;
+    global_max_particles_per_cell     = 0;
   }
 
 
@@ -217,18 +219,22 @@ namespace Particles
   void
   ParticleHandler<dim, spacedim>::update_cached_numbers()
   {
-    local_number_of_particles                          = 0;
+    number_of_locally_owned_particles                  = 0;
     types::particle_index local_max_particle_index     = 0;
     types::particle_index local_max_particles_per_cell = 0;
 
-    for (const auto &particles_in_cell : particles)
+    for (const auto &cell : triangulation->active_cell_iterators())
       {
+        const auto &particles_in_cell = particles[cell->active_cell_index()];
+
         const types::particle_index n_particles_in_cell =
           particles_in_cell.size();
 
         local_max_particles_per_cell =
           std::max(local_max_particles_per_cell, n_particles_in_cell);
-        local_number_of_particles += n_particles_in_cell;
+
+        if (cell->is_locally_owned())
+          number_of_locally_owned_particles += n_particles_in_cell;
 
         for (const auto &particle : particles_in_cell)
           {
@@ -243,7 +249,7 @@ namespace Particles
             &*triangulation))
       {
         global_number_of_particles = dealii::Utilities::MPI::sum(
-          local_number_of_particles,
+          number_of_locally_owned_particles,
           parallel_triangulation->get_communicator());
 
         if (global_number_of_particles == 0)
@@ -267,7 +273,7 @@ namespace Particles
       }
     else
       {
-        global_number_of_particles = local_number_of_particles;
+        global_number_of_particles = number_of_locally_owned_particles;
         next_free_particle_index =
           global_number_of_particles == 0 ? 0 : local_max_particle_index + 1;
         global_max_particles_per_cell = local_max_particles_per_cell;
@@ -368,6 +374,12 @@ namespace Particles
         property_pool->deregister_particle(handle);
       }
 
+    // need to reduce the cached number before deleting, because the iterator
+    // may be invalid after removing the particle even if only
+    // accessing the cell
+    if (particle->get_surrounding_cell()->is_locally_owned())
+      --number_of_locally_owned_particles;
+
     if (particles_on_cell.size() > 1)
       {
         particles_on_cell[particle->particle_index_within_cell] =
@@ -378,8 +390,6 @@ namespace Particles
       {
         particles_on_cell.clear();
       }
-
-    --local_number_of_particles;
   }
 
 
@@ -480,7 +490,7 @@ namespace Particles
 
     data = particle_it->read_particle_data_from_memory(data);
 
-    ++local_number_of_particles;
+    ++number_of_locally_owned_particles;
 
     return particle_it;
   }
@@ -519,7 +529,7 @@ namespace Particles
     if (properties.size() != 0)
       particle_it->set_properties(properties);
 
-    ++local_number_of_particles;
+    ++number_of_locally_owned_particles;
 
     return particle_it;
   }
@@ -894,7 +904,7 @@ namespace Particles
   types::particle_index
   ParticleHandler<dim, spacedim>::n_locally_owned_particles() const
   {
-    return local_number_of_particles;
+    return number_of_locally_owned_particles;
   }
 
 
@@ -1296,7 +1306,7 @@ namespace Particles
 
     // now make sure particle data is sorted in order of iteration
     std::vector<typename PropertyPool<dim, spacedim>::Handle> unsorted_handles;
-    unsorted_handles.reserve(local_number_of_particles);
+    unsorted_handles.reserve(property_pool->n_registered_slots());
 
     typename PropertyPool<dim, spacedim>::Handle sorted_handle = 0;
     for (auto &particles_in_cell : particles)
@@ -1354,8 +1364,7 @@ namespace Particles
       parallel_triangulation->ghost_owners();
     for (const auto ghost_owner : ghost_owners)
       ghost_particles_cache.ghost_particles_by_domain[ghost_owner].reserve(
-        static_cast<typename std::vector<particle_iterator>::size_type>(
-          local_number_of_particles * 0.25));
+        n_locally_owned_particles() / 4);
 
     const std::vector<std::set<unsigned int>> vertex_to_neighbor_subdomain =
       triangulation_cache->get_vertex_to_neighbor_subdomain();
index 1189a664b1982848d609d4b3ecc6a2d1e71f23b4..9a03f3195f62b150fe31e757f8810c6b583afe44 100644 (file)
@@ -167,6 +167,29 @@ namespace Particles
     return n_properties;
   }
 
+
+
+  template <int dim, int spacedim>
+  unsigned int
+  PropertyPool<dim, spacedim>::n_registered_slots() const
+  {
+    Assert(locations.size() == reference_locations.size(),
+           ExcMessage("Number of registered locations is not equal to number "
+                      "of registered reference locations."));
+
+    Assert(locations.size() == ids.size(),
+           ExcMessage("Number of registered locations is not equal to number "
+                      "of registered ids."));
+
+    Assert(locations.size() * n_properties == properties.size(),
+           ExcMessage("Number of registered locations is not equal to number "
+                      "of registered property slots."));
+
+    return locations.size() - currently_available_handles.size();
+  }
+
+
+
   template <int dim, int spacedim>
   void
   PropertyPool<dim, spacedim>::sort_memory_slots(
diff --git a/tests/particles/exchange_ghosts_01.cc b/tests/particles/exchange_ghosts_01.cc
new file mode 100644 (file)
index 0000000..a5dbe16
--- /dev/null
@@ -0,0 +1,124 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2021 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.
+//
+// ---------------------------------------------------------------------
+
+
+
+// like particle_handler_06, but tests that the particle handlers function
+// n_locally_owned_particles only counts locally owned particles.
+
+#include <deal.II/distributed/tria.h>
+
+#include <deal.II/fe/mapping_q.h>
+
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/grid_tools.h>
+
+#include <deal.II/particles/particle_handler.h>
+
+#include "../tests.h"
+
+template <int dim, int spacedim>
+void
+test()
+{
+  {
+    parallel::distributed::Triangulation<dim, spacedim> tr(MPI_COMM_WORLD);
+
+    GridGenerator::hyper_cube(tr);
+    tr.refine_global(2);
+    MappingQ<dim, spacedim> mapping(1);
+
+    // both processes create a particle handler with a particle in a
+    // position that is in a ghost cell to the other process
+    Particles::ParticleHandler<dim, spacedim> particle_handler(tr, mapping);
+
+    Point<spacedim> position;
+    Point<dim>      reference_position;
+
+    if (Utilities::MPI::this_mpi_process(tr.get_communicator()) == 0)
+      for (unsigned int i = 0; i < dim; ++i)
+        position(i) = 0.475;
+    else
+      for (unsigned int i = 0; i < dim; ++i)
+        position(i) = 0.525;
+
+    Particles::Particle<dim, spacedim> particle(
+      position,
+      reference_position,
+      Utilities::MPI::this_mpi_process(tr.get_communicator()));
+
+    // We give a local random cell hint to check that sorting and
+    // transferring ghost particles works.
+    typename Triangulation<dim, spacedim>::active_cell_iterator cell =
+      tr.begin_active();
+    while (!cell->is_locally_owned())
+      ++cell;
+
+    particle_handler.insert_particle(particle, cell);
+
+    particle_handler.sort_particles_into_subdomains_and_cells();
+
+    deallog << "Before ghost exchange: "
+            << particle_handler.n_locally_owned_particles()
+            << " locally owned particles on process "
+            << Utilities::MPI::this_mpi_process(tr.get_communicator())
+            << std::endl;
+
+    deallog << "Before ghost exchange: "
+            << particle_handler.get_property_pool().n_registered_slots()
+            << " stored particles on process "
+            << Utilities::MPI::this_mpi_process(tr.get_communicator())
+            << std::endl;
+
+    particle_handler.exchange_ghost_particles();
+
+    // Usually this update is unnecessary, because exchanging ghost particles
+    // does not change anything about the cache. However, there was a bug
+    // that ghost particles were counted as locally owned, so make sure
+    // this does not happen again.
+    particle_handler.update_cached_numbers();
+
+    deallog << "After ghost exchange: "
+            << particle_handler.n_locally_owned_particles()
+            << " locally owned particles on process "
+            << Utilities::MPI::this_mpi_process(tr.get_communicator())
+            << std::endl;
+
+    deallog << "After ghost exchange: "
+            << particle_handler.get_property_pool().n_registered_slots()
+            << " stored particles on process "
+            << Utilities::MPI::this_mpi_process(tr.get_communicator())
+            << std::endl;
+  }
+
+  deallog << "OK" << 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("3d/3d");
+  test<3, 3>();
+  deallog.pop();
+}
diff --git a/tests/particles/exchange_ghosts_01.with_p4est=true.mpirun=2.output b/tests/particles/exchange_ghosts_01.with_p4est=true.mpirun=2.output
new file mode 100644 (file)
index 0000000..8e729f8
--- /dev/null
@@ -0,0 +1,23 @@
+
+DEAL:0:2d/2d::Before ghost exchange: 1 locally owned particles on process 0
+DEAL:0:2d/2d::Before ghost exchange: 1 stored particles on process 0
+DEAL:0:2d/2d::After ghost exchange: 1 locally owned particles on process 0
+DEAL:0:2d/2d::After ghost exchange: 2 stored particles on process 0
+DEAL:0:2d/2d::OK
+DEAL:0:3d/3d::Before ghost exchange: 1 locally owned particles on process 0
+DEAL:0:3d/3d::Before ghost exchange: 1 stored particles on process 0
+DEAL:0:3d/3d::After ghost exchange: 1 locally owned particles on process 0
+DEAL:0:3d/3d::After ghost exchange: 2 stored particles on process 0
+DEAL:0:3d/3d::OK
+
+DEAL:1:2d/2d::Before ghost exchange: 1 locally owned particles on process 1
+DEAL:1:2d/2d::Before ghost exchange: 1 stored particles on process 1
+DEAL:1:2d/2d::After ghost exchange: 1 locally owned particles on process 1
+DEAL:1:2d/2d::After ghost exchange: 2 stored particles on process 1
+DEAL:1:2d/2d::OK
+DEAL:1:3d/3d::Before ghost exchange: 1 locally owned particles on process 1
+DEAL:1:3d/3d::Before ghost exchange: 1 stored particles on process 1
+DEAL:1:3d/3d::After ghost exchange: 1 locally owned particles on process 1
+DEAL:1:3d/3d::After ghost exchange: 2 stored particles on process 1
+DEAL:1:3d/3d::OK
+

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.