From: Luca Heltai Date: Thu, 21 Nov 2019 18:57:16 +0000 (+0100) Subject: Generate particles on support or quadrature points on other grid X-Git-Tag: v9.2.0-rc1~727^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=34a04816f1b997c418d6248212f40ca3f95bdd8b;p=dealii.git Generate particles on support or quadrature points on other grid Functions that generates particles at the location of the support points of a possibly non-matching grid or at the quadrature points of a possibly non-matching grid The total number of particles that is added to the particle_handler object is the number of dofs of the DoFHandler that is passed that are within the triangulation and whose components are within the ComponentMask. Co-authored-by: Bruno Blais Co-authored-by: Luca Heltai --- diff --git a/include/deal.II/particles/generators.h b/include/deal.II/particles/generators.h index bac34ab166..05f966d0ec 100644 --- a/include/deal.II/particles/generators.h +++ b/include/deal.II/particles/generators.h @@ -19,9 +19,12 @@ #include #include +#include #include +#include + #include #include @@ -162,6 +165,87 @@ namespace Particles StaticMappingQ1::mapping, const unsigned int random_number_seed = 5432); + + /** + * A function that generates particles at the locations of the support + * points of a DoFHandler, possibly based on a different Triangulation with + * respect to the one used to construct the ParticleHandler. + * The total number of particles that is added to the @p particle_handler object is + * the number of dofs of the DoFHandler that is passed that are within the + * triangulation and whose components are within the ComponentMask. + * This function uses insert_global_particles and consequently may induce + * considerable mpi communication overhead. + * + * @param[in] dof_handler A DOF handler that may live on another + * triangulation that is used to establsh the positions of the particles. + * + * @param[in] A vector that contains all the bounding boxes for all + * processors. This vector can be established by first using + * 'GridTools::compute_mesh_predicate_bounding_box()' and gathering all the + * bounding boxes using 'Utilities::MPI::all_gather(). + * + * @param[in,out] particle_handler The particle handler that will take + * ownership of the generated particles. The particles that are generated + * will be appended to the particles currently owned by the particle + * handler. + * + * @param[in] mapping An optional mapping object that is used to map + * the DOF locations. If no mapping is provided a MappingQ1 is assumed. + * + * @param[in] components Component mask that decides which subset of the + * support points of the dof_handler are used to generate the particles. + * + * @author Bruno Blais, Luca Heltai, 2019 + */ + template + void + dof_support_points(const DoFHandler &dof_handler, + const std::vector>> + & global_bounding_boxes, + ParticleHandler &particle_handler, + const Mapping & mapping = + StaticMappingQ1::mapping, + const ComponentMask &components = ComponentMask()); + + /** + * A function that generates particles at the locations of the quadrature + * points of a Triangulation. This Triangulation can be different + * from the one used to construct the ParticleHandler. + * The total number of particles that is added to the @p particle_handler object is the + * number of cells multiplied by the number of particle reference locations + * which are generally constructed using a quadrature. + * This function uses insert_global_particles and consequently may + * induce considerable mpi communication overhead. + * + * @param[in] triangulation The possibly non-matching triangulation which is + * used to insert the particles into the domain. + * + * @param[in] quadrature A quadrature whose reference location are used to + * insert the particles within the cells. + * + * @param[in] global_bounding_boxes A vector that contains all the bounding + * boxes for all processors. This vector can be established by first using + * 'GridTools::compute_mesh_predicate_bounding_box()' and gathering all the + * bounding boxes using 'Utilities::MPI::all_gather. + * + * @param[in,out] particle_handler The particle handler that will take + * ownership of the generated particles. + * + * @param[in] mapping An optional mapping object that is used to map + * the quadrature locations. If no mapping is provided a MappingQ1 is + * assumed. + * + * @author Bruno Blais, Luca Heltai, 2019 + */ + template + void + quadrature_points(const Triangulation &particle_tria, + const Quadrature & quadrature, + const std::vector>> + & global_bounding_boxes, + ParticleHandler &particle_handler, + const Mapping & mapping = + StaticMappingQ1::mapping); } // namespace Generators } // namespace Particles diff --git a/source/particles/generators.cc b/source/particles/generators.cc index 9619fa0e25..1d3713b255 100644 --- a/source/particles/generators.cc +++ b/source/particles/generators.cc @@ -18,9 +18,15 @@ #include #include +#include + #include #include +#include +#include +#include + #include DEAL_II_NAMESPACE_OPEN @@ -391,6 +397,75 @@ namespace Particles particle_handler.insert_particles(particles); } } + + template + void + dof_support_points(const DoFHandler &dof_handler, + const std::vector>> + & global_bounding_boxes, + ParticleHandler &particle_handler, + const Mapping & mapping, + const ComponentMask & components) + { + const auto &fe = dof_handler.get_fe(); + + // Take care of components + const ComponentMask mask = + (components.size() == 0 ? ComponentMask(fe.n_components(), true) : + components); + + std::map> support_points_map; + + DoFTools::map_dofs_to_support_points(mapping, + dof_handler, + support_points_map, + mask); + + // Generate the vector of points from the map + // Memory is reserved for efficiency reasons + std::vector> support_points_vec; + support_points_vec.reserve(support_points_map.size()); + for (auto const &element : support_points_map) + support_points_vec.push_back(element.second); + + particle_handler.insert_global_particles(support_points_vec, + global_bounding_boxes); + } + + + template + void + quadrature_points( + const Triangulation &particle_tria, + const Quadrature & quadrature, + // const std::vector> & particle_reference_locations, + const std::vector>> + & global_bounding_boxes, + ParticleHandler &particle_handler, + const Mapping & mapping) + { + const std::vector> &particle_reference_locations = + quadrature.get_points(); + std::vector> points_to_generate; + + // Loop through cells and gather gauss points + for (const auto &cell : particle_tria.active_cell_iterators()) + { + if (cell->is_locally_owned()) + { + for (const auto &reference_location : + particle_reference_locations) + { + const Point position_real = + mapping.transform_unit_to_real_cell(cell, + reference_location); + points_to_generate.push_back(position_real); + } + } + } + particle_handler.insert_global_particles(points_to_generate, + global_bounding_boxes); + } } // namespace Generators } // namespace Particles diff --git a/source/particles/generators.inst.in b/source/particles/generators.inst.in index bef4d7dc82..e1dfb5ec3d 100644 --- a/source/particles/generators.inst.in +++ b/source/particles/generators.inst.in @@ -51,6 +51,29 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS) &particle_handler, const Mapping &mapping, const unsigned int random_number_seed); + + template void + dof_support_points( + const DoFHandler + &particle_dof_handler, + const std::vector>> + &global_bounding_boxes, + ParticleHandler + &particle_handler, + const Mapping &mapping, + const ComponentMask &components); + + template void + quadrature_points( + const Triangulation + & particle_tria, + const Quadrature &quadrature, + const std::vector>> + &global_bounding_boxes, + ParticleHandler + &particle_handler, + const Mapping &mapping); + \} \} #endif diff --git a/tests/particles/generators_07.cc b/tests/particles/generators_07.cc new file mode 100644 index 0000000000..018dd3523e --- /dev/null +++ b/tests/particles/generators_07.cc @@ -0,0 +1,124 @@ +// --------------------------------------------------------------------- +// +// 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. +// +// --------------------------------------------------------------------- + +// Insert particles within an hyper_cube triangulation using the degrees of +// freedom associated with a non-matching hyper_cube triangulation and then +// check if the particles are correctly positioned + +#include +#include + +#include + +#include + +#include +#include +#include +#include + +#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); + + DoFHandler dof_handler(tr); + const FE_Nothing fe_nothing; + dof_handler.distribute_dofs(fe_nothing); + + const MappingQ mapping(1); + + Particles::ParticleHandler particle_handler(tr, mapping); + + parallel::distributed::Triangulation particles_tr( + MPI_COMM_WORLD); + GridGenerator::hyper_cube(particles_tr, 0.1, 0.9); + + // Generate the necessary bounding boxes for the generator + const auto my_bounding_box = GridTools::compute_mesh_predicate_bounding_box( + tr, IteratorFilters::LocallyOwnedCell()); + const auto global_bounding_boxes = + Utilities::MPI::all_gather(MPI_COMM_WORLD, my_bounding_box); + + DoFHandler particles_dof_handler(particles_tr); + const FE_Q particles_fe(1); + particles_dof_handler.distribute_dofs(particles_fe); + + + Particles::Generators::dof_support_points(particles_dof_handler, + global_bounding_boxes, + particle_handler); + + + { + deallog << "Locally owned active cells: " + << tr.n_locally_owned_active_cells() << std::endl; + + deallog << "Global particles: " << particle_handler.n_global_particles() + << std::endl; + + for (const auto &cell : tr.active_cell_iterators()) + { + if (cell->is_locally_owned()) + { + deallog << "Cell " << cell << " has " + << particle_handler.n_particles_in_cell(cell) + << " particles." << std::endl; + } + } + + for (const auto &particle : particle_handler) + { + deallog << "Particle index " << particle.get_id() << " is in cell " + << particle.get_surrounding_cell(tr) << std::endl; + deallog << "Particle location: " << particle.get_location() + << std::endl; + } + } + deallog << "OK" << std::endl; +} + + + +int +main(int argc, char *argv[]) +{ + Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1); + + MPILogInitAll init; + + 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/generators_07.with_p4est=true.mpirun=1.output b/tests/particles/generators_07.with_p4est=true.mpirun=1.output new file mode 100644 index 0000000000..7e79cffbc0 --- /dev/null +++ b/tests/particles/generators_07.with_p4est=true.mpirun=1.output @@ -0,0 +1,138 @@ + +DEAL:0:2d/2d::Locally owned active cells: 16 +DEAL:0:2d/2d::Global particles: 4 +DEAL:0:2d/2d::Cell 2.0 has 1 particles. +DEAL:0:2d/2d::Cell 2.1 has 0 particles. +DEAL:0:2d/2d::Cell 2.2 has 0 particles. +DEAL:0:2d/2d::Cell 2.3 has 0 particles. +DEAL:0:2d/2d::Cell 2.4 has 0 particles. +DEAL:0:2d/2d::Cell 2.5 has 1 particles. +DEAL:0:2d/2d::Cell 2.6 has 0 particles. +DEAL:0:2d/2d::Cell 2.7 has 0 particles. +DEAL:0:2d/2d::Cell 2.8 has 0 particles. +DEAL:0:2d/2d::Cell 2.9 has 0 particles. +DEAL:0:2d/2d::Cell 2.10 has 1 particles. +DEAL:0:2d/2d::Cell 2.11 has 0 particles. +DEAL:0:2d/2d::Cell 2.12 has 0 particles. +DEAL:0:2d/2d::Cell 2.13 has 0 particles. +DEAL:0:2d/2d::Cell 2.14 has 0 particles. +DEAL:0:2d/2d::Cell 2.15 has 1 particles. +DEAL:0:2d/2d::Particle index 0 is in cell 2.0 +DEAL:0:2d/2d::Particle location: 0.100000 0.100000 +DEAL:0:2d/2d::Particle index 1 is in cell 2.5 +DEAL:0:2d/2d::Particle location: 0.900000 0.100000 +DEAL:0:2d/2d::Particle index 2 is in cell 2.10 +DEAL:0:2d/2d::Particle location: 0.100000 0.900000 +DEAL:0:2d/2d::Particle index 3 is in cell 2.15 +DEAL:0:2d/2d::Particle location: 0.900000 0.900000 +DEAL:0:2d/2d::OK +DEAL:0:2d/3d::Locally owned active cells: 16 +DEAL:0:2d/3d::Global particles: 4 +DEAL:0:2d/3d::Cell 2.0 has 1 particles. +DEAL:0:2d/3d::Cell 2.1 has 0 particles. +DEAL:0:2d/3d::Cell 2.2 has 0 particles. +DEAL:0:2d/3d::Cell 2.3 has 0 particles. +DEAL:0:2d/3d::Cell 2.4 has 0 particles. +DEAL:0:2d/3d::Cell 2.5 has 1 particles. +DEAL:0:2d/3d::Cell 2.6 has 0 particles. +DEAL:0:2d/3d::Cell 2.7 has 0 particles. +DEAL:0:2d/3d::Cell 2.8 has 0 particles. +DEAL:0:2d/3d::Cell 2.9 has 0 particles. +DEAL:0:2d/3d::Cell 2.10 has 1 particles. +DEAL:0:2d/3d::Cell 2.11 has 0 particles. +DEAL:0:2d/3d::Cell 2.12 has 0 particles. +DEAL:0:2d/3d::Cell 2.13 has 0 particles. +DEAL:0:2d/3d::Cell 2.14 has 0 particles. +DEAL:0:2d/3d::Cell 2.15 has 1 particles. +DEAL:0:2d/3d::Particle index 0 is in cell 2.0 +DEAL:0:2d/3d::Particle location: 0.100000 0.100000 0.00000 +DEAL:0:2d/3d::Particle index 1 is in cell 2.5 +DEAL:0:2d/3d::Particle location: 0.900000 0.100000 0.00000 +DEAL:0:2d/3d::Particle index 2 is in cell 2.10 +DEAL:0:2d/3d::Particle location: 0.100000 0.900000 0.00000 +DEAL:0:2d/3d::Particle index 3 is in cell 2.15 +DEAL:0:2d/3d::Particle location: 0.900000 0.900000 0.00000 +DEAL:0:2d/3d::OK +DEAL:0:3d/3d::Locally owned active cells: 64 +DEAL:0:3d/3d::Global particles: 8 +DEAL:0:3d/3d::Cell 2.0 has 1 particles. +DEAL:0:3d/3d::Cell 2.1 has 0 particles. +DEAL:0:3d/3d::Cell 2.2 has 0 particles. +DEAL:0:3d/3d::Cell 2.3 has 0 particles. +DEAL:0:3d/3d::Cell 2.4 has 0 particles. +DEAL:0:3d/3d::Cell 2.5 has 0 particles. +DEAL:0:3d/3d::Cell 2.6 has 0 particles. +DEAL:0:3d/3d::Cell 2.7 has 0 particles. +DEAL:0:3d/3d::Cell 2.8 has 0 particles. +DEAL:0:3d/3d::Cell 2.9 has 1 particles. +DEAL:0:3d/3d::Cell 2.10 has 0 particles. +DEAL:0:3d/3d::Cell 2.11 has 0 particles. +DEAL:0:3d/3d::Cell 2.12 has 0 particles. +DEAL:0:3d/3d::Cell 2.13 has 0 particles. +DEAL:0:3d/3d::Cell 2.14 has 0 particles. +DEAL:0:3d/3d::Cell 2.15 has 0 particles. +DEAL:0:3d/3d::Cell 2.16 has 0 particles. +DEAL:0:3d/3d::Cell 2.17 has 0 particles. +DEAL:0:3d/3d::Cell 2.18 has 1 particles. +DEAL:0:3d/3d::Cell 2.19 has 0 particles. +DEAL:0:3d/3d::Cell 2.20 has 0 particles. +DEAL:0:3d/3d::Cell 2.21 has 0 particles. +DEAL:0:3d/3d::Cell 2.22 has 0 particles. +DEAL:0:3d/3d::Cell 2.23 has 0 particles. +DEAL:0:3d/3d::Cell 2.24 has 0 particles. +DEAL:0:3d/3d::Cell 2.25 has 0 particles. +DEAL:0:3d/3d::Cell 2.26 has 0 particles. +DEAL:0:3d/3d::Cell 2.27 has 1 particles. +DEAL:0:3d/3d::Cell 2.28 has 0 particles. +DEAL:0:3d/3d::Cell 2.29 has 0 particles. +DEAL:0:3d/3d::Cell 2.30 has 0 particles. +DEAL:0:3d/3d::Cell 2.31 has 0 particles. +DEAL:0:3d/3d::Cell 2.32 has 0 particles. +DEAL:0:3d/3d::Cell 2.33 has 0 particles. +DEAL:0:3d/3d::Cell 2.34 has 0 particles. +DEAL:0:3d/3d::Cell 2.35 has 0 particles. +DEAL:0:3d/3d::Cell 2.36 has 1 particles. +DEAL:0:3d/3d::Cell 2.37 has 0 particles. +DEAL:0:3d/3d::Cell 2.38 has 0 particles. +DEAL:0:3d/3d::Cell 2.39 has 0 particles. +DEAL:0:3d/3d::Cell 2.40 has 0 particles. +DEAL:0:3d/3d::Cell 2.41 has 0 particles. +DEAL:0:3d/3d::Cell 2.42 has 0 particles. +DEAL:0:3d/3d::Cell 2.43 has 0 particles. +DEAL:0:3d/3d::Cell 2.44 has 0 particles. +DEAL:0:3d/3d::Cell 2.45 has 1 particles. +DEAL:0:3d/3d::Cell 2.46 has 0 particles. +DEAL:0:3d/3d::Cell 2.47 has 0 particles. +DEAL:0:3d/3d::Cell 2.48 has 0 particles. +DEAL:0:3d/3d::Cell 2.49 has 0 particles. +DEAL:0:3d/3d::Cell 2.50 has 0 particles. +DEAL:0:3d/3d::Cell 2.51 has 0 particles. +DEAL:0:3d/3d::Cell 2.52 has 0 particles. +DEAL:0:3d/3d::Cell 2.53 has 0 particles. +DEAL:0:3d/3d::Cell 2.54 has 1 particles. +DEAL:0:3d/3d::Cell 2.55 has 0 particles. +DEAL:0:3d/3d::Cell 2.56 has 0 particles. +DEAL:0:3d/3d::Cell 2.57 has 0 particles. +DEAL:0:3d/3d::Cell 2.58 has 0 particles. +DEAL:0:3d/3d::Cell 2.59 has 0 particles. +DEAL:0:3d/3d::Cell 2.60 has 0 particles. +DEAL:0:3d/3d::Cell 2.61 has 0 particles. +DEAL:0:3d/3d::Cell 2.62 has 0 particles. +DEAL:0:3d/3d::Cell 2.63 has 1 particles. +DEAL:0:3d/3d::Particle index 0 is in cell 2.0 +DEAL:0:3d/3d::Particle location: 0.100000 0.100000 0.100000 +DEAL:0:3d/3d::Particle index 1 is in cell 2.9 +DEAL:0:3d/3d::Particle location: 0.900000 0.100000 0.100000 +DEAL:0:3d/3d::Particle index 2 is in cell 2.18 +DEAL:0:3d/3d::Particle location: 0.100000 0.900000 0.100000 +DEAL:0:3d/3d::Particle index 3 is in cell 2.27 +DEAL:0:3d/3d::Particle location: 0.900000 0.900000 0.100000 +DEAL:0:3d/3d::Particle index 4 is in cell 2.36 +DEAL:0:3d/3d::Particle location: 0.100000 0.100000 0.900000 +DEAL:0:3d/3d::Particle index 5 is in cell 2.45 +DEAL:0:3d/3d::Particle location: 0.900000 0.100000 0.900000 +DEAL:0:3d/3d::Particle index 6 is in cell 2.54 +DEAL:0:3d/3d::Particle location: 0.100000 0.900000 0.900000 +DEAL:0:3d/3d::Particle index 7 is in cell 2.63 +DEAL:0:3d/3d::Particle location: 0.900000 0.900000 0.900000 +DEAL:0:3d/3d::OK diff --git a/tests/particles/generators_07.with_p4est=true.mpirun=3.output b/tests/particles/generators_07.with_p4est=true.mpirun=3.output new file mode 100644 index 0000000000..b652c2d23d --- /dev/null +++ b/tests/particles/generators_07.with_p4est=true.mpirun=3.output @@ -0,0 +1,160 @@ + +DEAL:0:2d/2d::Locally owned active cells: 4 +DEAL:0:2d/2d::Global particles: 4 +DEAL:0:2d/2d::Cell 2.0 has 1 particles. +DEAL:0:2d/2d::Cell 2.1 has 0 particles. +DEAL:0:2d/2d::Cell 2.2 has 0 particles. +DEAL:0:2d/2d::Cell 2.3 has 0 particles. +DEAL:0:2d/2d::Particle index 0 is in cell 2.0 +DEAL:0:2d/2d::Particle location: 0.100000 0.100000 +DEAL:0:2d/2d::OK +DEAL:0:2d/3d::Locally owned active cells: 4 +DEAL:0:2d/3d::Global particles: 4 +DEAL:0:2d/3d::Cell 2.0 has 1 particles. +DEAL:0:2d/3d::Cell 2.1 has 0 particles. +DEAL:0:2d/3d::Cell 2.2 has 0 particles. +DEAL:0:2d/3d::Cell 2.3 has 0 particles. +DEAL:0:2d/3d::Particle index 0 is in cell 2.0 +DEAL:0:2d/3d::Particle location: 0.100000 0.100000 0.00000 +DEAL:0:2d/3d::OK +DEAL:0:3d/3d::Locally owned active cells: 24 +DEAL:0:3d/3d::Global particles: 8 +DEAL:0:3d/3d::Cell 2.0 has 1 particles. +DEAL:0:3d/3d::Cell 2.1 has 0 particles. +DEAL:0:3d/3d::Cell 2.2 has 0 particles. +DEAL:0:3d/3d::Cell 2.3 has 0 particles. +DEAL:0:3d/3d::Cell 2.4 has 0 particles. +DEAL:0:3d/3d::Cell 2.5 has 0 particles. +DEAL:0:3d/3d::Cell 2.6 has 0 particles. +DEAL:0:3d/3d::Cell 2.7 has 0 particles. +DEAL:0:3d/3d::Cell 2.8 has 0 particles. +DEAL:0:3d/3d::Cell 2.9 has 1 particles. +DEAL:0:3d/3d::Cell 2.10 has 0 particles. +DEAL:0:3d/3d::Cell 2.11 has 0 particles. +DEAL:0:3d/3d::Cell 2.12 has 0 particles. +DEAL:0:3d/3d::Cell 2.13 has 0 particles. +DEAL:0:3d/3d::Cell 2.14 has 0 particles. +DEAL:0:3d/3d::Cell 2.15 has 0 particles. +DEAL:0:3d/3d::Cell 2.16 has 0 particles. +DEAL:0:3d/3d::Cell 2.17 has 0 particles. +DEAL:0:3d/3d::Cell 2.18 has 1 particles. +DEAL:0:3d/3d::Cell 2.19 has 0 particles. +DEAL:0:3d/3d::Cell 2.20 has 0 particles. +DEAL:0:3d/3d::Cell 2.21 has 0 particles. +DEAL:0:3d/3d::Cell 2.22 has 0 particles. +DEAL:0:3d/3d::Cell 2.23 has 0 particles. +DEAL:0:3d/3d::Particle index 0 is in cell 2.0 +DEAL:0:3d/3d::Particle location: 0.100000 0.100000 0.100000 +DEAL:0:3d/3d::Particle index 1 is in cell 2.9 +DEAL:0:3d/3d::Particle location: 0.900000 0.100000 0.100000 +DEAL:0:3d/3d::Particle index 2 is in cell 2.18 +DEAL:0:3d/3d::Particle location: 0.100000 0.900000 0.100000 +DEAL:0:3d/3d::OK + +DEAL:1:2d/2d::Locally owned active cells: 8 +DEAL:1:2d/2d::Global particles: 4 +DEAL:1:2d/2d::Cell 2.4 has 0 particles. +DEAL:1:2d/2d::Cell 2.5 has 1 particles. +DEAL:1:2d/2d::Cell 2.6 has 0 particles. +DEAL:1:2d/2d::Cell 2.7 has 0 particles. +DEAL:1:2d/2d::Cell 2.8 has 0 particles. +DEAL:1:2d/2d::Cell 2.9 has 0 particles. +DEAL:1:2d/2d::Cell 2.10 has 1 particles. +DEAL:1:2d/2d::Cell 2.11 has 0 particles. +DEAL:1:2d/2d::Particle index 1 is in cell 2.5 +DEAL:1:2d/2d::Particle location: 0.900000 0.100000 +DEAL:1:2d/2d::Particle index 2 is in cell 2.10 +DEAL:1:2d/2d::Particle location: 0.100000 0.900000 +DEAL:1:2d/2d::OK +DEAL:1:2d/3d::Locally owned active cells: 8 +DEAL:1:2d/3d::Global particles: 4 +DEAL:1:2d/3d::Cell 2.4 has 0 particles. +DEAL:1:2d/3d::Cell 2.5 has 1 particles. +DEAL:1:2d/3d::Cell 2.6 has 0 particles. +DEAL:1:2d/3d::Cell 2.7 has 0 particles. +DEAL:1:2d/3d::Cell 2.8 has 0 particles. +DEAL:1:2d/3d::Cell 2.9 has 0 particles. +DEAL:1:2d/3d::Cell 2.10 has 1 particles. +DEAL:1:2d/3d::Cell 2.11 has 0 particles. +DEAL:1:2d/3d::Particle index 1 is in cell 2.5 +DEAL:1:2d/3d::Particle location: 0.900000 0.100000 0.00000 +DEAL:1:2d/3d::Particle index 2 is in cell 2.10 +DEAL:1:2d/3d::Particle location: 0.100000 0.900000 0.00000 +DEAL:1:2d/3d::OK +DEAL:1:3d/3d::Locally owned active cells: 16 +DEAL:1:3d/3d::Global particles: 8 +DEAL:1:3d/3d::Cell 2.24 has 0 particles. +DEAL:1:3d/3d::Cell 2.25 has 0 particles. +DEAL:1:3d/3d::Cell 2.26 has 0 particles. +DEAL:1:3d/3d::Cell 2.27 has 1 particles. +DEAL:1:3d/3d::Cell 2.28 has 0 particles. +DEAL:1:3d/3d::Cell 2.29 has 0 particles. +DEAL:1:3d/3d::Cell 2.30 has 0 particles. +DEAL:1:3d/3d::Cell 2.31 has 0 particles. +DEAL:1:3d/3d::Cell 2.32 has 0 particles. +DEAL:1:3d/3d::Cell 2.33 has 0 particles. +DEAL:1:3d/3d::Cell 2.34 has 0 particles. +DEAL:1:3d/3d::Cell 2.35 has 0 particles. +DEAL:1:3d/3d::Cell 2.36 has 1 particles. +DEAL:1:3d/3d::Cell 2.37 has 0 particles. +DEAL:1:3d/3d::Cell 2.38 has 0 particles. +DEAL:1:3d/3d::Cell 2.39 has 0 particles. +DEAL:1:3d/3d::Particle index 3 is in cell 2.27 +DEAL:1:3d/3d::Particle location: 0.900000 0.900000 0.100000 +DEAL:1:3d/3d::Particle index 4 is in cell 2.36 +DEAL:1:3d/3d::Particle location: 0.100000 0.100000 0.900000 +DEAL:1:3d/3d::OK + + +DEAL:2:2d/2d::Locally owned active cells: 4 +DEAL:2:2d/2d::Global particles: 4 +DEAL:2:2d/2d::Cell 2.12 has 0 particles. +DEAL:2:2d/2d::Cell 2.13 has 0 particles. +DEAL:2:2d/2d::Cell 2.14 has 0 particles. +DEAL:2:2d/2d::Cell 2.15 has 1 particles. +DEAL:2:2d/2d::Particle index 3 is in cell 2.15 +DEAL:2:2d/2d::Particle location: 0.900000 0.900000 +DEAL:2:2d/2d::OK +DEAL:2:2d/3d::Locally owned active cells: 4 +DEAL:2:2d/3d::Global particles: 4 +DEAL:2:2d/3d::Cell 2.12 has 0 particles. +DEAL:2:2d/3d::Cell 2.13 has 0 particles. +DEAL:2:2d/3d::Cell 2.14 has 0 particles. +DEAL:2:2d/3d::Cell 2.15 has 1 particles. +DEAL:2:2d/3d::Particle index 3 is in cell 2.15 +DEAL:2:2d/3d::Particle location: 0.900000 0.900000 0.00000 +DEAL:2:2d/3d::OK +DEAL:2:3d/3d::Locally owned active cells: 24 +DEAL:2:3d/3d::Global particles: 8 +DEAL:2:3d/3d::Cell 2.40 has 0 particles. +DEAL:2:3d/3d::Cell 2.41 has 0 particles. +DEAL:2:3d/3d::Cell 2.42 has 0 particles. +DEAL:2:3d/3d::Cell 2.43 has 0 particles. +DEAL:2:3d/3d::Cell 2.44 has 0 particles. +DEAL:2:3d/3d::Cell 2.45 has 1 particles. +DEAL:2:3d/3d::Cell 2.46 has 0 particles. +DEAL:2:3d/3d::Cell 2.47 has 0 particles. +DEAL:2:3d/3d::Cell 2.48 has 0 particles. +DEAL:2:3d/3d::Cell 2.49 has 0 particles. +DEAL:2:3d/3d::Cell 2.50 has 0 particles. +DEAL:2:3d/3d::Cell 2.51 has 0 particles. +DEAL:2:3d/3d::Cell 2.52 has 0 particles. +DEAL:2:3d/3d::Cell 2.53 has 0 particles. +DEAL:2:3d/3d::Cell 2.54 has 1 particles. +DEAL:2:3d/3d::Cell 2.55 has 0 particles. +DEAL:2:3d/3d::Cell 2.56 has 0 particles. +DEAL:2:3d/3d::Cell 2.57 has 0 particles. +DEAL:2:3d/3d::Cell 2.58 has 0 particles. +DEAL:2:3d/3d::Cell 2.59 has 0 particles. +DEAL:2:3d/3d::Cell 2.60 has 0 particles. +DEAL:2:3d/3d::Cell 2.61 has 0 particles. +DEAL:2:3d/3d::Cell 2.62 has 0 particles. +DEAL:2:3d/3d::Cell 2.63 has 1 particles. +DEAL:2:3d/3d::Particle index 5 is in cell 2.45 +DEAL:2:3d/3d::Particle location: 0.900000 0.100000 0.900000 +DEAL:2:3d/3d::Particle index 6 is in cell 2.54 +DEAL:2:3d/3d::Particle location: 0.100000 0.900000 0.900000 +DEAL:2:3d/3d::Particle index 7 is in cell 2.63 +DEAL:2:3d/3d::Particle location: 0.900000 0.900000 0.900000 +DEAL:2:3d/3d::OK + diff --git a/tests/particles/generators_08.cc b/tests/particles/generators_08.cc new file mode 100644 index 0000000000..8f037638d8 --- /dev/null +++ b/tests/particles/generators_08.cc @@ -0,0 +1,121 @@ +// --------------------------------------------------------------------- +// +// 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. +// +// --------------------------------------------------------------------- + +// Insert particles within an hyper_cube triangulation using a Q1 quadrature +// defined on a non-matching hyper_cube triangulation and then check if the +// particles are correctly positioned + + +#include + +#include + +#include + +#include +#include +#include +#include + +#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); + + DoFHandler dof_handler(tr); + const FE_Nothing fe_nothing; + dof_handler.distribute_dofs(fe_nothing); + + const MappingQ mapping(1); + + Particles::ParticleHandler particle_handler(tr, mapping); + + parallel::distributed::Triangulation particles_tr( + MPI_COMM_WORLD); + GridGenerator::hyper_cube(particles_tr, 0.1, 0.9); + + const QGauss quadrature(2); + + // Generate the necessary bounding boxes for the generator + const auto my_bounding_box = GridTools::compute_mesh_predicate_bounding_box( + tr, IteratorFilters::LocallyOwnedCell()); + const auto global_bounding_boxes = + Utilities::MPI::all_gather(MPI_COMM_WORLD, my_bounding_box); + + Particles::Generators::quadrature_points(particles_tr, + quadrature, + global_bounding_boxes, + particle_handler); + + { + deallog << "Locally owned active cells: " + << tr.n_locally_owned_active_cells() << std::endl; + + deallog << "Global particles: " << particle_handler.n_global_particles() + << std::endl; + + for (const auto &cell : tr.active_cell_iterators()) + { + if (cell->is_locally_owned()) + { + deallog << "Cell " << cell << " has " + << particle_handler.n_particles_in_cell(cell) + << " particles." << std::endl; + } + } + + for (const auto &particle : particle_handler) + { + deallog << "Particle index " << particle.get_id() << " is in cell " + << particle.get_surrounding_cell(tr) << std::endl; + deallog << "Particle location: " << particle.get_location() + << std::endl; + } + } + deallog << "OK" << std::endl; +} + + + +int +main(int argc, char *argv[]) +{ + Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1); + + MPILogInitAll init; + + 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/generators_08.with_p4est=true.mpirun=1.output b/tests/particles/generators_08.with_p4est=true.mpirun=1.output new file mode 100644 index 0000000000..d81a723c65 --- /dev/null +++ b/tests/particles/generators_08.with_p4est=true.mpirun=1.output @@ -0,0 +1,138 @@ + +DEAL:0:2d/2d::Locally owned active cells: 16 +DEAL:0:2d/2d::Global particles: 4 +DEAL:0:2d/2d::Cell 2.0 has 0 particles. +DEAL:0:2d/2d::Cell 2.1 has 0 particles. +DEAL:0:2d/2d::Cell 2.2 has 0 particles. +DEAL:0:2d/2d::Cell 2.3 has 1 particles. +DEAL:0:2d/2d::Cell 2.4 has 0 particles. +DEAL:0:2d/2d::Cell 2.5 has 0 particles. +DEAL:0:2d/2d::Cell 2.6 has 1 particles. +DEAL:0:2d/2d::Cell 2.7 has 0 particles. +DEAL:0:2d/2d::Cell 2.8 has 0 particles. +DEAL:0:2d/2d::Cell 2.9 has 1 particles. +DEAL:0:2d/2d::Cell 2.10 has 0 particles. +DEAL:0:2d/2d::Cell 2.11 has 0 particles. +DEAL:0:2d/2d::Cell 2.12 has 1 particles. +DEAL:0:2d/2d::Cell 2.13 has 0 particles. +DEAL:0:2d/2d::Cell 2.14 has 0 particles. +DEAL:0:2d/2d::Cell 2.15 has 0 particles. +DEAL:0:2d/2d::Particle index 0 is in cell 2.3 +DEAL:0:2d/2d::Particle location: 0.269060 0.269060 +DEAL:0:2d/2d::Particle index 1 is in cell 2.6 +DEAL:0:2d/2d::Particle location: 0.730940 0.269060 +DEAL:0:2d/2d::Particle index 2 is in cell 2.9 +DEAL:0:2d/2d::Particle location: 0.269060 0.730940 +DEAL:0:2d/2d::Particle index 3 is in cell 2.12 +DEAL:0:2d/2d::Particle location: 0.730940 0.730940 +DEAL:0:2d/2d::OK +DEAL:0:2d/3d::Locally owned active cells: 16 +DEAL:0:2d/3d::Global particles: 4 +DEAL:0:2d/3d::Cell 2.0 has 0 particles. +DEAL:0:2d/3d::Cell 2.1 has 0 particles. +DEAL:0:2d/3d::Cell 2.2 has 0 particles. +DEAL:0:2d/3d::Cell 2.3 has 1 particles. +DEAL:0:2d/3d::Cell 2.4 has 0 particles. +DEAL:0:2d/3d::Cell 2.5 has 0 particles. +DEAL:0:2d/3d::Cell 2.6 has 1 particles. +DEAL:0:2d/3d::Cell 2.7 has 0 particles. +DEAL:0:2d/3d::Cell 2.8 has 0 particles. +DEAL:0:2d/3d::Cell 2.9 has 1 particles. +DEAL:0:2d/3d::Cell 2.10 has 0 particles. +DEAL:0:2d/3d::Cell 2.11 has 0 particles. +DEAL:0:2d/3d::Cell 2.12 has 1 particles. +DEAL:0:2d/3d::Cell 2.13 has 0 particles. +DEAL:0:2d/3d::Cell 2.14 has 0 particles. +DEAL:0:2d/3d::Cell 2.15 has 0 particles. +DEAL:0:2d/3d::Particle index 0 is in cell 2.3 +DEAL:0:2d/3d::Particle location: 0.269060 0.269060 0.00000 +DEAL:0:2d/3d::Particle index 1 is in cell 2.6 +DEAL:0:2d/3d::Particle location: 0.730940 0.269060 0.00000 +DEAL:0:2d/3d::Particle index 2 is in cell 2.9 +DEAL:0:2d/3d::Particle location: 0.269060 0.730940 0.00000 +DEAL:0:2d/3d::Particle index 3 is in cell 2.12 +DEAL:0:2d/3d::Particle location: 0.730940 0.730940 0.00000 +DEAL:0:2d/3d::OK +DEAL:0:3d/3d::Locally owned active cells: 64 +DEAL:0:3d/3d::Global particles: 8 +DEAL:0:3d/3d::Cell 2.0 has 0 particles. +DEAL:0:3d/3d::Cell 2.1 has 0 particles. +DEAL:0:3d/3d::Cell 2.2 has 0 particles. +DEAL:0:3d/3d::Cell 2.3 has 0 particles. +DEAL:0:3d/3d::Cell 2.4 has 0 particles. +DEAL:0:3d/3d::Cell 2.5 has 0 particles. +DEAL:0:3d/3d::Cell 2.6 has 0 particles. +DEAL:0:3d/3d::Cell 2.7 has 1 particles. +DEAL:0:3d/3d::Cell 2.8 has 0 particles. +DEAL:0:3d/3d::Cell 2.9 has 0 particles. +DEAL:0:3d/3d::Cell 2.10 has 0 particles. +DEAL:0:3d/3d::Cell 2.11 has 0 particles. +DEAL:0:3d/3d::Cell 2.12 has 0 particles. +DEAL:0:3d/3d::Cell 2.13 has 0 particles. +DEAL:0:3d/3d::Cell 2.14 has 1 particles. +DEAL:0:3d/3d::Cell 2.15 has 0 particles. +DEAL:0:3d/3d::Cell 2.16 has 0 particles. +DEAL:0:3d/3d::Cell 2.17 has 0 particles. +DEAL:0:3d/3d::Cell 2.18 has 0 particles. +DEAL:0:3d/3d::Cell 2.19 has 0 particles. +DEAL:0:3d/3d::Cell 2.20 has 0 particles. +DEAL:0:3d/3d::Cell 2.21 has 1 particles. +DEAL:0:3d/3d::Cell 2.22 has 0 particles. +DEAL:0:3d/3d::Cell 2.23 has 0 particles. +DEAL:0:3d/3d::Cell 2.24 has 0 particles. +DEAL:0:3d/3d::Cell 2.25 has 0 particles. +DEAL:0:3d/3d::Cell 2.26 has 0 particles. +DEAL:0:3d/3d::Cell 2.27 has 0 particles. +DEAL:0:3d/3d::Cell 2.28 has 1 particles. +DEAL:0:3d/3d::Cell 2.29 has 0 particles. +DEAL:0:3d/3d::Cell 2.30 has 0 particles. +DEAL:0:3d/3d::Cell 2.31 has 0 particles. +DEAL:0:3d/3d::Cell 2.32 has 0 particles. +DEAL:0:3d/3d::Cell 2.33 has 0 particles. +DEAL:0:3d/3d::Cell 2.34 has 0 particles. +DEAL:0:3d/3d::Cell 2.35 has 1 particles. +DEAL:0:3d/3d::Cell 2.36 has 0 particles. +DEAL:0:3d/3d::Cell 2.37 has 0 particles. +DEAL:0:3d/3d::Cell 2.38 has 0 particles. +DEAL:0:3d/3d::Cell 2.39 has 0 particles. +DEAL:0:3d/3d::Cell 2.40 has 0 particles. +DEAL:0:3d/3d::Cell 2.41 has 0 particles. +DEAL:0:3d/3d::Cell 2.42 has 1 particles. +DEAL:0:3d/3d::Cell 2.43 has 0 particles. +DEAL:0:3d/3d::Cell 2.44 has 0 particles. +DEAL:0:3d/3d::Cell 2.45 has 0 particles. +DEAL:0:3d/3d::Cell 2.46 has 0 particles. +DEAL:0:3d/3d::Cell 2.47 has 0 particles. +DEAL:0:3d/3d::Cell 2.48 has 0 particles. +DEAL:0:3d/3d::Cell 2.49 has 1 particles. +DEAL:0:3d/3d::Cell 2.50 has 0 particles. +DEAL:0:3d/3d::Cell 2.51 has 0 particles. +DEAL:0:3d/3d::Cell 2.52 has 0 particles. +DEAL:0:3d/3d::Cell 2.53 has 0 particles. +DEAL:0:3d/3d::Cell 2.54 has 0 particles. +DEAL:0:3d/3d::Cell 2.55 has 0 particles. +DEAL:0:3d/3d::Cell 2.56 has 1 particles. +DEAL:0:3d/3d::Cell 2.57 has 0 particles. +DEAL:0:3d/3d::Cell 2.58 has 0 particles. +DEAL:0:3d/3d::Cell 2.59 has 0 particles. +DEAL:0:3d/3d::Cell 2.60 has 0 particles. +DEAL:0:3d/3d::Cell 2.61 has 0 particles. +DEAL:0:3d/3d::Cell 2.62 has 0 particles. +DEAL:0:3d/3d::Cell 2.63 has 0 particles. +DEAL:0:3d/3d::Particle index 0 is in cell 2.7 +DEAL:0:3d/3d::Particle location: 0.269060 0.269060 0.269060 +DEAL:0:3d/3d::Particle index 1 is in cell 2.14 +DEAL:0:3d/3d::Particle location: 0.730940 0.269060 0.269060 +DEAL:0:3d/3d::Particle index 2 is in cell 2.21 +DEAL:0:3d/3d::Particle location: 0.269060 0.730940 0.269060 +DEAL:0:3d/3d::Particle index 3 is in cell 2.28 +DEAL:0:3d/3d::Particle location: 0.730940 0.730940 0.269060 +DEAL:0:3d/3d::Particle index 4 is in cell 2.35 +DEAL:0:3d/3d::Particle location: 0.269060 0.269060 0.730940 +DEAL:0:3d/3d::Particle index 5 is in cell 2.42 +DEAL:0:3d/3d::Particle location: 0.730940 0.269060 0.730940 +DEAL:0:3d/3d::Particle index 6 is in cell 2.49 +DEAL:0:3d/3d::Particle location: 0.269060 0.730940 0.730940 +DEAL:0:3d/3d::Particle index 7 is in cell 2.56 +DEAL:0:3d/3d::Particle location: 0.730940 0.730940 0.730940 +DEAL:0:3d/3d::OK diff --git a/tests/particles/generators_08.with_p4est=true.mpirun=3.output b/tests/particles/generators_08.with_p4est=true.mpirun=3.output new file mode 100644 index 0000000000..770b7ca3e5 --- /dev/null +++ b/tests/particles/generators_08.with_p4est=true.mpirun=3.output @@ -0,0 +1,160 @@ + +DEAL:0:2d/2d::Locally owned active cells: 4 +DEAL:0:2d/2d::Global particles: 4 +DEAL:0:2d/2d::Cell 2.0 has 0 particles. +DEAL:0:2d/2d::Cell 2.1 has 0 particles. +DEAL:0:2d/2d::Cell 2.2 has 0 particles. +DEAL:0:2d/2d::Cell 2.3 has 1 particles. +DEAL:0:2d/2d::Particle index 0 is in cell 2.3 +DEAL:0:2d/2d::Particle location: 0.269060 0.269060 +DEAL:0:2d/2d::OK +DEAL:0:2d/3d::Locally owned active cells: 4 +DEAL:0:2d/3d::Global particles: 4 +DEAL:0:2d/3d::Cell 2.0 has 0 particles. +DEAL:0:2d/3d::Cell 2.1 has 0 particles. +DEAL:0:2d/3d::Cell 2.2 has 0 particles. +DEAL:0:2d/3d::Cell 2.3 has 1 particles. +DEAL:0:2d/3d::Particle index 0 is in cell 2.3 +DEAL:0:2d/3d::Particle location: 0.269060 0.269060 0.00000 +DEAL:0:2d/3d::OK +DEAL:0:3d/3d::Locally owned active cells: 24 +DEAL:0:3d/3d::Global particles: 8 +DEAL:0:3d/3d::Cell 2.0 has 0 particles. +DEAL:0:3d/3d::Cell 2.1 has 0 particles. +DEAL:0:3d/3d::Cell 2.2 has 0 particles. +DEAL:0:3d/3d::Cell 2.3 has 0 particles. +DEAL:0:3d/3d::Cell 2.4 has 0 particles. +DEAL:0:3d/3d::Cell 2.5 has 0 particles. +DEAL:0:3d/3d::Cell 2.6 has 0 particles. +DEAL:0:3d/3d::Cell 2.7 has 1 particles. +DEAL:0:3d/3d::Cell 2.8 has 0 particles. +DEAL:0:3d/3d::Cell 2.9 has 0 particles. +DEAL:0:3d/3d::Cell 2.10 has 0 particles. +DEAL:0:3d/3d::Cell 2.11 has 0 particles. +DEAL:0:3d/3d::Cell 2.12 has 0 particles. +DEAL:0:3d/3d::Cell 2.13 has 0 particles. +DEAL:0:3d/3d::Cell 2.14 has 1 particles. +DEAL:0:3d/3d::Cell 2.15 has 0 particles. +DEAL:0:3d/3d::Cell 2.16 has 0 particles. +DEAL:0:3d/3d::Cell 2.17 has 0 particles. +DEAL:0:3d/3d::Cell 2.18 has 0 particles. +DEAL:0:3d/3d::Cell 2.19 has 0 particles. +DEAL:0:3d/3d::Cell 2.20 has 0 particles. +DEAL:0:3d/3d::Cell 2.21 has 1 particles. +DEAL:0:3d/3d::Cell 2.22 has 0 particles. +DEAL:0:3d/3d::Cell 2.23 has 0 particles. +DEAL:0:3d/3d::Particle index 0 is in cell 2.7 +DEAL:0:3d/3d::Particle location: 0.269060 0.269060 0.269060 +DEAL:0:3d/3d::Particle index 1 is in cell 2.14 +DEAL:0:3d/3d::Particle location: 0.730940 0.269060 0.269060 +DEAL:0:3d/3d::Particle index 2 is in cell 2.21 +DEAL:0:3d/3d::Particle location: 0.269060 0.730940 0.269060 +DEAL:0:3d/3d::OK + +DEAL:1:2d/2d::Locally owned active cells: 8 +DEAL:1:2d/2d::Global particles: 4 +DEAL:1:2d/2d::Cell 2.4 has 0 particles. +DEAL:1:2d/2d::Cell 2.5 has 0 particles. +DEAL:1:2d/2d::Cell 2.6 has 1 particles. +DEAL:1:2d/2d::Cell 2.7 has 0 particles. +DEAL:1:2d/2d::Cell 2.8 has 0 particles. +DEAL:1:2d/2d::Cell 2.9 has 1 particles. +DEAL:1:2d/2d::Cell 2.10 has 0 particles. +DEAL:1:2d/2d::Cell 2.11 has 0 particles. +DEAL:1:2d/2d::Particle index 1 is in cell 2.6 +DEAL:1:2d/2d::Particle location: 0.730940 0.269060 +DEAL:1:2d/2d::Particle index 2 is in cell 2.9 +DEAL:1:2d/2d::Particle location: 0.269060 0.730940 +DEAL:1:2d/2d::OK +DEAL:1:2d/3d::Locally owned active cells: 8 +DEAL:1:2d/3d::Global particles: 4 +DEAL:1:2d/3d::Cell 2.4 has 0 particles. +DEAL:1:2d/3d::Cell 2.5 has 0 particles. +DEAL:1:2d/3d::Cell 2.6 has 1 particles. +DEAL:1:2d/3d::Cell 2.7 has 0 particles. +DEAL:1:2d/3d::Cell 2.8 has 0 particles. +DEAL:1:2d/3d::Cell 2.9 has 1 particles. +DEAL:1:2d/3d::Cell 2.10 has 0 particles. +DEAL:1:2d/3d::Cell 2.11 has 0 particles. +DEAL:1:2d/3d::Particle index 1 is in cell 2.6 +DEAL:1:2d/3d::Particle location: 0.730940 0.269060 0.00000 +DEAL:1:2d/3d::Particle index 2 is in cell 2.9 +DEAL:1:2d/3d::Particle location: 0.269060 0.730940 0.00000 +DEAL:1:2d/3d::OK +DEAL:1:3d/3d::Locally owned active cells: 16 +DEAL:1:3d/3d::Global particles: 8 +DEAL:1:3d/3d::Cell 2.24 has 0 particles. +DEAL:1:3d/3d::Cell 2.25 has 0 particles. +DEAL:1:3d/3d::Cell 2.26 has 0 particles. +DEAL:1:3d/3d::Cell 2.27 has 0 particles. +DEAL:1:3d/3d::Cell 2.28 has 1 particles. +DEAL:1:3d/3d::Cell 2.29 has 0 particles. +DEAL:1:3d/3d::Cell 2.30 has 0 particles. +DEAL:1:3d/3d::Cell 2.31 has 0 particles. +DEAL:1:3d/3d::Cell 2.32 has 0 particles. +DEAL:1:3d/3d::Cell 2.33 has 0 particles. +DEAL:1:3d/3d::Cell 2.34 has 0 particles. +DEAL:1:3d/3d::Cell 2.35 has 1 particles. +DEAL:1:3d/3d::Cell 2.36 has 0 particles. +DEAL:1:3d/3d::Cell 2.37 has 0 particles. +DEAL:1:3d/3d::Cell 2.38 has 0 particles. +DEAL:1:3d/3d::Cell 2.39 has 0 particles. +DEAL:1:3d/3d::Particle index 3 is in cell 2.28 +DEAL:1:3d/3d::Particle location: 0.730940 0.730940 0.269060 +DEAL:1:3d/3d::Particle index 4 is in cell 2.35 +DEAL:1:3d/3d::Particle location: 0.269060 0.269060 0.730940 +DEAL:1:3d/3d::OK + + +DEAL:2:2d/2d::Locally owned active cells: 4 +DEAL:2:2d/2d::Global particles: 4 +DEAL:2:2d/2d::Cell 2.12 has 1 particles. +DEAL:2:2d/2d::Cell 2.13 has 0 particles. +DEAL:2:2d/2d::Cell 2.14 has 0 particles. +DEAL:2:2d/2d::Cell 2.15 has 0 particles. +DEAL:2:2d/2d::Particle index 3 is in cell 2.12 +DEAL:2:2d/2d::Particle location: 0.730940 0.730940 +DEAL:2:2d/2d::OK +DEAL:2:2d/3d::Locally owned active cells: 4 +DEAL:2:2d/3d::Global particles: 4 +DEAL:2:2d/3d::Cell 2.12 has 1 particles. +DEAL:2:2d/3d::Cell 2.13 has 0 particles. +DEAL:2:2d/3d::Cell 2.14 has 0 particles. +DEAL:2:2d/3d::Cell 2.15 has 0 particles. +DEAL:2:2d/3d::Particle index 3 is in cell 2.12 +DEAL:2:2d/3d::Particle location: 0.730940 0.730940 0.00000 +DEAL:2:2d/3d::OK +DEAL:2:3d/3d::Locally owned active cells: 24 +DEAL:2:3d/3d::Global particles: 8 +DEAL:2:3d/3d::Cell 2.40 has 0 particles. +DEAL:2:3d/3d::Cell 2.41 has 0 particles. +DEAL:2:3d/3d::Cell 2.42 has 1 particles. +DEAL:2:3d/3d::Cell 2.43 has 0 particles. +DEAL:2:3d/3d::Cell 2.44 has 0 particles. +DEAL:2:3d/3d::Cell 2.45 has 0 particles. +DEAL:2:3d/3d::Cell 2.46 has 0 particles. +DEAL:2:3d/3d::Cell 2.47 has 0 particles. +DEAL:2:3d/3d::Cell 2.48 has 0 particles. +DEAL:2:3d/3d::Cell 2.49 has 1 particles. +DEAL:2:3d/3d::Cell 2.50 has 0 particles. +DEAL:2:3d/3d::Cell 2.51 has 0 particles. +DEAL:2:3d/3d::Cell 2.52 has 0 particles. +DEAL:2:3d/3d::Cell 2.53 has 0 particles. +DEAL:2:3d/3d::Cell 2.54 has 0 particles. +DEAL:2:3d/3d::Cell 2.55 has 0 particles. +DEAL:2:3d/3d::Cell 2.56 has 1 particles. +DEAL:2:3d/3d::Cell 2.57 has 0 particles. +DEAL:2:3d/3d::Cell 2.58 has 0 particles. +DEAL:2:3d/3d::Cell 2.59 has 0 particles. +DEAL:2:3d/3d::Cell 2.60 has 0 particles. +DEAL:2:3d/3d::Cell 2.61 has 0 particles. +DEAL:2:3d/3d::Cell 2.62 has 0 particles. +DEAL:2:3d/3d::Cell 2.63 has 0 particles. +DEAL:2:3d/3d::Particle index 5 is in cell 2.42 +DEAL:2:3d/3d::Particle location: 0.730940 0.269060 0.730940 +DEAL:2:3d/3d::Particle index 6 is in cell 2.49 +DEAL:2:3d/3d::Particle location: 0.269060 0.730940 0.730940 +DEAL:2:3d/3d::Particle index 7 is in cell 2.56 +DEAL:2:3d/3d::Particle location: 0.730940 0.730940 0.730940 +DEAL:2:3d/3d::OK +