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
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
, 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()
, 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()
// 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;
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;
}
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)
{
&*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)
}
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;
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] =
{
particles_on_cell.clear();
}
-
- --local_number_of_particles;
}
data = particle_it->read_particle_data_from_memory(data);
- ++local_number_of_particles;
+ ++number_of_locally_owned_particles;
return particle_it;
}
if (properties.size() != 0)
particle_it->set_properties(properties);
- ++local_number_of_particles;
+ ++number_of_locally_owned_particles;
return particle_it;
}
types::particle_index
ParticleHandler<dim, spacedim>::n_locally_owned_particles() const
{
- return local_number_of_particles;
+ return number_of_locally_owned_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)
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();
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(
--- /dev/null
+// ---------------------------------------------------------------------
+//
+// 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();
+}
--- /dev/null
+
+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
+