From 80ec2e0b47070686089d87b6468246dc8acc6ee1 Mon Sep 17 00:00:00 2001 From: Rene Gassmoeller Date: Tue, 20 Feb 2024 03:00:15 -0500 Subject: [PATCH] Modify usage and include in tests --- examples/step-68/step-68.cc | 7 +--- include/deal.II/particles/particle.h | 15 -------- include/deal.II/particles/particle_accessor.h | 37 +++++++++---------- include/deal.II/particles/particle_handler.h | 13 ++++--- include/deal.II/particles/property_pool.h | 30 --------------- source/particles/particle_handler.cc | 8 ++-- tests/particles/step-68.cc | 7 +--- tests/performance/timing_step_68.cc | 3 +- tests/simplex/step-68.cc | 3 +- 9 files changed, 36 insertions(+), 87 deletions(-) diff --git a/examples/step-68/step-68.cc b/examples/step-68/step-68.cc index 149892bc5f..12214e9661 100644 --- a/examples/step-68/step-68.cc +++ b/examples/step-68/step-68.cc @@ -466,7 +466,7 @@ namespace Step68 { // We calculate the velocity of the particles using their current // location. - Point particle_location = particle.get_location(); + Point &particle_location = particle.get_location(); velocity.vector_value(particle_location, particle_velocity); // This updates the position of the particles and sets the old position @@ -474,8 +474,6 @@ namespace Step68 for (int d = 0; d < dim; ++d) particle_location[d] += particle_velocity[d] * dt; - particle.set_location(particle_location); - // We store the processor id (a scalar) and the particle velocity (a // vector) in the particle properties. In this example, this is done // purely for visualization purposes. @@ -534,11 +532,10 @@ namespace Step68 for (unsigned int particle_index = 0; particle != pic.end(); ++particle, ++particle_index) { - Point particle_location = particle->get_location(); + Point &particle_location = particle->get_location(); const Tensor<1, dim> &particle_velocity = evaluator.get_value(particle_index); particle_location += particle_velocity * dt; - particle->set_location(particle_location); // Again, we store the particle velocity and the processor id in the // particle properties for visualization purposes. diff --git a/include/deal.II/particles/particle.h b/include/deal.II/particles/particle.h index 600d671293..1d60e894a4 100644 --- a/include/deal.II/particles/particle.h +++ b/include/deal.II/particles/particle.h @@ -287,12 +287,6 @@ namespace Particles const Point & get_reference_location() const; - /** - * Return the reference location of this particle in its current cell. - */ - Point & - get_reference_location(); - /** * Return the ID number of this particle. The ID of a particle is intended * to be a property that is globally unique even in parallel computations @@ -578,15 +572,6 @@ namespace Particles - template - inline Point & - Particle::get_reference_location() - { - return property_pool->get_reference_location(property_pool_handle); - } - - - template inline types::particle_index Particle::get_id() const diff --git a/include/deal.II/particles/particle_accessor.h b/include/deal.II/particles/particle_accessor.h index 1bce933fe6..267d5c10cb 100644 --- a/include/deal.II/particles/particle_accessor.h +++ b/include/deal.II/particles/particle_accessor.h @@ -149,7 +149,7 @@ namespace Particles set_location(const Point &new_location); /** - * Get the location of this particle. + * Get read-access to the location of this particle. * * @return The location of this particle. */ @@ -157,7 +157,23 @@ namespace Particles get_location() const; /** - * Get the location of this particle. + * Get read- and write-access to the location of this particle. + * Note that changing the location does not check + * whether this is a valid location in the simulation domain. + * + * @note In parallel programs, the ParticleHandler class stores particles + * on both the locally owned cells, as well as on ghost cells. The + * particles on the latter are *copies* of particles owned on other + * processors, and should therefore be treated in the same way as + * ghost entries in + * @ref GlossGhostedVector "vectors with ghost elements" + * or + * @ref GlossGhostCell "ghost cells": + * In both cases, one should + * treat the ghost elements or cells as `const` objects that shouldn't + * be modified even if the objects allow for calls that modify + * properties. Rather, properties should only be modified on processors + * that actually *own* the particle. * * @return The location of this particle. */ @@ -193,12 +209,6 @@ namespace Particles const Point & get_reference_location() const; - /** - * Return the reference location of this particle in its current cell. - */ - Point & - get_reference_location(); - /** * Set the ID number of this particle. */ @@ -698,17 +708,6 @@ namespace Particles - template - inline Point & - ParticleAccessor::get_reference_location() - { - Assert(state() == IteratorState::valid, ExcInternalError()); - - return property_pool->get_reference_location(get_handle()); - } - - - template inline types::particle_index ParticleAccessor::get_id() const diff --git a/include/deal.II/particles/particle_handler.h b/include/deal.II/particles/particle_handler.h index 4a5c98d20d..672a116c37 100644 --- a/include/deal.II/particles/particle_handler.h +++ b/include/deal.II/particles/particle_handler.h @@ -1341,12 +1341,15 @@ namespace Particles get_next_free_particle_index() * spacedim); for (auto &p : *this) { - Point new_point(displace_particles ? p.get_location() : - Point()); + Point &position = p.get_location(); const auto id = p.get_id(); - for (unsigned int i = 0; i < spacedim; ++i) - new_point[i] += input_vector[id * spacedim + i]; - p.set_location(new_point); + + if (displace_particles) + for (unsigned int i = 0; i < spacedim; ++i) + position[i] += input_vector[id * spacedim + i]; + else + for (unsigned int i = 0; i < spacedim; ++i) + position[i] = input_vector[id * spacedim + i]; } sort_particles_into_subdomains_and_cells(); } diff --git a/include/deal.II/particles/property_pool.h b/include/deal.II/particles/property_pool.h index 7131149956..5e4643553d 100644 --- a/include/deal.II/particles/property_pool.h +++ b/include/deal.II/particles/property_pool.h @@ -177,13 +177,6 @@ namespace Particles const Point & get_reference_location(const Handle handle) const; - /** - * Return a writeable reference to the reference location of a particle - * identified by the given `handle`. - */ - Point & - get_reference_location(const Handle handle); - /** * Set the reference location of a particle identified by the given * `handle`. @@ -391,29 +384,6 @@ namespace Particles - template - inline Point & - PropertyPool::get_reference_location(const Handle handle) - { - const std::vector::size_type data_index = - (handle != invalid_handle) ? handle : 0; - - // Ideally we would need to assert that 'handle' has not been deallocated - // by searching through 'currently_available_handles'. However, this - // is expensive and this function is performance critical, so instead - // just check against the array range, and rely on the fact - // that handles are invalidated when handed over to - // deallocate_properties_array(). - Assert(data_index <= reference_locations.size() - 1, - ExcMessage("Invalid location handle. This can happen if the " - "handle was duplicated and then one copy was deallocated " - "before trying to access the properties.")); - - return reference_locations[data_index]; - } - - - template inline void PropertyPool::set_reference_location( diff --git a/source/particles/particle_handler.cc b/source/particles/particle_handler.cc index df6d147959..0d47946dc6 100644 --- a/source/particles/particle_handler.cc +++ b/source/particles/particle_handler.cc @@ -1127,10 +1127,11 @@ namespace Particles unsigned int i = 0; for (auto it = begin(); it != end(); ++it, ++i) { + Point &location = it->get_location(); if (displace_particles) - it->set_location(it->get_location() + new_positions[i]); + location += new_positions[i]; else - it->set_location(new_positions[i]); + location = new_positions[i]; } sort_particles_into_subdomains_and_cells(); } @@ -1150,7 +1151,7 @@ namespace Particles Vector new_position(spacedim); for (auto &particle : *this) { - Point particle_location = particle.get_location(); + Point &particle_location = particle.get_location(); function.vector_value(particle_location, new_position); if (displace_particles) for (unsigned int d = 0; d < spacedim; ++d) @@ -1158,7 +1159,6 @@ namespace Particles else for (unsigned int d = 0; d < spacedim; ++d) particle_location[d] = new_position[d]; - particle.set_location(particle_location); } sort_particles_into_subdomains_and_cells(); } diff --git a/tests/particles/step-68.cc b/tests/particles/step-68.cc index 572ad60767..24c1c41bd7 100644 --- a/tests/particles/step-68.cc +++ b/tests/particles/step-68.cc @@ -440,7 +440,7 @@ namespace Step68 { // We calculate the velocity of the particles using their current // location. - Point particle_location = particle.get_location(); + Point &particle_location = particle.get_location(); velocity.vector_value(particle_location, particle_velocity); // This updates the position of the particles and sets the old position @@ -448,8 +448,6 @@ namespace Step68 for (int d = 0; d < dim; ++d) particle_location[d] += particle_velocity[d] * dt; - particle.set_location(particle_location); - // We store the processor id (a scalar) and the particle velocity (a // vector) in the particle properties. In this example, this is done // purely for visualization purposes. @@ -513,10 +511,9 @@ namespace Step68 local_dof_values[j]; } - Point particle_location = particle->get_location(); + Point &particle_location = particle->get_location(); for (int d = 0; d < dim; ++d) particle_location[d] += particle_velocity[d] * dt; - p.set_location(particle_location); // Again, we store the particle velocity and the processor id in the // particle properties for visualization purposes. diff --git a/tests/performance/timing_step_68.cc b/tests/performance/timing_step_68.cc index 4150ab19f7..0caa517841 100644 --- a/tests/performance/timing_step_68.cc +++ b/tests/performance/timing_step_68.cc @@ -317,11 +317,10 @@ namespace Step68 for (unsigned int particle_index = 0; particle != pic.end(); ++particle, ++particle_index) { - Point particle_location = particle->get_location(); + Point &particle_location = particle->get_location(); const Tensor<1, dim> &particle_velocity = evaluator.get_value(particle_index); particle_location += particle_velocity * dt; - particle->set_location(particle_location); ArrayView properties = particle->get_properties(); for (int d = 0; d < dim; ++d) diff --git a/tests/simplex/step-68.cc b/tests/simplex/step-68.cc index ee645b6d26..524aac8696 100644 --- a/tests/simplex/step-68.cc +++ b/tests/simplex/step-68.cc @@ -377,10 +377,9 @@ namespace Step68 local_dof_values[j]; } - Point particle_location = particle->get_location(); + Point &particle_location = particle->get_location(); for (int d = 0; d < dim; ++d) particle_location[d] += particle_velocity[d] * dt; - p.set_location(particle_location); // Again, we store the particle velocity and the processor id in the // particle properties for visualization purposes. -- 2.39.5