From ccfa0aac59a92ecd079cc3084d31d4fdb06ab84b Mon Sep 17 00:00:00 2001 From: Martin Kronbichler Date: Sun, 4 Jul 2021 22:48:20 +0200 Subject: [PATCH] Test case for new function --- tests/particles/local_particle_index.cc | 137 ++++++++++++++++++ ...icle_index.with_p4est=true.mpirun=3.output | 71 +++++++++ 2 files changed, 208 insertions(+) create mode 100644 tests/particles/local_particle_index.cc create mode 100644 tests/particles/local_particle_index.with_p4est=true.mpirun=3.output diff --git a/tests/particles/local_particle_index.cc b/tests/particles/local_particle_index.cc new file mode 100644 index 0000000000..a45df3cd21 --- /dev/null +++ b/tests/particles/local_particle_index.cc @@ -0,0 +1,137 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2020 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 ParticleAccessor::local_particle_index + +#include + +#include + +#include +#include +#include + +#include +#include +#include + +#include +#include + +#include +#include + +#include "../tests.h" + + + +template +void +test() +{ + // create similar mesh as step-68 + parallel::distributed::Triangulation background_triangulation( + MPI_COMM_WORLD); + + GridGenerator::hyper_cube(background_triangulation, 0, 1); + background_triangulation.refine_global(6 - dim); + + const MappingQGeneric mapping(1); + + Particles::ParticleHandler particle_handler; + particle_handler.initialize(background_triangulation, mapping, 1); + + Point center; + for (unsigned int d = 0; d < dim; ++d) + center[d] = 0.5; + + const double outer_radius = 0.15; + const double inner_radius = 0.01; + + parallel::distributed::Triangulation particle_triangulation( + MPI_COMM_WORLD); + + GridGenerator::hyper_shell( + particle_triangulation, center, inner_radius, outer_radius, 6); + particle_triangulation.refine_global(5 - dim); + + const auto my_bounding_box = GridTools::compute_mesh_predicate_bounding_box( + background_triangulation, IteratorFilters::LocallyOwnedCell()); + const auto global_bounding_boxes = + Utilities::MPI::all_gather(MPI_COMM_WORLD, my_bounding_box); + + std::vector> properties( + particle_triangulation.n_locally_owned_active_cells(), + std::vector(1, 0.)); + Particles::Generators::quadrature_points(particle_triangulation, + QMidpoint(), + global_bounding_boxes, + particle_handler, + mapping, + properties); + + for (unsigned int i = 0; i < 2; ++i) + { + // Write out some information about the local particles + deallog << "Number of global particles: " + << particle_handler.n_global_particles() << std::endl; + deallog << "Number of locally owned particles: " + << particle_handler.n_locally_owned_particles() << std::endl; + deallog << "Maximal local particle index: " + << particle_handler.get_max_local_particle_index() << std::endl; + + std::vector my_particle_numbers( + particle_handler.get_max_local_particle_index(), + numbers::invalid_dof_index); + + // set the id numbers for the particles + unsigned int count = 0; + for (auto particle = particle_handler.begin(); + particle != particle_handler.end(); + ++particle, ++count) + my_particle_numbers[particle->get_local_index()] = particle->get_id(); + + for (auto particle = particle_handler.begin_ghost(); + particle != particle_handler.end_ghost(); + ++particle, ++count) + my_particle_numbers[particle->get_local_index()] = particle->get_id(); + + deallog << "Set information of " << count << " particles" << std::endl; + + unsigned int n_unset = 0; + for (const types::particle_index id : my_particle_numbers) + if (id == numbers::invalid_dof_index) + ++n_unset; + + deallog << "Number of particles not set: " << n_unset << std::endl; + + particle_handler.exchange_ghost_particles(); + } + deallog << "OK" << std::endl; +} + + + +int +main(int argc, char **argv) +{ + Utilities::MPI::MPI_InitFinalize mpi(argc, argv, 1); + MPILogInitAll log; + + test<2>(); + test<3>(); +} diff --git a/tests/particles/local_particle_index.with_p4est=true.mpirun=3.output b/tests/particles/local_particle_index.with_p4est=true.mpirun=3.output new file mode 100644 index 0000000000..9b32d47d62 --- /dev/null +++ b/tests/particles/local_particle_index.with_p4est=true.mpirun=3.output @@ -0,0 +1,71 @@ + +DEAL:0::Number of global particles: 384 +DEAL:0::Number of locally owned particles: 96 +DEAL:0::Maximal local particle index: 96 +DEAL:0::Set information of 96 particles +DEAL:0::Number of particles not set: 0 +DEAL:0::Number of global particles: 384 +DEAL:0::Number of locally owned particles: 96 +DEAL:0::Maximal local particle index: 263 +DEAL:0::Set information of 262 particles +DEAL:0::Number of particles not set: 1 +DEAL:0::OK +DEAL:0::Number of global particles: 384 +DEAL:0::Number of locally owned particles: 96 +DEAL:0::Maximal local particle index: 96 +DEAL:0::Set information of 96 particles +DEAL:0::Number of particles not set: 0 +DEAL:0::Number of global particles: 384 +DEAL:0::Number of locally owned particles: 96 +DEAL:0::Maximal local particle index: 385 +DEAL:0::Set information of 384 particles +DEAL:0::Number of particles not set: 1 +DEAL:0::OK + +DEAL:1::Number of global particles: 384 +DEAL:1::Number of locally owned particles: 192 +DEAL:1::Maximal local particle index: 192 +DEAL:1::Set information of 192 particles +DEAL:1::Number of particles not set: 0 +DEAL:1::Number of global particles: 384 +DEAL:1::Number of locally owned particles: 192 +DEAL:1::Maximal local particle index: 365 +DEAL:1::Set information of 364 particles +DEAL:1::Number of particles not set: 1 +DEAL:1::OK +DEAL:1::Number of global particles: 384 +DEAL:1::Number of locally owned particles: 192 +DEAL:1::Maximal local particle index: 192 +DEAL:1::Set information of 192 particles +DEAL:1::Number of particles not set: 0 +DEAL:1::Number of global particles: 384 +DEAL:1::Number of locally owned particles: 192 +DEAL:1::Maximal local particle index: 385 +DEAL:1::Set information of 384 particles +DEAL:1::Number of particles not set: 1 +DEAL:1::OK + + +DEAL:2::Number of global particles: 384 +DEAL:2::Number of locally owned particles: 96 +DEAL:2::Maximal local particle index: 96 +DEAL:2::Set information of 96 particles +DEAL:2::Number of particles not set: 0 +DEAL:2::Number of global particles: 384 +DEAL:2::Number of locally owned particles: 96 +DEAL:2::Maximal local particle index: 263 +DEAL:2::Set information of 262 particles +DEAL:2::Number of particles not set: 1 +DEAL:2::OK +DEAL:2::Number of global particles: 384 +DEAL:2::Number of locally owned particles: 96 +DEAL:2::Maximal local particle index: 96 +DEAL:2::Set information of 96 particles +DEAL:2::Number of particles not set: 0 +DEAL:2::Number of global particles: 384 +DEAL:2::Number of locally owned particles: 96 +DEAL:2::Maximal local particle index: 385 +DEAL:2::Set information of 384 particles +DEAL:2::Number of particles not set: 1 +DEAL:2::OK + -- 2.39.5