From 670fe7e416ff89bd9431619230f56a18ccc4796c Mon Sep 17 00:00:00 2001 From: Rene Gassmoeller Date: Mon, 21 Jun 2021 17:52:48 -0400 Subject: [PATCH] Fix bug in remove particle --- source/particles/particle_handler.cc | 12 +-- tests/particles/remove_particle_01.cc | 108 ++++++++++++++++++++++ tests/particles/remove_particle_01.output | 53 +++++++++++ 3 files changed, 167 insertions(+), 6 deletions(-) create mode 100644 tests/particles/remove_particle_01.cc create mode 100644 tests/particles/remove_particle_01.output diff --git a/source/particles/particle_handler.cc b/source/particles/particle_handler.cc index 81fcefaf00..2269b73c74 100644 --- a/source/particles/particle_handler.cc +++ b/source/particles/particle_handler.cc @@ -357,7 +357,7 @@ namespace Particles ParticleHandler::remove_particle( const ParticleHandler::particle_iterator &particle) { - auto &container = *(particle->particles_on_cell); + auto &particles_on_cell = *(particle->particles_on_cell); // if the particle has an invalid handle (e.g. because it has // been duplicated before calling this function) do not try @@ -368,15 +368,15 @@ namespace Particles property_pool->deregister_particle(handle); } - if (container.size() > 1) + if (particles_on_cell.size() > 1) { - container[particle->particle_index_within_cell] = - std::move(container.back()); - container.resize(particles.size() - 1); + particles_on_cell[particle->particle_index_within_cell] = + std::move(particles_on_cell.back()); + particles_on_cell.resize(particles_on_cell.size() - 1); } else { - container.clear(); + particles_on_cell.clear(); } --local_number_of_particles; diff --git a/tests/particles/remove_particle_01.cc b/tests/particles/remove_particle_01.cc new file mode 100644 index 0000000000..edf8334f50 --- /dev/null +++ b/tests/particles/remove_particle_01.cc @@ -0,0 +1,108 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2017 - 2019 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 the removal of a single particle works as expected + +#include + +#include + +#include +#include + +#include +#include + +#include "../tests.h" + +template +void +test() +{ + { + Triangulation tr; + + GridGenerator::hyper_cube(tr); + tr.refine_global(1); + MappingQ mapping(1); + + Particles::ParticleHandler particle_handler(tr, mapping); + + std::vector> particle_reference_locations(3, Point()); + + for (unsigned int i = 0; i < dim; ++i) + { + particle_reference_locations[0](i) = 0.25; + particle_reference_locations[1](i) = 0.5; + particle_reference_locations[2](i) = 0.75; + } + + Particles::Generators::regular_reference_locations( + tr, particle_reference_locations, particle_handler); + + deallog << "Particle number: " << particle_handler.n_global_particles() + << std::endl; + + for (const auto &cell : tr.active_cell_iterators()) + { + const unsigned int n_particles_to_remove = + cell->active_cell_index() % 4; + + for (unsigned int i = 0; i < n_particles_to_remove; ++i) + { + const unsigned int particle_index_to_remove = + Testing::rand() % particle_handler.n_particles_in_cell(cell); + auto particle_to_remove = + particle_handler.particles_in_cell(cell).begin(); + std::advance(particle_to_remove, particle_index_to_remove); + + deallog << "Removing particle index: " + << particle_to_remove->get_id() + << ". Advanced by: " << particle_index_to_remove + << std::endl; + + particle_handler.remove_particle(particle_to_remove); + } + + deallog << "Cell index: " << cell->active_cell_index() + << ". N Particles: " + << particle_handler.n_particles_in_cell(cell) << std::endl; + } + + for (const auto &particle : particle_handler) + { + deallog << "Particle id: " << particle.get_id() << std::endl; + } + } + + deallog << "OK" << std::endl; +} + + + +int +main(int argc, char *argv[]) +{ + initlog(); + + deallog.push("2d/2d"); + test<2, 2>(); + deallog.pop(); + deallog.push("3d/3d"); + test<3, 3>(); + deallog.pop(); +} diff --git a/tests/particles/remove_particle_01.output b/tests/particles/remove_particle_01.output new file mode 100644 index 0000000000..af546c9f0c --- /dev/null +++ b/tests/particles/remove_particle_01.output @@ -0,0 +1,53 @@ + +DEAL:2d/2d::Particle number: 12 +DEAL:2d/2d::Cell index: 0. N Particles: 3 +DEAL:2d/2d::Removing particle index: 4. Advanced by: 1 +DEAL:2d/2d::Cell index: 1. N Particles: 2 +DEAL:2d/2d::Removing particle index: 7. Advanced by: 1 +DEAL:2d/2d::Removing particle index: 8. Advanced by: 1 +DEAL:2d/2d::Cell index: 2. N Particles: 1 +DEAL:2d/2d::Removing particle index: 10. Advanced by: 1 +DEAL:2d/2d::Removing particle index: 11. Advanced by: 1 +DEAL:2d/2d::Removing particle index: 9. Advanced by: 0 +DEAL:2d/2d::Cell index: 3. N Particles: 0 +DEAL:2d/2d::Particle id: 0 +DEAL:2d/2d::Particle id: 1 +DEAL:2d/2d::Particle id: 2 +DEAL:2d/2d::Particle id: 3 +DEAL:2d/2d::Particle id: 5 +DEAL:2d/2d::Particle id: 6 +DEAL:2d/2d::OK +DEAL:3d/3d::Particle number: 24 +DEAL:3d/3d::Cell index: 0. N Particles: 3 +DEAL:3d/3d::Removing particle index: 4. Advanced by: 1 +DEAL:3d/3d::Cell index: 1. N Particles: 2 +DEAL:3d/3d::Removing particle index: 6. Advanced by: 0 +DEAL:3d/3d::Removing particle index: 7. Advanced by: 1 +DEAL:3d/3d::Cell index: 2. N Particles: 1 +DEAL:3d/3d::Removing particle index: 10. Advanced by: 1 +DEAL:3d/3d::Removing particle index: 9. Advanced by: 0 +DEAL:3d/3d::Removing particle index: 11. Advanced by: 0 +DEAL:3d/3d::Cell index: 3. N Particles: 0 +DEAL:3d/3d::Cell index: 4. N Particles: 3 +DEAL:3d/3d::Removing particle index: 17. Advanced by: 2 +DEAL:3d/3d::Cell index: 5. N Particles: 2 +DEAL:3d/3d::Removing particle index: 19. Advanced by: 1 +DEAL:3d/3d::Removing particle index: 20. Advanced by: 1 +DEAL:3d/3d::Cell index: 6. N Particles: 1 +DEAL:3d/3d::Removing particle index: 22. Advanced by: 1 +DEAL:3d/3d::Removing particle index: 21. Advanced by: 0 +DEAL:3d/3d::Removing particle index: 23. Advanced by: 0 +DEAL:3d/3d::Cell index: 7. N Particles: 0 +DEAL:3d/3d::Particle id: 0 +DEAL:3d/3d::Particle id: 1 +DEAL:3d/3d::Particle id: 2 +DEAL:3d/3d::Particle id: 3 +DEAL:3d/3d::Particle id: 5 +DEAL:3d/3d::Particle id: 8 +DEAL:3d/3d::Particle id: 12 +DEAL:3d/3d::Particle id: 13 +DEAL:3d/3d::Particle id: 14 +DEAL:3d/3d::Particle id: 15 +DEAL:3d/3d::Particle id: 16 +DEAL:3d/3d::Particle id: 18 +DEAL:3d/3d::OK -- 2.39.5