From: Wolfgang Bangerth Date: Thu, 24 Aug 2023 23:35:26 +0000 (-0600) Subject: Provide AffineConstraints with the set of locally owned DoFs. X-Git-Tag: relicensing~496^2~1 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=c55de31c3a776813d0d1a08a92da122f8ea3d60b;p=dealii.git Provide AffineConstraints with the set of locally owned DoFs. --- diff --git a/include/deal.II/lac/affine_constraints.h b/include/deal.II/lac/affine_constraints.h index 635f82fd64..7426cfa9eb 100644 --- a/include/deal.II/lac/affine_constraints.h +++ b/include/deal.II/lac/affine_constraints.h @@ -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); /** - * 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 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 inline AffineConstraints::AffineConstraints() - : AffineConstraints(IndexSet()) + : AffineConstraints(IndexSet(), IndexSet()) {} template inline AffineConstraints::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::AffineConstraints( +template +inline AffineConstraints::AffineConstraints( + const IndexSet &locally_stored_constraints) + : AffineConstraints(locally_stored_constraints, + locally_stored_constraints) +{} + + + template inline AffineConstraints::AffineConstraints( const AffineConstraints &affine_constraints) diff --git a/include/deal.II/lac/affine_constraints.templates.h b/include/deal.II/lac/affine_constraints.templates.h index 55e2d3e0a6..729a86b781 100644 --- a/include/deal.II/lac/affine_constraints.templates.h +++ b/include/deal.II/lac/affine_constraints.templates.h @@ -995,29 +995,31 @@ AffineConstraints::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 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 &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::clear() lines_cache.swap(tmp); } + locally_owned_dofs = {}; + needed_elements_for_distribute = {}; + sorted = false; } @@ -1101,16 +1106,43 @@ AffineConstraints::clear() template void -AffineConstraints::reinit(const IndexSet &local_constraints) +AffineConstraints::reinit() { - local_lines = local_constraints; + reinit(IndexSet(), IndexSet()); +} + + + +template +void +AffineConstraints::reinit(const IndexSet &locally_stored_constraints) +{ + reinit(locally_stored_constraints, locally_stored_constraints); +} + + + +template +void +AffineConstraints::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::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;