From f6ef8977de547c05a46eeba891ae03ebbc79f34e Mon Sep 17 00:00:00 2001 From: Wolfgang Bangerth Date: Thu, 26 Jan 2023 15:45:20 -0700 Subject: [PATCH] Better document CommunicationPatternBase and derived classes. --- .../deal.II/base/communication_pattern_base.h | 45 ++++++++++++++++--- .../base/mpi_noncontiguous_partitioner.h | 26 +++++++---- include/deal.II/base/partitioner.h | 8 +++- 3 files changed, 62 insertions(+), 17 deletions(-) diff --git a/include/deal.II/base/communication_pattern_base.h b/include/deal.II/base/communication_pattern_base.h index 94c27e16a1..fb3c24c8f3 100644 --- a/include/deal.II/base/communication_pattern_base.h +++ b/include/deal.II/base/communication_pattern_base.h @@ -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; /** diff --git a/include/deal.II/base/mpi_noncontiguous_partitioner.h b/include/deal.II/base/mpi_noncontiguous_partitioner.h index 776436167f..e0927bbf00 100644 --- a/include/deal.II/base/mpi_noncontiguous_partitioner.h +++ b/include/deal.II/base/mpi_noncontiguous_partitioner.h @@ -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 &indices_locally_owned, - const std::vector &indices_ghost, + reinit(const std::vector &locally_owned_indices, + const std::vector &ghost_indices, const MPI_Comm & communicator); private: diff --git a/include/deal.II/base/partitioner.h b/include/deal.II/base/partitioner.h index c1ec508ff7..8a66433d7e 100644 --- a/include/deal.II/base/partitioner.h +++ b/include/deal.II/base/partitioner.h @@ -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 -- 2.39.5