]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Provide AffineConstraints with the set of locally owned DoFs.
authorWolfgang Bangerth <bangerth@colostate.edu>
Thu, 24 Aug 2023 23:35:26 +0000 (17:35 -0600)
committerWolfgang Bangerth <bangerth@colostate.edu>
Fri, 15 Sep 2023 00:19:21 +0000 (18:19 -0600)
include/deal.II/lac/affine_constraints.h
include/deal.II/lac/affine_constraints.templates.h

index 635f82fd6481e07205b3c0140a5c2125a0864ef6..7426cfa9eb350bed07578fda41f42f1e667670df 100644 (file)
@@ -567,8 +567,8 @@ public:
   AffineConstraints();
 
   /**
-   * Constructor. The supplied IndexSet defines which indices might be
-   * constrained inside this AffineConstraints container. In a calculation
+   * Constructor. The supplied IndexSet defines for which indices
+   * this object will store constraints. In a calculation
    * with a DoFHandler object based on parallel::distributed::Triangulation
    * or parallel::shared::Triangulation, one should use the set of locally
    * relevant DoFs (see
@@ -576,9 +576,53 @@ public:
    *
    * The given IndexSet allows the AffineConstraints container to save
    * memory by just not caring about degrees of freedom that are not of
-   * importance to the current processor.
+   * importance to the current processor. In contrast, in parallel
+   * computations, if you do not provide such an index set (here, or
+   * using the reinit() function that takes such an argument), the
+   * current object will allocate memory proportional to the *total*
+   * number of degrees of freedom (accumulated over all processes),
+   * which is clearly wasteful and not efficient -- and should be considered
+   * a bug.
+   *
+   * @note This constructor is equivalent to calling the following one with
+   *   both of its arguments equal to the index set provided here. This
+   *   is not wrong, but inefficient. Use the following constructor instead.
    */
-  explicit AffineConstraints(const IndexSet &local_constraints);
+  explicit AffineConstraints(const IndexSet &locally_stored_constraints);
+
+  /**
+   * Constructor. The supplied IndexSet objects define which of the indices
+   * the object may encounter are locally owned, and for which indices
+   * this object will store constraints. In a calculation
+   * with a DoFHandler object based on parallel::distributed::Triangulation
+   * or parallel::shared::Triangulation, the first of these arguments
+   * should be the set of locally owned indices (see
+   * @ref GlossLocallyOwnedDof)
+   * and the second argument should use the set of locally
+   * relevant DoFs (see
+   * @ref GlossLocallyRelevantDof).
+   *
+   * The @p locally_owned_dofs argument indicates which elements of
+   * vectors that will be passed to the member functions of this
+   * object are locally owned (which is generally equal to the set of
+   * locally owned degrees of freedom -- thus the name of the
+   * argument). Knowledge of this set makes it possible to up front,
+   * in close(), determine which vector elements need to be imported
+   * when you later call the distribute() or related functions.
+   *
+   * The @p locally_stored_constraints IndexSet allows the AffineConstraints
+   * container to save memory by just not caring about degrees of freedom
+   * that are not of importance to the current processor. In contrast,
+   * in parallel computations, if you do not provide such an index set
+   * (here, or using the reinit() function that takes such an argument), the
+   * current object will allocate memory proportional to the *total*
+   * number of degrees of freedom (accumulated over all processes),
+   * which is clearly wasteful and not efficient -- and should be considered
+   * a bug. As mentioned above, one typically wants to provide the set of
+   * locally relevant degrees of freedom for this argument.
+   */
+  AffineConstraints(const IndexSet &locally_owned_dofs,
+                    const IndexSet &locally_stored_constraints);
 
   /**
    * Copy constructor
@@ -621,13 +665,36 @@ public:
   copy_from(const AffineConstraints<other_number> &other);
 
   /**
-   * clear() the AffineConstraints object and supply an IndexSet with lines
-   * that may be constrained. This function is only relevant in the
-   * distributed case to supply a different IndexSet. Otherwise this routine
-   * is equivalent to calling clear(). See the constructor for details.
+   * Reset the AffineConstraints object. This results in an object that
+   * is as if it had been default-constructed. See the discussion in the
+   * documentation of the various constructors of this class for whether
+   * this is what you want -- in particular, for parallel computations,
+   * you do not want to call this function, but rather one of the other
+   * reinit() functions that take index sets.
    */
   void
-  reinit(const IndexSet &local_constraints = IndexSet());
+  reinit();
+
+  /**
+   * clear() the AffineConstraints object and supply an IndexSet that describes
+   * for which degrees of freedom this object can store constraints. See
+   * the discussion in the documentation of the constructor of this class
+   * that takes a single index set as argument.
+   */
+  void
+  reinit(const IndexSet &locally_stored_constraints);
+
+  /**
+   * clear() the AffineConstraints object and supply both an IndexSet
+   * object that describes which degrees of freedom the current process
+   * owns, and an IndexSet that describes
+   * for which degrees of freedom this object can store constraints. See
+   * the discussion in the documentation of the constructor of this class
+   * that takes two index sets as argument for more information.
+   */
+  void
+  reinit(const IndexSet &locally_owned_dofs,
+         const IndexSet &locally_stored_constraints);
 
   /**
    * Determines if we can store a constraint for the given @p line_n. This
@@ -1947,6 +2014,14 @@ private:
    */
   std::vector<size_type> lines_cache;
 
+  /**
+   * This IndexSet denotes the set of locally owned DoFs (or, more correctly,
+   * the locally owned vector elements when operating on parallel vectors)
+   * and is used in the distribute() function when determining for which
+   * degrees of freedom to import dependent DoFs.
+   */
+  IndexSet locally_owned_dofs;
+
   /**
    * This IndexSet is used to limit the lines to save in the AffineConstraints
    * to a subset. This is necessary, because the lines_cache vector would
@@ -2048,18 +2123,24 @@ private:
 
 template <typename number>
 inline AffineConstraints<number>::AffineConstraints()
-  : AffineConstraints<number>(IndexSet())
+  : AffineConstraints<number>(IndexSet(), IndexSet())
 {}
 
 
 
 template <typename number>
 inline AffineConstraints<number>::AffineConstraints(
-  const IndexSet &local_constraints)
+  const IndexSet &locally_owned_dofs,
+  const IndexSet &locally_stored_constraints)
   : lines()
-  , local_lines(local_constraints)
+  , locally_owned_dofs(locally_owned_dofs)
+  , local_lines(locally_stored_constraints)
   , sorted(false)
 {
+  Assert(locally_owned_dofs.is_subset_of(locally_stored_constraints),
+         ExcMessage("The set of locally stored constraints needs to be a "
+                    "superset of the locally owned DoFs."));
+
   // make sure the IndexSet is compressed. Otherwise this can lead to crashes
   // that are hard to find (only happen in release mode).
   // see tests/mpi/affine_constraints_crash_01
@@ -2068,6 +2149,15 @@ inline AffineConstraints<number>::AffineConstraints(
 
 
 
+template <typename number>
+inline AffineConstraints<number>::AffineConstraints(
+  const IndexSet &locally_stored_constraints)
+  : AffineConstraints<number>(locally_stored_constraints,
+                              locally_stored_constraints)
+{}
+
+
+
 template <typename number>
 inline AffineConstraints<number>::AffineConstraints(
   const AffineConstraints &affine_constraints)
index 55e2d3e0a609c58af5b512a5fb467b63cc535b79..729a86b781627c884060c2f6479b4e7342a65483 100644 (file)
@@ -995,29 +995,31 @@ AffineConstraints<number>::close()
 
   // Now also compute which indices we need in distribute().
   //
-  // This processor owns only part of the vector. one may think that
-  // every processor should be able to simply communicate those elements
-  // it owns and for which it knows that they act as sources to constrained
-  // DoFs to the owner of these DoFs. This would lead to a scheme where all
-  // we need to do is to add some local elements to (possibly non-local)
-  // ones and then call compress().
+  // This processor owns only part of the vectors we will work on. One may think
+  // that every processor should be able to simply communicate those elements it
+  // owns and for which it knows that they act as sources to constrained DoFs to
+  // the owner of these DoFs. This would lead to a scheme where all we need to
+  // do is to add some local elements to (possibly non-local) ones and then call
+  // compress().
   //
-  // Alas, this scheme does not work as evidenced by the disaster of bug
-  // #51, see http://code.google.com/p/dealii/issues/detail?id=51 and the
-  // reversion of one attempt that implements this in r29662. Rather, we
-  // need to get a vector that has all the *sources* or constraints we
-  // own locally, possibly as ghost vector elements, then read from them,
-  // and finally throw away the ghosted vector. Implement this in the
-  // following.
-  needed_elements_for_distribute = local_lines;
-
-  if (needed_elements_for_distribute != complete_index_set(local_lines.size()))
+  // Alas, this scheme does not work as evidenced by the disaster of
+  // bug #51 (originally stored in the Google Code repository, but now
+  // unfortunately no longer available because that platform has gone
+  // away) and reversion of one attempt that implements this in svn
+  // r29662. Rather, we need to get a vector that has all the
+  // *sources* or constraints we own locally, possibly as ghost vector
+  // elements, then read from them, and finally throw away the ghosted
+  // vector. Implement this in the following.
+  needed_elements_for_distribute = locally_owned_dofs;
+
+  if (needed_elements_for_distribute !=
+      complete_index_set(locally_owned_dofs.size()))
     {
       std::vector<types::global_dof_index> additional_elements;
       for (const ConstraintLine &line : lines)
-        if (local_lines.is_element(line.index))
+        if (locally_owned_dofs.is_element(line.index))
           for (const std::pair<size_type, number> &entry : line.entries)
-            if (!local_lines.is_element(entry.first))
+            if (!locally_owned_dofs.is_element(entry.first))
               additional_elements.emplace_back(entry.first);
       std::sort(additional_elements.begin(), additional_elements.end());
       needed_elements_for_distribute.add_indices(additional_elements.begin(),
@@ -1094,6 +1096,9 @@ AffineConstraints<number>::clear()
     lines_cache.swap(tmp);
   }
 
+  locally_owned_dofs             = {};
+  needed_elements_for_distribute = {};
+
   sorted = false;
 }
 
@@ -1101,16 +1106,43 @@ AffineConstraints<number>::clear()
 
 template <typename number>
 void
-AffineConstraints<number>::reinit(const IndexSet &local_constraints)
+AffineConstraints<number>::reinit()
 {
-  local_lines = local_constraints;
+  reinit(IndexSet(), IndexSet());
+}
+
+
+
+template <typename number>
+void
+AffineConstraints<number>::reinit(const IndexSet &locally_stored_constraints)
+{
+  reinit(locally_stored_constraints, locally_stored_constraints);
+}
+
+
+
+template <typename number>
+void
+AffineConstraints<number>::reinit(const IndexSet &locally_owned_dofs,
+                                  const IndexSet &locally_stored_constraints)
+{
+  // First clear previous content
+  clear();
+
+  // Then set the objects that describe the index sets of DoFs we care about:
+  Assert(locally_owned_dofs.is_subset_of(locally_stored_constraints),
+         ExcMessage("The set of locally stored constraints needs to be a "
+                    "superset of the locally owned DoFs."));
+
+  this->locally_owned_dofs = locally_owned_dofs;
+  this->local_lines        = locally_stored_constraints;
 
   // make sure the IndexSet is compressed. Otherwise this can lead to crashes
   // that are hard to find (only happen in release mode).
   // see tests/mpi/affine_constraints_crash_01
-  local_lines.compress();
-
-  clear();
+  this->locally_owned_dofs.compress();
+  this->local_lines.compress();
 }
 
 
@@ -2567,9 +2599,10 @@ AffineConstraints<number>::distribute(VectorType &vec) const
       // of the locally-owned indices. But you never know what people
       // pass as arguments...
 #ifdef DEBUG
-      for (const auto i : vec_owned_elements)
-        Assert(needed_elements_for_distribute.is_element(i),
-               ExcInternalError());
+      if (needed_elements_for_distribute != IndexSet())
+        for (const auto i : vec_owned_elements)
+          Assert(needed_elements_for_distribute.is_element(i),
+                 ExcInternalError());
 #endif
 
       VectorType ghosted_vector;

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.