]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Warn against the use of certain functions.
authorWolfgang Bangerth <bangerth@colostate.edu>
Mon, 29 Jun 2020 01:03:39 +0000 (19:03 -0600)
committerWolfgang Bangerth <bangerth@colostate.edu>
Mon, 29 Jun 2020 01:04:01 +0000 (19:04 -0600)
include/deal.II/particles/particle_handler.h

index a06a03b27c517bd41488ac13f41a26706b9fd9f0..4038da69bdcecd9018c0072519ee47c8a4fe94a2 100644 (file)
@@ -177,6 +177,26 @@ namespace Particles
 
     /**
      * Return the number of particles that live on the given cell.
+     *
+     * @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. Counting how many
+     *   steps this took will then give you the number you are looking for,
+     *   at a cost of ${\cal O)(\log N)$ when accumulated over all cells.
+     *   This is the approach used in step-70, for example. The approach is
+     *   also detailed in the "Possibilities for extensions section"
+     *   of step-19.
      */
     types::particle_index
     n_particles_in_cell(
@@ -190,6 +210,25 @@ 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(
@@ -202,6 +241,25 @@ 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(

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.