From 0db41097da126be56b264de5ed0ca3b84b973bd3 Mon Sep 17 00:00:00 2001 From: Wolfgang Bangerth Date: Thu, 5 Oct 2023 14:52:18 -0600 Subject: [PATCH] Add AffineConstraints::add_view(). --- include/deal.II/lac/affine_constraints.h | 51 +++++++++++++++ .../lac/affine_constraints.templates.h | 62 +++++++++++++++++++ 2 files changed, 113 insertions(+) diff --git a/include/deal.II/lac/affine_constraints.h b/include/deal.II/lac/affine_constraints.h index e683d2d8c4..abc995b0cb 100644 --- a/include/deal.II/lac/affine_constraints.h +++ b/include/deal.II/lac/affine_constraints.h @@ -942,6 +942,57 @@ public: void shift(const size_type offset); + /** + * This function provides a "view" into a constraint + * object. Specifically, given a "mask" index set that describes + * which constraints we are interested in, it returns an + * AffineConstraints object that contains only those constraints + * that correspond to degrees of freedom that are listed in the + * mask, with indices shifted so that they correspond to the + * position within the mask. This process is the same as how + * IndexSet::get_view() computes the shifted indices. The function + * is typically used to extract from an AffineConstraints object + * corresponding to a DoFHandler only those constraints that + * correspond to a specific variable (say, to the velocity in a + * Stokes system) so that the resulting AffineConstraints object can + * be applied to a single block of a block vector of solutions; in + * this case, the mask would be the index set of velocity degrees of + * freedom, as a subset of all degrees of freedom. + * + * This function can only work if the degrees of freedom selected by + * the mask are constrained only against other degrees of freedom + * that are listed in the mask. In the example above, this means + * that constraints for the selected velocity degrees of freedom are + * only against other velocity degrees of freedom, but not against + * any pressure degrees of freedom. If that is not so, an assertion + * will be triggered. + * + * A typical case where this function is useful is as follows. Say, + * you have a block linear system in which you have blocks + * corresponding to variables $(u,p,T,c)$ (which you can think of as + * velocity, pressure, temperature, and chemical composition -- or + * whatever other kind of problem you are currently considering in + * your own work). Let's assume we have developed a linear solver or + * preconditioner that first solves the coupled $u$-$T$ system, and + * once that is done, solves the $p$-$c$ system. In this case, it is + * often useful to set up block vectors with only two components + * corresponding to the $u$ and $T$ components, and later for only + * the $p$-$c$ components of the solution. As part of this, you will + * want to apply constraints (using the distribute() function of this + * class) to only the 2-block vector, but for this you need to obtain + * an AffineConstraints object that represents only those constraints + * that correspond to the variables in question, and in the order + * in which they appear in the 2-block vector rather than in global + * 4-block vectors. This function allows you to extract such an + * object corresponding to a subset of constraints by applying a mask + * to the global constraints object that corresponds to the + * variables we're currently interested in. For the $u$-$T$ system, + * we need a mask that contains all indices of $u$ degrees of + * freedom as well as all indices of $T$ degrees of freedom. + */ + AffineConstraints + get_view(const IndexSet &mask) const; + /** * Clear all entries of this matrix. Reset the flag determining whether new * entries are accepted or not. diff --git a/include/deal.II/lac/affine_constraints.templates.h b/include/deal.II/lac/affine_constraints.templates.h index c00d803f46..c33b2511e9 100644 --- a/include/deal.II/lac/affine_constraints.templates.h +++ b/include/deal.II/lac/affine_constraints.templates.h @@ -1082,6 +1082,68 @@ AffineConstraints::shift(const size_type offset) +template +AffineConstraints +AffineConstraints::get_view(const IndexSet &mask) const +{ + Assert(sorted == true, ExcMatrixNotClosed()); + + AffineConstraints view; + + // If index sets were associated with the current object, take views + // of those as well: + if (locally_owned_dofs != IndexSet()) + { + Assert(locally_owned_dofs.size() == mask.size(), + ExcMessage("This operation only makes sense if the index " + "space described by the mask is the same as " + "the index space associated with the " + "AffineConstraints object into which you take " + "a view.")); + Assert(local_lines.size() == mask.size(), + ExcMessage("This operation only makes sense if the index " + "space described by the mask is the same as " + "the index space associated with the " + "AffineConstraints object into which you take " + "a view.")); + + view.reinit(locally_owned_dofs.get_view(mask), + local_lines.get_view(mask)); + } + + for (const ConstraintLine &line : lines) + if (mask.is_element(line.index)) + { + const size_type row = mask.index_within_set(line.index); + view.add_line(row); + for (const std::pair &entry : line.entries) + { + Assert( + mask.is_element(entry.first), + ExcMessage( + "In creating a view of an AffineConstraints " + "object, the constraint on degree of freedom " + + std::to_string(line.index) + " (which corresponds to the " + + std::to_string(row) + + "th degree of freedom selected in the mask) " + "is constrained against degree of freedom " + + std::to_string(entry.first) + + ", but this degree of freedom is not listed in the mask and " + "consequently cannot be transcribed into the index space " + "of the output object.")); + view.add_entry(row, + mask.index_within_set(entry.first), + entry.second); + } + view.set_inhomogeneity(row, line.inhomogeneity); + } + + view.close(); + return view; +} + + + template void AffineConstraints::clear() -- 2.39.5