From: Matthias Maier Date: Fri, 25 May 2018 22:16:13 +0000 (-0500) Subject: RCM: Update lac, part 2: constrained_linear_operator X-Git-Tag: v9.1.0-rc1~1067^2~22 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=a3cba727c8ce7432c23ced77f5aab0370055c0ea;p=dealii.git RCM: Update lac, part 2: constrained_linear_operator --- diff --git a/include/deal.II/lac/constrained_linear_operator.h b/include/deal.II/lac/constrained_linear_operator.h index 0560f20a3a..96eccd4fc8 100644 --- a/include/deal.II/lac/constrained_linear_operator.h +++ b/include/deal.II/lac/constrained_linear_operator.h @@ -16,7 +16,7 @@ #ifndef dealii_constrained_linear_operator_h #define dealii_constrained_linear_operator_h -#include +#include #include #include @@ -31,16 +31,16 @@ DEAL_II_NAMESPACE_OPEN /** - * This function takes a ConstraintMatrix @p constraint_matrix and an operator - * exemplar @p exemplar (this exemplar is usually a linear operator that - * describes the system matrix - it is only used to create domain and range - * vectors of appropriate sizes, its action vmult is never used). A - * LinearOperator object associated with the "homogeneous action" of the - * underlying ConstraintMatrix object is returned: + * This function takes an AffineConstraints object @p constraints and + * an operator exemplar @p exemplar (this exemplar is usually a linear + * operator that describes the system matrix - it is only used to create + * domain and range vectors of appropriate sizes, its action vmult + * is never used). A LinearOperator object associated with the "homogeneous + * action" of the underlying AffineConstraints object is returned: * * Applying the LinearOperator object on a vector u results in a * vector v that stores the result of calling - * ConstraintMatrix::distribute() on u - with one important + * AffineConstraints::distribute() on u - with one important * difference: inhomogeneities are not applied, but always treated as 0 * instead. * @@ -62,12 +62,12 @@ DEAL_II_NAMESPACE_OPEN template LinearOperator distribute_constraints_linear_operator( - const ConstraintMatrix & constraint_matrix, - const LinearOperator &exemplar) + const AffineConstraints &constraints, + const LinearOperator & exemplar) { LinearOperator return_op = exemplar; - return_op.vmult_add = [&constraint_matrix](Range &v, const Domain &u) { + return_op.vmult_add = [&constraints](Range &v, const Domain &u) { Assert(!dealii::PointerComparison::equal(&v, &u), dealii::ExcMessage("The domain and range vectors must be different " "storage locations")); @@ -77,7 +77,7 @@ distribute_constraints_linear_operator( v += u; const auto &locally_owned_elements = v.locally_owned_elements(); - for (const auto &line : constraint_matrix.get_lines()) + for (const auto &line : constraints.get_lines()) { const auto i = line.index; if (locally_owned_elements.is_element(i)) @@ -95,7 +95,7 @@ distribute_constraints_linear_operator( v.compress(VectorOperation::add); }; - return_op.Tvmult_add = [&constraint_matrix](Domain &v, const Range &u) { + return_op.Tvmult_add = [&constraints](Domain &v, const Range &u) { Assert(!dealii::PointerComparison::equal(&v, &u), dealii::ExcMessage("The domain and range vectors must be different " "storage locations")); @@ -105,7 +105,7 @@ distribute_constraints_linear_operator( v += u; const auto &locally_owned_elements = v.locally_owned_elements(); - for (const auto &line : constraint_matrix.get_lines()) + for (const auto &line : constraints.get_lines()) { const auto i = line.index; @@ -145,7 +145,7 @@ distribute_constraints_linear_operator( /** - * Given a ConstraintMatrix @p constraint_matrix and an operator exemplar @p + * Given a AffineConstraints @p constraints and an operator exemplar @p * exemplar, return a LinearOperator that is the projection to the subspace of * constrained degrees of freedom, i.e. all entries of the result vector that * correspond to unconstrained degrees of freedom are set to zero. @@ -158,14 +158,14 @@ distribute_constraints_linear_operator( template LinearOperator project_to_constrained_linear_operator( - const ConstraintMatrix & constraint_matrix, - const LinearOperator &exemplar) + const AffineConstraints &constraints, + const LinearOperator & exemplar) { LinearOperator return_op = exemplar; - return_op.vmult_add = [&constraint_matrix](Range &v, const Domain &u) { + return_op.vmult_add = [&constraints](Range &v, const Domain &u) { const auto &locally_owned_elements = v.locally_owned_elements(); - for (const auto &line : constraint_matrix.get_lines()) + for (const auto &line : constraints.get_lines()) { const auto i = line.index; if (locally_owned_elements.is_element(i)) @@ -177,9 +177,9 @@ project_to_constrained_linear_operator( v.compress(VectorOperation::add); }; - return_op.Tvmult_add = [&constraint_matrix](Domain &v, const Range &u) { + return_op.Tvmult_add = [&constraints](Domain &v, const Range &u) { const auto &locally_owned_elements = v.locally_owned_elements(); - for (const auto &line : constraint_matrix.get_lines()) + for (const auto &line : constraints.get_lines()) { const auto i = line.index; if (locally_owned_elements.is_element(i)) @@ -210,7 +210,7 @@ project_to_constrained_linear_operator( /** - * Given a ConstraintMatrix object @p constraint_matrix and a LinearOperator + * Given a AffineConstraints object @p constraints and a LinearOperator * @p linop, this function creates a LinearOperator object consisting of the * composition of three operations and a regularization: * @code @@ -218,9 +218,9 @@ project_to_constrained_linear_operator( * @endcode * with * @code - * C = distribute_constraints_linear_operator(constraint_matrix, linop); + * C = distribute_constraints_linear_operator(constraints, linop); * Ct = transpose_operator(C); - * Id_c = project_to_constrained_linear_operator(constraint_matrix, linop); + * Id_c = project_to_constrained_linear_operator(constraints, linop); * @endcode * and Id_c is the projection to the subspace consisting of all * vector entries associated with constrained degrees of freedom. @@ -248,20 +248,19 @@ project_to_constrained_linear_operator( */ template LinearOperator -constrained_linear_operator(const ConstraintMatrix &constraint_matrix, - const LinearOperator &linop) +constrained_linear_operator( + const AffineConstraints &constraints, + const LinearOperator & linop) { - const auto C = - distribute_constraints_linear_operator(constraint_matrix, linop); - const auto Ct = transpose_operator(C); - const auto Id_c = - project_to_constrained_linear_operator(constraint_matrix, linop); + const auto C = distribute_constraints_linear_operator(constraints, linop); + const auto Ct = transpose_operator(C); + const auto Id_c = project_to_constrained_linear_operator(constraints, linop); return Ct * linop * C + Id_c; } /** - * Given a ConstraintMatrix object @p constraint_matrix, a LinearOperator @p + * Given a AffineConstraints object @p constraints, a LinearOperator @p * linop and a right-hand side @p right_hand_side, this function creates a * PackagedOperation that stores the following computation: * @code @@ -269,7 +268,7 @@ constrained_linear_operator(const ConstraintMatrix &constraint_matrix, * @endcode * with * @code - * C = distribute_constraints_linear_operator(constraint_matrix, linop); + * C = distribute_constraints_linear_operator(constraints, linop); * Ct = transpose_operator(C); * @endcode * @@ -296,27 +295,26 @@ constrained_linear_operator(const ConstraintMatrix &constraint_matrix, */ template PackagedOperation -constrained_right_hand_side(const ConstraintMatrix &constraint_matrix, - const LinearOperator &linop, - const Range &right_hand_side) +constrained_right_hand_side( + const AffineConstraints &constraints, + const LinearOperator & linop, + const Range & right_hand_side) { PackagedOperation return_comp; return_comp.reinit_vector = linop.reinit_range_vector; - return_comp.apply_add = - [&constraint_matrix, &linop, &right_hand_side](Range &v) { - const auto C = - distribute_constraints_linear_operator(constraint_matrix, linop); - const auto Ct = transpose_operator(C); + return_comp.apply_add = [&constraints, &linop, &right_hand_side](Range &v) { + const auto C = distribute_constraints_linear_operator(constraints, linop); + const auto Ct = transpose_operator(C); - GrowingVectorMemory vector_memory; - typename VectorMemory::Pointer k(vector_memory); - linop.reinit_domain_vector(*k, /*bool fast=*/false); - constraint_matrix.distribute(*k); + GrowingVectorMemory vector_memory; + typename VectorMemory::Pointer k(vector_memory); + linop.reinit_domain_vector(*k, /*bool fast=*/false); + constraints.distribute(*k); - v += Ct * (right_hand_side - linop * *k); - }; + v += Ct * (right_hand_side - linop * *k); + }; // lambda capture expressions are a C++14 feature... const auto apply_add = return_comp.apply_add;