From 56596e3a126f1ca185065933558b597452880af1 Mon Sep 17 00:00:00 2001 From: Bruno Date: Mon, 2 Nov 2020 20:56:10 -0500 Subject: [PATCH] Added unit test Improved documentation --- include/deal.II/particles/particle_handler.h | 5 +- source/particles/particle_handler.cc | 22 ++- tests/particles/particle_handler_19.cc | 157 ++++++++++++++++++ ...handler_19.with_p4est=true.mpirun=2.output | 39 +++++ 4 files changed, 210 insertions(+), 13 deletions(-) create mode 100644 tests/particles/particle_handler_19.cc create mode 100644 tests/particles/particle_handler_19.with_p4est=true.mpirun=2.output diff --git a/include/deal.II/particles/particle_handler.h b/include/deal.II/particles/particle_handler.h index 897adad997..bce703bb59 100644 --- a/include/deal.II/particles/particle_handler.h +++ b/include/deal.II/particles/particle_handler.h @@ -709,7 +709,10 @@ namespace Particles /** * Update all particles that live in cells that are ghost cells to - * other processes. + * other processes. In this context, update means to update the + * location and the properties of the ghost particles assuming that + * the ghost particles have not changed cells. Consquently, this will + * not update the reference location of the particles. */ void update_ghost_particles(); diff --git a/source/particles/particle_handler.cc b/source/particles/particle_handler.cc index f0e665db42..56094fb7c4 100644 --- a/source/particles/particle_handler.cc +++ b/source/particles/particle_handler.cc @@ -1233,14 +1233,12 @@ namespace Particles const auto parallel_triangulation = dynamic_cast *>( &*triangulation); - if (parallel_triangulation != nullptr) + if (parallel_triangulation == nullptr || + dealii::Utilities::MPI::n_mpi_processes( + parallel_triangulation->get_communicator()) == 1) { - if (dealii::Utilities::MPI::n_mpi_processes( - parallel_triangulation->get_communicator()) == 1) - return; + return; } - else - return; #ifdef DEAL_II_WITH_MPI // First clear the current ghost_particle information @@ -1707,8 +1705,8 @@ namespace Particles switch (status) { - case parallel::distributed::Triangulation::CELL_PERSIST: { + case parallel::distributed::Triangulation::CELL_PERSIST: + { auto position_hint = particles.end(); for (const auto &particle : loaded_particles_on_cell) { @@ -1727,8 +1725,8 @@ namespace Particles } break; - case parallel::distributed::Triangulation::CELL_COARSEN: { + case parallel::distributed::Triangulation::CELL_COARSEN: + { typename std::multimap>::iterator position_hint = particles.end(); @@ -1753,8 +1751,8 @@ namespace Particles } break; - case parallel::distributed::Triangulation::CELL_REFINE: { + case parallel::distributed::Triangulation::CELL_REFINE: + { std::vector< typename std::multimap>::iterator> diff --git a/tests/particles/particle_handler_19.cc b/tests/particles/particle_handler_19.cc new file mode 100644 index 0000000000..a0933e5c7e --- /dev/null +++ b/tests/particles/particle_handler_19.cc @@ -0,0 +1,157 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2017 - 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. +// +// --------------------------------------------------------------------- + + + +// like particle_handler_06, but after exchanging the ghost particles +// modifies the location of the particles and updates them through +// the update_ghosts mechanism. This test also adds a single property +// to the particles. This property is also modified and then updated +// through the update_ghosts mechanism. + +#include + +#include + +#include +#include + +#include + +#include "../tests.h" + +template +void +test() +{ + { + parallel::distributed::Triangulation tr(MPI_COMM_WORLD); + + GridGenerator::hyper_cube(tr); + tr.refine_global(2); + MappingQ mapping(1); + + Particles::ParticleHandler particle_handler(tr, mapping, 1); + + Point position; + Point 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 particle( + position, + reference_position, + Utilities::MPI::this_mpi_process(tr.get_communicator())); + typename Triangulation::active_cell_iterator cell = + tr.begin_active(); + particle_handler.insert_particle(particle, cell); + + particle_handler.sort_particles_into_subdomains_and_cells(); + + // Set the properties of the particle to be a unique number + for (auto particle = particle_handler.begin(); + particle != particle_handler.end(); + ++particle) + { + particle->get_properties()[0] = + 10 + Utilities::MPI::this_mpi_process(tr.get_communicator()); + } + + + particle_handler.exchange_ghost_particles(); + + for (auto particle = particle_handler.begin(); + particle != particle_handler.end(); + ++particle) + deallog << "Particle id : " << particle->get_id() + << " location : " << particle->get_location() + << " property : " << particle->get_properties()[0] + << " is local on process : " + << Utilities::MPI::this_mpi_process(tr.get_communicator()) + << std::endl; + + for (auto particle = particle_handler.begin_ghost(); + particle != particle_handler.end_ghost(); + ++particle) + deallog << "Particle id : " << particle->get_id() + << " location : " << particle->get_location() + << " property : " << particle->get_properties()[0] + << " is ghost on process : " + << Utilities::MPI::this_mpi_process(tr.get_communicator()) + << std::endl; + + deallog << "Modifying particles positions and properties" << std::endl; + + // Modify the location of a single particle on processor 0 and update the + // ghosts + for (auto particle = particle_handler.begin(); + particle != particle_handler.end(); + ++particle) + { + auto location = particle->get_location(); + location[0] += 0.1; + particle->get_properties()[0] += 10; + particle->set_location(location); + } + particle_handler.update_ghost_particles(); + + for (auto particle = particle_handler.begin(); + particle != particle_handler.end(); + ++particle) + deallog << "Particle id : " << particle->get_id() + << " location : " << particle->get_location() + << " property : " << particle->get_properties()[0] + << " is local on process : " + << Utilities::MPI::this_mpi_process(tr.get_communicator()) + << std::endl; + + for (auto particle = particle_handler.begin_ghost(); + particle != particle_handler.end_ghost(); + ++particle) + deallog << "Particle id : " << particle->get_id() + << " location : " << particle->get_location() + << " property : " << particle->get_properties()[0] + << " is ghost 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("2d/3d"); + test<2, 3>(); + deallog.pop(); + deallog.push("3d/3d"); + test<3, 3>(); + deallog.pop(); +} diff --git a/tests/particles/particle_handler_19.with_p4est=true.mpirun=2.output b/tests/particles/particle_handler_19.with_p4est=true.mpirun=2.output new file mode 100644 index 0000000000..91fc304f0f --- /dev/null +++ b/tests/particles/particle_handler_19.with_p4est=true.mpirun=2.output @@ -0,0 +1,39 @@ + +DEAL:0:2d/2d::Particle id : 0 location : 0.475000 0.475000 property : 10.0000 is local on process : 0 +DEAL:0:2d/2d::Particle id : 1 location : 0.525000 0.525000 property : 11.0000 is ghost on process : 0 +DEAL:0:2d/2d::Modifying particles positions and properties +DEAL:0:2d/2d::Particle id : 0 location : 0.575000 0.475000 property : 20.0000 is local on process : 0 +DEAL:0:2d/2d::Particle id : 1 location : 0.625000 0.525000 property : 21.0000 is ghost on process : 0 +DEAL:0:2d/2d::OK +DEAL:0:2d/3d::Particle id : 0 location : 0.475000 0.475000 0.00000 property : 10.0000 is local on process : 0 +DEAL:0:2d/3d::Particle id : 1 location : 0.525000 0.525000 0.00000 property : 11.0000 is ghost on process : 0 +DEAL:0:2d/3d::Modifying particles positions and properties +DEAL:0:2d/3d::Particle id : 0 location : 0.575000 0.475000 0.00000 property : 20.0000 is local on process : 0 +DEAL:0:2d/3d::Particle id : 1 location : 0.625000 0.525000 0.00000 property : 21.0000 is ghost on process : 0 +DEAL:0:2d/3d::OK +DEAL:0:3d/3d::Particle id : 0 location : 0.475000 0.475000 0.475000 property : 10.0000 is local on process : 0 +DEAL:0:3d/3d::Particle id : 1 location : 0.525000 0.525000 0.525000 property : 11.0000 is ghost on process : 0 +DEAL:0:3d/3d::Modifying particles positions and properties +DEAL:0:3d/3d::Particle id : 0 location : 0.575000 0.475000 0.475000 property : 20.0000 is local on process : 0 +DEAL:0:3d/3d::Particle id : 1 location : 0.625000 0.525000 0.525000 property : 21.0000 is ghost on process : 0 +DEAL:0:3d/3d::OK + +DEAL:1:2d/2d::Particle id : 1 location : 0.525000 0.525000 property : 11.0000 is local on process : 1 +DEAL:1:2d/2d::Particle id : 0 location : 0.475000 0.475000 property : 10.0000 is ghost on process : 1 +DEAL:1:2d/2d::Modifying particles positions and properties +DEAL:1:2d/2d::Particle id : 1 location : 0.625000 0.525000 property : 21.0000 is local on process : 1 +DEAL:1:2d/2d::Particle id : 0 location : 0.575000 0.475000 property : 20.0000 is ghost on process : 1 +DEAL:1:2d/2d::OK +DEAL:1:2d/3d::Particle id : 1 location : 0.525000 0.525000 0.00000 property : 11.0000 is local on process : 1 +DEAL:1:2d/3d::Particle id : 0 location : 0.475000 0.475000 0.00000 property : 10.0000 is ghost on process : 1 +DEAL:1:2d/3d::Modifying particles positions and properties +DEAL:1:2d/3d::Particle id : 1 location : 0.625000 0.525000 0.00000 property : 21.0000 is local on process : 1 +DEAL:1:2d/3d::Particle id : 0 location : 0.575000 0.475000 0.00000 property : 20.0000 is ghost on process : 1 +DEAL:1:2d/3d::OK +DEAL:1:3d/3d::Particle id : 1 location : 0.525000 0.525000 0.525000 property : 11.0000 is local on process : 1 +DEAL:1:3d/3d::Particle id : 0 location : 0.475000 0.475000 0.475000 property : 10.0000 is ghost on process : 1 +DEAL:1:3d/3d::Modifying particles positions and properties +DEAL:1:3d/3d::Particle id : 1 location : 0.625000 0.525000 0.525000 property : 21.0000 is local on process : 1 +DEAL:1:3d/3d::Particle id : 0 location : 0.575000 0.475000 0.475000 property : 20.0000 is ghost on process : 1 +DEAL:1:3d/3d::OK + -- 2.39.5