From: Timo Heister Date: Thu, 7 Aug 2014 15:39:57 +0000 (-0400) Subject: patch by Alexander Grayver for condense() X-Git-Tag: v8.2.0-rc1~208^2~2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=15cc263e75b533dc24c9cf0437e83f28feca4db6;p=dealii.git patch by Alexander Grayver for condense() --- diff --git a/include/deal.II/lac/constraint_matrix.h b/include/deal.II/lac/constraint_matrix.h index 75b91077b7..282fa05f95 100644 --- a/include/deal.II/lac/constraint_matrix.h +++ b/include/deal.II/lac/constraint_matrix.h @@ -641,26 +641,6 @@ public: template void condense (BlockSparseMatrix &matrix) const; - /** - * Condense the given vector @p uncondensed into @p condensed. It is the - * user's responsibility to guarantee that all entries of @p condensed be - * zero. Note that this function does not take any inhomogeneity into - * account and throws an exception in case there are any - * inhomogeneities. Use the function using both a matrix and vector for that - * case. - * - * The @p VectorType may be a Vector, Vector, - * BlockVector<...>, a PETSc or Trilinos vector wrapper class, or - * any other type having the same interface. - * - * @deprecated The functions converting an uncondensed matrix into - * its condensed form are deprecated. Use the functions doing the - * in-place condensation leaving the size of the linear system unchanged. - */ - template - void condense (const VectorType &uncondensed, - VectorType &condensed) const DEAL_II_DEPRECATED; - /** * Condense the given vector in-place. The @p VectorType may be a * Vector, Vector, BlockVector<...>, a PETSc or @@ -669,10 +649,23 @@ public: * account and throws an exception in case there are any * inhomogeneities. Use the function using both a matrix and vector for that * case. + * Note: this function does not work for MPI vectors because it would require + * ghosted vector in that case, however all ghosted vectors are read-only objects. + * Instead, use the analogue with two input arguments. */ template void condense (VectorType &vec) const; + /** + * Same as above, but accepts two input arguments. If called in parallel, one + * vector is supposed to be a ghosted vector and the other vector without ghost + * elements. The function copies and condenses values from @p vec_ghosted into + * @p vec. + */ + template + void condense (const VectorType &vec_ghosted, + VectorType &vec) const; + /** * Condense a given matrix and a given vector. The associated matrix struct * should be condensed and compressed. It is the user's responsibility to diff --git a/include/deal.II/lac/constraint_matrix.templates.h b/include/deal.II/lac/constraint_matrix.templates.h index a1a97f6a95..d54eca03e6 100644 --- a/include/deal.II/lac/constraint_matrix.templates.h +++ b/include/deal.II/lac/constraint_matrix.templates.h @@ -73,82 +73,12 @@ ConstraintMatrix::condense (BlockSparseMatrix &uncondensed) const template void -ConstraintMatrix::condense (const VectorType &uncondensed, - VectorType &condensed) const +ConstraintMatrix::condense (const VectorType &vec_ghosted, + VectorType &vec) const { Assert (sorted == true, ExcMatrixNotClosed()); - AssertDimension (condensed.size()+n_constraints(), uncondensed.size()); - - // store for each line of the - // vector its new line number after - // compression. If the shift is -1, - // this line will be condensed away - std::vector new_line; - - new_line.reserve (uncondensed.size()); - - std::vector::const_iterator next_constraint = lines.begin(); - size_type shift = 0; - size_type n_rows = uncondensed.size(); - - if (next_constraint == lines.end()) - // if no constraint is to be handled - for (size_type row=0; row!=n_rows; ++row) - new_line.push_back (row); - else - for (size_type row=0; row!=n_rows; ++row) - if (row == next_constraint->line) - { - // this line is constrained - new_line.push_back (-1); - // note that @p lines is ordered - ++shift; - ++next_constraint; - if (next_constraint == lines.end()) - // nothing more to do; finish rest - // of loop - { - for (size_type i=row+1; ientries.size(); ++q) - condensed(new_line[next_constraint->entries[q].first]) - += - uncondensed(row) * next_constraint->entries[q].second; - - ++next_constraint; - }; -} - - - -template -void -ConstraintMatrix::condense (VectorType &vec) const -{ - Assert (sorted == true, ExcMatrixNotClosed()); + vec = vec_ghosted; // distribute all entries, and set them to zero. do so in // two loops because in the first one we need to add to elements @@ -166,7 +96,7 @@ ConstraintMatrix::condense (VectorType &vec) const ExcMessage ("Inhomogeneous constraint cannot be condensed " "without any matrix specified.")); - const typename VectorType::value_type old_value = vec(constraint_line->line); + const typename VectorType::value_type old_value = vec_ghosted(constraint_line->line); for (size_type q=0; q!=constraint_line->entries.size(); ++q) if (vec.in_local_range(constraint_line->entries[q].first) == true) vec(constraint_line->entries[q].first) @@ -188,6 +118,15 @@ ConstraintMatrix::condense (VectorType &vec) const +template +void +ConstraintMatrix::condense (VectorType &vec) const +{ + condense(vec, vec); +} + + + template void ConstraintMatrix::condense (const SparseMatrix &uncondensed,