From 341ea1c70823bf67e21e14cb9b4466f7b205c226 Mon Sep 17 00:00:00 2001 From: Rene Gassmoeller Date: Fri, 11 Jun 2021 19:08:04 -0400 Subject: [PATCH] Update documentation. Fix missing function. --- examples/step-19/doc/results.dox | 71 ---------------- include/deal.II/particles/particle_accessor.h | 18 ++++ include/deal.II/particles/particle_handler.h | 85 ++++++++----------- source/particles/particle_handler.cc | 2 + 4 files changed, 57 insertions(+), 119 deletions(-) diff --git a/examples/step-19/doc/results.dox b/examples/step-19/doc/results.dox index 190dcdfaa3..448f15495f 100644 --- a/examples/step-19/doc/results.dox +++ b/examples/step-19/doc/results.dox @@ -159,77 +159,6 @@ electrode.

Possibilities for extensions

-

Avoiding a performance bottleneck with particles

- -The `assemble_system()`, `move_particles()`, and `update_timestep_size()` -functions all call Particles::ParticleHandler::particles_in_cell() and -Particles::ParticleHandler::n_particles_in_cell() that query information -about the particles located on the current cell. While this is convenient, -it's also inefficient. To understand why this is so, one needs to know -how particles are stored in Particles::ParticleHandler: namely, in a -data structure in which particles are ordered in some kind of linear -fashion sorted by the cell they are on. Consequently, in order to find -the particles associated with a given cell, these functions need to -search for the first (and possibly last) particle on a given cell -- -an effort that costs ${\cal O}(\log N)$ operations where $N$ is the -number of particles. But this is repeated on every cell; assuming that -for large computations, the number of cells and particles are roughly -proportional, the accumulated cost of these function calls is then -${\cal O}(N \log N)$ and consequently larger than the ${\cal O}(N)$ -cost that we should shoot for with all parts of a program. - -We can make this cheaper, though. First, instead of calling -Particles::ParticleHandler::n_particles_in_cell(), we might first call -Particles::ParticleHandler::particles_in_cell() and then compute the -number of particles on a cell by just computing the distance of the last -to the first particle on the current cell: -@code - const typename Particles::ParticleHandler::particle_iterator_range - particles_in_cell = particle_handler.particles_in_cell(cell); - const unsigned int - n_particles_in_cell = std::distance (particles_in_cell.begin(), - particles_in_cell.end()); -@endcode -The first of these calls is of course still ${\cal O}(\log N)$, -but at least the second call only takes a compute time proportional to -the number of particles on the current cell and so, when accumulated -over all cells, has a cost of ${\cal O}(N)$. - -But we can even get rid of the first of these calls with some proper algorithm -design. That's because particles are ordered in the same way as cells, and so -we can just walk them as we move along on the cells. The following outline -of an algorithm does this: -@code - auto begin_particle_on_cell = particle_handler.begin(); - for (const auto &cell : dof_handler.active_cell_iterators()) - { - unsigned int n_particles_on_cell = 0; - auto end_particle_on_cell = begin_particle_on_cell; - while (end_particle_on_cell->get_surrounding_cell(triangulation) - == cell) - { - ++n_particles_on_cell; - ++end_particle_on_cell; - } - - ...now operate on the range of particles from begin_particle_on_cell - to end_particle_on_cell, all of which are known to be on the current - cell...; - - // Move the begin iterator forward so that it points to the first - // particle on the next cell - begin_particle_on_cell = end_particle_on_cell; - } -@endcode - -In this code, we touch every cell exactly once and we never have to search -the big data structure for the first or last particle on each cell. As a -consequence, the algorithm costs a total of ${\cal O}(N)$ for a complete -sweep of all particles and all cells. - -It would not be very difficult to implement this scheme for all three of the -functions in this program that have this issue. -

More statistics about electrons

diff --git a/include/deal.II/particles/particle_accessor.h b/include/deal.II/particles/particle_accessor.h index 073f3a69c0..e035a3941e 100644 --- a/include/deal.II/particles/particle_accessor.h +++ b/include/deal.II/particles/particle_accessor.h @@ -193,6 +193,14 @@ namespace Particles const ArrayView get_properties() const; + /** + * @deprecated: This function is only kept for backward compatibility + * and has no meaning any more. ParticleAccessors always use the + * property pool of the owning particle handler. + */ + void + set_property_pool(PropertyPool &property_pool); + /** * Return the size in bytes this particle occupies if all of its data is * serialized (i.e. the number of bytes that is written by the write_data @@ -644,6 +652,16 @@ namespace Particles + template + inline void + ParticleAccessor::set_property_pool( + PropertyPool & /*new_property_pool*/) + { + /* nothing to do */ + } + + + template inline const typename Triangulation::cell_iterator & ParticleAccessor::get_surrounding_cell() const diff --git a/include/deal.II/particles/particle_handler.h b/include/deal.II/particles/particle_handler.h index 6e1175f46c..03edfd0ef0 100644 --- a/include/deal.II/particles/particle_handler.h +++ b/include/deal.II/particles/particle_handler.h @@ -51,9 +51,20 @@ namespace Particles * is designed in a similar way as the triangulation class. In particular, * we call particles in the domain of the local process local particles, * and particles that belong to neighbor processes and live in the ghost cells - * around the locally owned domain "ghost particles". + * around the locally owned domain "ghost particles". The class also includes + * functionality that is similar to the DoFHandler() class (it knows which + * particles live on which cells) and the SolutionTransfer() class (it know + * how to transfer particles between cells and subdomains). * - * This class is used in step-70. + * @note: While the class can be used in any kind of triangulation, transfer + * of particles during mesh refinement is currently only implemented for + * distributed triangulations. You can still use the class for serial + * triangulations, but you cannot change the mesh while particles + * exist inside the particle handler. + * + * For examples on how to use this class to track particles, store properties + * on particles, and let the properties on the particles influence the + * finite-element solution see step-19, step-68, and step-70. * * @ingroup Particle */ @@ -84,7 +95,7 @@ namespace Particles /** * Constructor that initializes the particle handler with - * a given triangulation and mapping. Since particles are stored in + * a given triangulation and mapping. Since particles are stored with * respect to their surrounding cells this information is necessary to * correctly organize the particle collection. * This constructor is equivalent to calling the default constructor and @@ -100,8 +111,8 @@ namespace Particles virtual ~ParticleHandler(); /** - * Initialize the particle handler. This function does not clear the - * internal data structures, it just sets the triangulation and the + * Initialize the particle handler. This function does clear the + * internal data structures, and sets the triangulation and the * mapping to be used. */ void @@ -133,7 +144,7 @@ namespace Particles copy_from(const ParticleHandler &particle_handler); /** - * Clear all particle related data. + * Clear all particle related data but keep the handler initialized. */ void clear(); @@ -158,13 +169,13 @@ namespace Particles update_cached_numbers(); /** - * Return an iterator to the first particle. + * Return an iterator to the first locally owned particle. */ particle_iterator begin() const; /** - * Return an iterator to the first particle. + * Return an iterator to the first locally owned particle. */ particle_iterator begin(); @@ -220,25 +231,6 @@ namespace Particles * * The number of elements in the returned range equals what the * n_particles_in_cell() function returns. - * - * @note While this function is used in step-19, it is not an efficient - * function to use if the number of particles is large. That is because - * to find the particles that are located in one cell costs - * ${\cal O}(\log N)$ where $N$ is the number of overall particles. Since - * you will likely do this for every cell, and assuming that the number - * of particles and the number of cells are roughly proportional, - * you end up with an ${\cal O}(N \log N)$ algorithm. A better approach - * is to use the fact that internally, particles are arranged in the - * order of the active cells they are in. In other words, if you iterate - * over all particles, you will encounter them in the same order as - * you walk over the active cells. You can exploit this by keeping an - * iterator to the first particle of the first cell, and when you move - * to the next cell, you increment the particle iterator as well until - * you find a particle located on that next cell. This is the approach - * used in step-70, for example, and has an overall cost of - * ${\cal O}(\log N)$ when accumulated over all cells. The approach is - * also detailed in the "Possibilities for extensions section" - * of step-19. */ particle_iterator_range particles_in_cell( @@ -258,18 +250,17 @@ namespace Particles const; /** - * Remove a particle pointed to by the iterator. Afterwards the iterator - * will point to one of the remaining particles in the cell. Note that - * particle iterators that point to other particles in the same cell as @p particle - * may not longer be valid after this function call. + * Remove a particle pointed to by the iterator. Note that @p particle + * and all iterators that point to other particles in the same cell + * as @p particle will be invalidated during this call. */ void remove_particle(const particle_iterator &particle); /** * Remove a vector of particles indicated by the particle iterators. - * The iterators and all other iterators are invalidated during the - * function call. + * The iterators and all other particle iterators are invalidated + * during the function call. */ void remove_particles(const std::vector &particles); @@ -595,7 +586,7 @@ namespace Particles &load_callback); /** - * Return the total number of particles that were managed by this class + * Return the total number of particles that are managed by this class * the last time the update_cached_numbers() function was called. * The actual number of particles may have changed since then if * particles have been added or removed. @@ -827,14 +818,14 @@ namespace Particles * properties. Since particles reference the property pool, the * latter has to be destroyed *after* the particles are destroyed. * This is achieved by making sure the `property_pool` member variable - * precedes the declaration of the `particles` and `ghost_particles` - * members. + * precedes the declaration of the `particles` + * member variable. */ std::unique_ptr> property_pool; /** * Set of particles currently living in the local domain including ghost - * cells, organized by the active cell of the cell they are in. + * cells, organized by the active cell index of the cell they are in. */ particle_container particles; @@ -933,10 +924,9 @@ namespace Particles * @param [in] particles_to_send All particles that should be sent and * their new subdomain_ids are in this map. * - * @param [in,out] received_particles Vector that stores all received - * particles. Note that it is not required nor checked that the list - * is empty, received particles are simply attached to the end of - * the vector. + * @param [in,out] received_particles Particle container that stores all received + * particles. Note that it is not required nor checked that the container + * is empty, received particles are simply inserted into the container. * * @param [in] new_cells_for_particles Optional vector of cell * iterators with the same structure as @p particles_to_send. If this @@ -968,12 +958,11 @@ namespace Particles const bool enable_cache = false); /** - * Transfer particles position and properties assuming that + * Transfer ghost particles position and properties assuming that * the particles have not changed cells. This routine uses the - * GhostParticlePartitioner as a caching structure to update the particles. + * GhostParticlePartitioner as a caching structure to know which particles + * are ghost to other processes, and where they need to be send. * It inherently assumes that particles cannot have changed cell. - * All updated particles will be appended to the - * @p received_particles container. * * @param [in] particles_to_send All particles for which information * should be sent and their new subdomain_ids are in this map. @@ -981,7 +970,7 @@ namespace Particles * @param [in,out] received_particles A map with all received * particles. Note that it is not required nor checked that the container * is empty, received particles are simply inserted into - * the map. + * the container. * */ void @@ -1012,8 +1001,8 @@ namespace Particles const typename Triangulation::CellStatus status) const; /** - * Called by listener functions after a refinement step. The local map - * of particles has to be read from the triangulation user_pointer. + * Called by listener functions after a refinement step to receive + * particles and insert them into the particle container. */ void load_particles( diff --git a/source/particles/particle_handler.cc b/source/particles/particle_handler.cc index 78ae3a3816..41de53f38d 100644 --- a/source/particles/particle_handler.cc +++ b/source/particles/particle_handler.cc @@ -116,6 +116,8 @@ namespace Particles const Mapping & new_mapping, const unsigned int n_properties) { + clear(); + triangulation = &new_triangulation; mapping = &new_mapping; -- 2.39.5