From: Luca Heltai Date: Tue, 14 Apr 2020 22:53:53 +0000 (+0200) Subject: Fix set particle positions. X-Git-Tag: v9.2.0-rc1~231^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=8a1077f8ff9e53df8fc6e69ac4751fd3da3f05a4;p=dealii.git Fix set particle positions. --- diff --git a/include/deal.II/particles/particle_handler.h b/include/deal.II/particles/particle_handler.h index fa9822991c..24404e0cd2 100644 --- a/include/deal.II/particles/particle_handler.h +++ b/include/deal.II/particles/particle_handler.h @@ -334,7 +334,8 @@ namespace Particles * @authors Luca Heltai, Bruno Blais, 2019. */ template - void + typename boost::disable_if< + std::is_convertible *>>::type set_particle_positions(const VectorType &input_vector, const bool displace_particles = true); @@ -794,7 +795,8 @@ namespace Particles template template - void + typename boost::disable_if< + std::is_convertible *>>::type ParticleHandler::set_particle_positions( const VectorType &input_vector, const bool displace_particles) diff --git a/tests/particles/set_particle_positions_03.cc b/tests/particles/set_particle_positions_03.cc new file mode 100644 index 0000000000..4be165967d --- /dev/null +++ b/tests/particles/set_particle_positions_03.cc @@ -0,0 +1,120 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2017 - 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. +// +// --------------------------------------------------------------------- + +// Test insert_global_particles and set_particles_positions using a +// Function + +#include +#include + +#include + +#include + +#include +#include +#include + +#include + +#include "../tests.h" + +template +class Positions : public Function +{ +public: + Positions(unsigned int n_components) + : Function(n_components) + {} + + virtual double + value(const Point &p, unsigned int component = 0) const + { + if (component == 0) + return p[component] + 0.2; + else + return p[component]; + } + +private: +}; + +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); + + const unsigned int n_points = 2; + const unsigned int my_cpu = Utilities::MPI::this_mpi_process(MPI_COMM_WORLD); + const unsigned int n_cpus = Utilities::MPI::n_mpi_processes(MPI_COMM_WORLD); + Testing::srand(my_cpu + 1); + + // Create the bounding boxes for the global insertion + auto my_bounding_box = GridTools::compute_mesh_predicate_bounding_box( + tr, IteratorFilters::LocallyOwnedCell()); + + auto global_bounding_boxes = + Utilities::MPI::all_gather(MPI_COMM_WORLD, my_bounding_box); + + std::vector> points(n_points); + for (auto &p : points) + p = random_point(0.1, 0.2); + + auto cpu_to_index = + particle_handler.insert_global_particles(points, global_bounding_boxes); + + for (const auto &particle : particle_handler) + { + deallog << "Particle id : " << particle.get_id(); + deallog << " Particle location: " << particle.get_location() << std::endl; + } + + Positions new_positions(spacedim); + + particle_handler.set_particle_positions(new_positions, false); + + deallog << "Displacement" << std::endl; + for (const auto &particle : particle_handler) + { + deallog << "Particle id : " << particle.get_id(); + deallog << " Particle location: " << particle.get_location() << 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 << "---" << std::endl; + deallog.push("3d/3d"); + test<3, 3>(); + deallog.pop(); +} diff --git a/tests/particles/set_particle_positions_03.with_p4est=true.mpirun=1.output b/tests/particles/set_particle_positions_03.with_p4est=true.mpirun=1.output new file mode 100644 index 0000000000..15553841ad --- /dev/null +++ b/tests/particles/set_particle_positions_03.with_p4est=true.mpirun=1.output @@ -0,0 +1,12 @@ + +DEAL:0:2d/2d::Particle id : 0 Particle location: 0.184019 0.139438 +DEAL:0:2d/2d::Particle id : 1 Particle location: 0.178310 0.179844 +DEAL:0:2d/2d::Displacement +DEAL:0:2d/2d::Particle id : 0 Particle location: 0.384019 0.139438 +DEAL:0:2d/2d::Particle id : 1 Particle location: 0.378310 0.179844 +DEAL:0::--- +DEAL:0:3d/3d::Particle id : 0 Particle location: 0.184019 0.139438 0.178310 +DEAL:0:3d/3d::Particle id : 1 Particle location: 0.179844 0.191165 0.119755 +DEAL:0:3d/3d::Displacement +DEAL:0:3d/3d::Particle id : 0 Particle location: 0.384019 0.139438 0.178310 +DEAL:0:3d/3d::Particle id : 1 Particle location: 0.379844 0.191165 0.119755