From b17d42e8f5e045835a6cd963344008018544a9d5 Mon Sep 17 00:00:00 2001 From: Rene Gassmoeller Date: Wed, 14 Aug 2019 13:44:34 -0700 Subject: [PATCH] Implement regular particle generation --- .../deal.II/particles/particle_generator.h | 52 +++++++++++ source/particles/CMakeLists.txt | 2 + source/particles/particle_generator.cc | 91 +++++++++++++++++++ source/particles/particle_generator.inst.in | 37 ++++++++ tests/particles/particle_generator_01.cc | 82 +++++++++++++++++ ...rticle_generator_01.with_p4est=true.output | 13 +++ 6 files changed, 277 insertions(+) create mode 100644 include/deal.II/particles/particle_generator.h create mode 100644 source/particles/particle_generator.cc create mode 100644 source/particles/particle_generator.inst.in create mode 100644 tests/particles/particle_generator_01.cc create mode 100644 tests/particles/particle_generator_01.with_p4est=true.output diff --git a/include/deal.II/particles/particle_generator.h b/include/deal.II/particles/particle_generator.h new file mode 100644 index 0000000000..e4cab73ff9 --- /dev/null +++ b/include/deal.II/particles/particle_generator.h @@ -0,0 +1,52 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 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. +// +// --------------------------------------------------------------------- + +#ifndef dealii_particles_particle_generator_h +#define dealii_particles_particle_generator_h + +#include + +#include +#include + +#include + +DEAL_II_NAMESPACE_OPEN + +namespace Particles +{ + /** + * A namespace that contains all classes that are related to the particle + * generation. + */ + namespace Generator + { + /** + * A function that generates particles in every cell at specified @p particle_reference_locations. + */ + template + void + regular_reference_locations( + const parallel::distributed::Triangulation &triangulation, + const std::vector> & particle_reference_locations, + ParticleHandler &particle_handler, + const Mapping & mapping = MappingQ1()); + + } // namespace Generator +} // namespace Particles + +DEAL_II_NAMESPACE_CLOSE + +#endif diff --git a/source/particles/CMakeLists.txt b/source/particles/CMakeLists.txt index d911d3ef1c..0bbe6baaba 100644 --- a/source/particles/CMakeLists.txt +++ b/source/particles/CMakeLists.txt @@ -21,6 +21,7 @@ SET(_src particle_accessor.cc particle_iterator.cc particle_handler.cc + particle_generator.cc property_pool.cc ) @@ -29,6 +30,7 @@ SET(_inst particle_accessor.inst.in particle_iterator.inst.in particle_handler.inst.in + particle_generator.inst.in ) FILE(GLOB _header diff --git a/source/particles/particle_generator.cc b/source/particles/particle_generator.cc new file mode 100644 index 0000000000..de923d07ec --- /dev/null +++ b/source/particles/particle_generator.cc @@ -0,0 +1,91 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 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. +// +// --------------------------------------------------------------------- + +#include + +#include + +DEAL_II_NAMESPACE_OPEN + +namespace Particles +{ + namespace Generator + { + template + void + regular_reference_locations( + const parallel::distributed::Triangulation &triangulation, + const std::vector> & particle_reference_locations, + ParticleHandler &particle_handler, + const Mapping & mapping) + { + types::particle_index n_particles_to_generate = + triangulation.n_locally_owned_active_cells() * + particle_reference_locations.size(); + + // The local start index is the number of all particles + // generated on lower MPI ranks. + types::particle_index local_start_index = 0; + MPI_Scan(&n_particles_to_generate, + &local_start_index, + 1, + DEAL_II_PARTICLE_INDEX_MPI_TYPE, + MPI_SUM, + triangulation.get_communicator()); + local_start_index -= n_particles_to_generate; + + types::particle_index particle_index = local_start_index; + + typename Triangulation::active_cell_iterator + cell = triangulation.begin_active(), + endc = triangulation.end(); + + for (; cell != endc; cell++) + { + if (cell->is_locally_owned()) + { + for (typename std::vector>::const_iterator + itr_particles_in_unit_cell = + particle_reference_locations.begin(); + itr_particles_in_unit_cell != + particle_reference_locations.end(); + itr_particles_in_unit_cell++) + { + const Point position_real = + mapping.transform_unit_to_real_cell( + cell, *itr_particles_in_unit_cell); + + const Particle particle( + position_real, *itr_particles_in_unit_cell, particle_index); + particle_handler.insert_particle(particle, cell); + ++particle_index; + } + } + } + + particle_handler.update_cached_numbers(); + } + } // namespace Generator +} // namespace Particles + + + +DEAL_II_NAMESPACE_CLOSE + +DEAL_II_NAMESPACE_OPEN + +#include "particle_generator.inst" + +DEAL_II_NAMESPACE_CLOSE diff --git a/source/particles/particle_generator.inst.in b/source/particles/particle_generator.inst.in new file mode 100644 index 0000000000..c642d282fb --- /dev/null +++ b/source/particles/particle_generator.inst.in @@ -0,0 +1,37 @@ +// --------------------------------------------------------------------- +// +// 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. +// +// --------------------------------------------------------------------- + + +for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS) + { +#if deal_II_dimension <= deal_II_space_dimension + namespace Particles + \{ + namespace Generator + \{ + template void + regular_reference_locations( + const parallel::distributed::Triangulation + &triangulation, + const std::vector> + &particle_reference_locations, + ParticleHandler + &particle_handler, + const Mapping &mapping); + \} + \} +#endif + } diff --git a/tests/particles/particle_generator_01.cc b/tests/particles/particle_generator_01.cc new file mode 100644 index 0000000000..4e1265caa8 --- /dev/null +++ b/tests/particles/particle_generator_01.cc @@ -0,0 +1,82 @@ +// --------------------------------------------------------------------- +// +// 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. +// +// --------------------------------------------------------------------- + + + +// check the creation and destruction of particle within the particle handler +// class using a particle generator. + +#include + +#include + +#include +#include + +#include +#include + +#include "../tests.h" + +template +void +test() +{ + { + parallel::distributed::Triangulation tr(MPI_COMM_WORLD); + + GridGenerator::hyper_cube(tr); + MappingQ mapping(1); + + Particles::ParticleHandler particle_handler(tr, mapping); + + std::vector> particle_reference_locations(1, Point()); + Particles::Generator::regular_reference_locations( + tr, particle_reference_locations, particle_handler); + + deallog << "Particle number: " << particle_handler.n_global_particles() + << std::endl; + + for (const auto &particle : particle_handler) + { + deallog << "Particle location: " << particle.get_location() + << std::endl; + deallog << "Particle reference location: " + << particle.get_reference_location() << std::endl; + } + } + + deallog << "OK" << std::endl; +} + + + +int +main(int argc, char *argv[]) +{ + Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1); + + initlog(); + + 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_generator_01.with_p4est=true.output b/tests/particles/particle_generator_01.with_p4est=true.output new file mode 100644 index 0000000000..89addf53f3 --- /dev/null +++ b/tests/particles/particle_generator_01.with_p4est=true.output @@ -0,0 +1,13 @@ + +DEAL:2d/2d::Particle number: 1 +DEAL:2d/2d::Particle location: 0.00000 0.00000 +DEAL:2d/2d::Particle reference location: 0.00000 0.00000 +DEAL:2d/2d::OK +DEAL:2d/3d::Particle number: 1 +DEAL:2d/3d::Particle location: 0.00000 0.00000 0.00000 +DEAL:2d/3d::Particle reference location: 0.00000 0.00000 +DEAL:2d/3d::OK +DEAL:3d/3d::Particle number: 1 +DEAL:3d/3d::Particle location: 0.00000 0.00000 0.00000 +DEAL:3d/3d::Particle reference location: 0.00000 0.00000 0.00000 +DEAL:3d/3d::OK -- 2.39.5