]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Better document CommunicationPatternBase and derived classes. 14740/head
authorWolfgang Bangerth <bangerth@colostate.edu>
Thu, 26 Jan 2023 22:45:20 +0000 (15:45 -0700)
committerWolfgang Bangerth <bangerth@colostate.edu>
Thu, 26 Jan 2023 22:45:20 +0000 (15:45 -0700)
include/deal.II/base/communication_pattern_base.h
include/deal.II/base/mpi_noncontiguous_partitioner.h
include/deal.II/base/partitioner.h

index 94c27e16a199a798f69023116b98db8f8e2c9f48..fb3c24c8f3ff8b33f43c09692936b6a808fd5cf6 100644 (file)
@@ -38,6 +38,31 @@ namespace Utilities
      * from the data that needs to be communicated. The goal is to reuse the
      * same communication pattern for different containers. This is similar to
      * the way SparseMatrix and SparsityPattern works.
+     *
+     * Conceptually, this class operates under the assumption that data
+     * are stored in one linearly indexed array of which every MPI process
+     * stores some elements (or possibly none). In practice it does of course
+     * not matter whether the elements are stored in contiguous arrays; the
+     * point is simply that a single integer index can be used to access a
+     * specific element. The elements of this logical array are (or at least
+     * may be) stored on different processes in a parallel MPI universe.
+     *
+     * In this world view, every process has a set of elements and their
+     * indices in the array form the "locally owned indices". Next, every
+     * process will as part of executing an algorithm require access to some
+     * of the elements stored elsewhere; we call the indices of these elements
+     * the "ghost indices", in analogy to how vectors and triangulations
+     * partition vector elements and mesh cells into locally
+     * owned ones and ghost ones (along, of course, with cells and ghosts stored
+     * on other processes that the current process simply does not care about
+     * and consequently needs not know anything about).
+     *
+     * The point of this class (and its implementations in derived classes) is
+     * to set up communication infrastructure whereby one process can inquire
+     * efficiently about the "ghost elements" stored on other processes, and
+     * to send those locally owned elements to those other processes that
+     * require knowledge of their value because they list these elements among
+     * their respective "ghost indices".
      */
     class CommunicationPatternBase
     {
@@ -48,15 +73,21 @@ namespace Utilities
       virtual ~CommunicationPatternBase() = default;
 
       /**
-       * Reinitialize the communication pattern. The first argument
-       * `vector_space_vector_index_set` is the index set associated to a
-       * VectorSpaceVector object. The second argument
-       * `read_write_vector_index_set` is the index set associated to a
-       * ReadWriteVector object.
+       * Reinitialize the communication pattern.
+       *
+       * @param[in] locally_owned_indices The set of indices of elements
+       *   in the array mentioned in the class documentation that are
+       *   stored on the current process.
+       * @param[in] ghost_indices The set of indices of elements in the
+       *   array mentioned in the class documentation that the current
+       *   process will need to be able to import.
+       * @param[in] communicator The MPI communicator used to describe the
+       *   entire set of processes that participate in the storage and
+       *   access to elements of the array.
        */
       virtual void
-      reinit(const IndexSet &vector_space_vector_index_set,
-             const IndexSet &read_write_vector_index_set,
+      reinit(const IndexSet &locally_owned_indices,
+             const IndexSet &ghost_indices,
              const MPI_Comm &communicator) = 0;
 
       /**
index 776436167f6fe3b3eb473208c94be1b8c1215808..e0927bbf004e67ded9bf9edaa5ca5bebfc332bee 100644 (file)
@@ -33,7 +33,16 @@ namespace Utilities
   {
     /**
      * A flexible Partitioner class, which does not impose restrictions
-     * regarding the order of the underlying index sets.
+     * regarding the order of the underlying index sets. In other words,
+     * this class implements the interface of the
+     * Utilities::MPI::CommunicationPatternBase base class with no
+     * assumption that every process stores a contiguous part of the
+     * array of objects, but that indeed the locally owned indices
+     * can be an arbitrary subset of all indices of elements of the array
+     * to which they refer.
+     *
+     * If you want to store only contiguous parts of these arrays on
+     * each process, take a look at Utilities::MPI::Partitioner.
      */
     class NoncontiguousPartitioner
       : public Utilities::MPI::CommunicationPatternBase
@@ -185,20 +194,19 @@ namespace Utilities
       const MPI_Comm &
       get_mpi_communicator() const override;
 
-      /**
-       * Initialize the inner data structures.
-       */
       void
-      reinit(const IndexSet &indexset_locally_owned,
-             const IndexSet &indexset_ghost,
+      reinit(const IndexSet &locally_owned_indices,
+             const IndexSet &ghost_indices,
              const MPI_Comm &communicator) override;
 
       /**
-       * Initialize the inner data structures.
+       * Initialize the inner data structures using explicit sets of
+       * indices. See the documentation of the other reinit() function for
+       * what the function does.
        */
       void
-      reinit(const std::vector<types::global_dof_index> &indices_locally_owned,
-             const std::vector<types::global_dof_index> &indices_ghost,
+      reinit(const std::vector<types::global_dof_index> &locally_owned_indices,
+             const std::vector<types::global_dof_index> &ghost_indices,
              const MPI_Comm &                            communicator);
 
     private:
index c1ec508ff771486f98b257f0d75450a97dad196e..8a66433d7e0156c1b57bca4678bed15410e744d6 100644 (file)
@@ -41,7 +41,13 @@ namespace Utilities
      * fact, any linear data structure) among processors using MPI.
      *
      * The partitioner stores the global vector size and the locally owned
-     * range as a half-open interval [@p lower, @p upper) on each process.
+     * range as a half-open interval [@p lower, @p upper) on each process. In
+     * other words, it assumes that every process stores a contiguous subset
+     * of the array mentioned in the documentation of the base class
+     * Utilities::MPI::CommunicationPatternBase. (If you want to store
+     * non-contiguous parts of these arrays on each process, take a look
+     * at Utilities::MPI::NoncontiguousPartitioner.)
+     *
      * Furthermore, it includes a structure for the point-to-point communication
      * patterns. It allows the inclusion of ghost indices (i.e. indices that a
      * current processor needs to have access to, but are owned by another

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.