From 9c87e785a2104ffab316468dabe97db173edb190 Mon Sep 17 00:00:00 2001 From: Luca Heltai Date: Thu, 21 Nov 2019 19:52:06 +0100 Subject: [PATCH] Add simple output for particles. This class manages the DataOut of a Particle Handler It currently only supports witing the particle position and their ID Co-authored-by: Bruno Blais --- include/deal.II/particles/particle_handler.h | 73 ++++++++++++- tests/particles/particle_output_01.cc | 102 ++++++++++++++++++ ..._output_01.with_p4est=true.mpirun=1.output | 43 ++++++++ 3 files changed, 217 insertions(+), 1 deletion(-) create mode 100644 tests/particles/particle_output_01.cc create mode 100644 tests/particles/particle_output_01.with_p4est=true.mpirun=1.output diff --git a/include/deal.II/particles/particle_handler.h b/include/deal.II/particles/particle_handler.h index cac550d4e0..d11c9565c0 100644 --- a/include/deal.II/particles/particle_handler.h +++ b/include/deal.II/particles/particle_handler.h @@ -524,7 +524,78 @@ namespace Particles &data_range); }; - /* ---------------------- inline and template functions ------------------ */ + /** + * This class manages the DataOut of a Particle Handler + * It currently only supports witing the particle position + * and their ID + * + * @ingroup Particle + * + * @author : Bruno Blais, Luca Heltai 2019 + */ + + template + class ParticleOutput : public dealii::DataOutInterface<0, spacedim> + { + public: + ParticleOutput() = default; + + void + build_patches(const Particles::ParticleHandler &particles) + { + dataset_names.reserve(1); + dataset_names.emplace_back("id"); + patches.resize(particles.n_locally_owned_particles()); + + auto particle = particles.begin(); + for (unsigned int i = 0; particle != particles.end(); ++particle, ++i) + { + patches[i].vertices[0] = particle->get_location(); + patches[i].patch_index = i; + patches[i].n_subdivisions = 1; + patches[i].data.reinit(dataset_names.size(), 1); + + patches[i].data(0, 0) = particle->get_id(); + } + } + + + + ~ParticleOutput() = default; + + protected: + /** + * Implementation of the corresponding function of the base class. + */ + virtual const std::vector> & + get_patches() const + { + return patches; + } + + /** + * Implementation of the corresponding function of the base class. + */ + virtual std::vector + get_dataset_names() const + { + return dataset_names; + } + + private: + std::vector> patches; + + /** + * A list of field names for all data components stored in patches. + */ + std::vector dataset_names; + }; + + + /* ---------------------- inline and template functions ------------------ + */ + + template template diff --git a/tests/particles/particle_output_01.cc b/tests/particles/particle_output_01.cc new file mode 100644 index 0000000000..f5bac30a3a --- /dev/null +++ b/tests/particles/particle_output_01.cc @@ -0,0 +1,102 @@ +// --------------------------------------------------------------------- +// +// 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. +// +// --------------------------------------------------------------------- + + + +// Create a particle handler then generate an output +// Tests the Visualization class of the particle handler + +#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); + + // both processes create a particle handler, but only the first creates + // particles + Particles::ParticleHandler particle_handler(tr, mapping); + std::vector> position(2); + std::vector> reference_position(2); + + for (unsigned int i = 0; i < dim; ++i) + { + position[0](i) = 0.125; + position[1](i) = 0.525; + } + + Particles::Particle particle1(position[0], + reference_position[0], + 0); + Particles::Particle particle2(position[1], + reference_position[1], + 1); + + typename Triangulation::active_cell_iterator cell1(&tr, + 2, + 0); + typename Triangulation::active_cell_iterator cell2(&tr, + 2, + 0); + + particle_handler.insert_particle(particle1, cell1); + particle_handler.insert_particle(particle2, cell2); + + Particles::ParticleOutput particle_viz; + particle_viz.build_patches(particle_handler); + particle_viz.write_gnuplot(deallog.get_file_stream()); + + for (const auto &particle : particle_handler) + deallog << "Particle id " << particle.get_id() << " is in cell " + << particle.get_surrounding_cell(tr) << 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_output_01.with_p4est=true.mpirun=1.output b/tests/particles/particle_output_01.with_p4est=true.mpirun=1.output new file mode 100644 index 0000000000..5f6694ce5d --- /dev/null +++ b/tests/particles/particle_output_01.with_p4est=true.mpirun=1.output @@ -0,0 +1,43 @@ + +# This file was generated by the deal.II library. + + +# +# For a description of the GNUPLOT format see the GNUPLOT manual. +# +# +0.125000 0.125000 0.00000 + +0.525000 0.525000 1.00000 + +DEAL:0:2d/2d::Particle id 0 is in cell 2.0 +DEAL:0:2d/2d::Particle id 1 is in cell 2.0 +DEAL:0:2d/2d::OK +# This file was generated by the deal.II library. + + +# +# For a description of the GNUPLOT format see the GNUPLOT manual. +# +# +0.125000 0.125000 0.00000 0.00000 + +0.525000 0.525000 0.00000 1.00000 + +DEAL:0:2d/3d::Particle id 0 is in cell 2.0 +DEAL:0:2d/3d::Particle id 1 is in cell 2.0 +DEAL:0:2d/3d::OK +# This file was generated by the deal.II library. + + +# +# For a description of the GNUPLOT format see the GNUPLOT manual. +# +# +0.125000 0.125000 0.125000 0.00000 + +0.525000 0.525000 0.525000 1.00000 + +DEAL:0:3d/3d::Particle id 0 is in cell 2.0 +DEAL:0:3d/3d::Particle id 1 is in cell 2.0 +DEAL:0:3d/3d::OK -- 2.39.5