From eb08e39d712dbf1924768aaf2bd6ed847baf5e8d Mon Sep 17 00:00:00 2001 From: Chayapol Chaoveeraprasit Date: Wed, 24 Jul 2024 17:08:00 +0300 Subject: [PATCH] added clear() and Tvmult for SparseVanka undo copyright notices & fix wording issues Add a changelog entry fix formatting of sparse_vanka.h fix typos add spacing --- doc/news/changes/minor/20240814chaycha | 4 ++ include/deal.II/lac/sparse_vanka.h | 35 ++++++----- include/deal.II/lac/sparse_vanka.templates.h | 61 ++++++++++++++++++-- source/lac/sparse_vanka.cc | 12 ++++ 4 files changed, 94 insertions(+), 18 deletions(-) create mode 100644 doc/news/changes/minor/20240814chaycha diff --git a/doc/news/changes/minor/20240814chaycha b/doc/news/changes/minor/20240814chaycha new file mode 100644 index 0000000000..e945e345ba --- /dev/null +++ b/doc/news/changes/minor/20240814chaycha @@ -0,0 +1,4 @@ +New: Added SparseVanka::Tvmult() and SparseVanka::clear() +SparseVanka can now be passed to MGSmootherPrecondition to be used as a multigrid smoother. +
+(Chayapol Chaoveeraprasit, 2024/08/14) \ No newline at end of file diff --git a/include/deal.II/lac/sparse_vanka.h b/include/deal.II/lac/sparse_vanka.h index 2d6a3e2b0d..63ede4c36f 100644 --- a/include/deal.II/lac/sparse_vanka.h +++ b/include/deal.II/lac/sparse_vanka.h @@ -15,8 +15,6 @@ #ifndef dealii_sparse_vanka_h #define dealii_sparse_vanka_h - - #include #include @@ -216,6 +214,12 @@ public: void Tvmult(Vector &dst, const Vector &src) const; + /** + * Clear all memory. + */ + void + clear(); + /** * Return the dimension of the codomain (or range) space. Note that the * matrix is of dimension $m \times n$. @@ -240,9 +244,11 @@ protected: /** * Apply the inverses corresponding to those degrees of freedom that have a * @p true value in @p dof_mask, to the @p src vector and move the result - * into @p dst. Actually, only values for allowed indices are written to @p - * dst, so the application of this function only does what is announced in - * the general documentation if the given mask sets all values to zero + * into @p dst. The transpose of the inverse is applied instead if + * @p transpose equals true. Actually, only values for allowed indices are + * written to @p dst, so the application of this function only does what is + * announced in the general documentation if the given mask sets all values + * to zero. * * The reason for providing the mask anyway is that in derived classes we * may want to apply the preconditioner to parts of the matrix only, in @@ -258,10 +264,12 @@ protected: * The @p vmult of this class of course calls this function with a null * pointer */ + template void apply_preconditioner(Vector &dst, const Vector &src, + const bool transpose = false, const std::vector *const dof_mask = nullptr) const; /** @@ -511,6 +519,14 @@ public: void vmult(Vector &dst, const Vector &src) const; + + /** + * Apply the transpose preconditioner. + */ + template + void + Tvmult(Vector &dst, const Vector &src) const; + /** * Determine an estimate for the memory consumption (in bytes) of this * object. @@ -564,15 +580,6 @@ SparseVanka::n() const return _n; } -template -template -inline void -SparseVanka::Tvmult(Vector & /*dst*/, - const Vector & /*src*/) const -{ - AssertThrow(false, ExcNotImplemented()); -} - #endif // DOXYGEN DEAL_II_NAMESPACE_CLOSE diff --git a/include/deal.II/lac/sparse_vanka.templates.h b/include/deal.II/lac/sparse_vanka.templates.h index 168419546e..34fda24b17 100644 --- a/include/deal.II/lac/sparse_vanka.templates.h +++ b/include/deal.II/lac/sparse_vanka.templates.h @@ -215,9 +215,42 @@ SparseVanka::vmult(Vector &dst, dst = 0; // then pass on to the function // that actually does the work - apply_preconditioner(dst, src); + apply_preconditioner(dst, src, false); } +template +template +void +SparseVanka::Tvmult(Vector &dst, + const Vector &src) const +{ + Assert(matrix != nullptr, ExcNotInitialized()); + Assert(selected != nullptr, ExcNotInitialized()); + + // first set output vector to zero + dst = 0; + // then pass on to the function + // that actually does the work + apply_preconditioner(dst, src, true); +} + + +template +void +SparseVanka::clear() +{ // Clear the inverses vector and deallocate memory + inverses.clear(); + + // Reset the matrix pointer + matrix = nullptr; + + // Reset the selected pointer + selected = nullptr; + + // Reset the sizes + _m = 0; + _n = 0; +} template template @@ -225,6 +258,7 @@ void SparseVanka::apply_preconditioner( Vector &dst, const Vector &src, + const bool transpose, const std::vector *const dof_mask) const { Assert(dst.size() == src.size(), @@ -336,7 +370,10 @@ SparseVanka::apply_preconditioner( } // apply preconditioner - inverses[row]->vmult(x, b); + if (transpose) + inverses[row]->Tvmult(x, b); + else + inverses[row]->vmult(x, b); // Distribute new values for (std::map::const_iterator is = @@ -567,8 +604,24 @@ SparseBlockVanka::vmult(Vector &dst, Threads::TaskGroup<> tasks; for (unsigned int block = 0; block < n_blocks; ++block) - tasks += Threads::new_task([&, block]() { - this->apply_preconditioner(dst, src, &dof_masks[block]); + tasks += Threads::new_task([&, block] { + this->apply_preconditioner(dst, src, false, &dof_masks[block]); + }); + tasks.join_all(); +} + +template +template +void +SparseBlockVanka::Tvmult(Vector &dst, + const Vector &src) const +{ + dst = 0; + + Threads::TaskGroup<> tasks; + for (unsigned int block = 0; block < n_blocks; ++block) + tasks += Threads::new_task([&, block] { + this->apply_preconditioner(dst, src, true, &dof_masks[block]); }); tasks.join_all(); } diff --git a/source/lac/sparse_vanka.cc b/source/lac/sparse_vanka.cc index 9cb74bb9cb..d49ba83de8 100644 --- a/source/lac/sparse_vanka.cc +++ b/source/lac/sparse_vanka.cc @@ -27,6 +27,12 @@ SparseVanka::vmult(Vector &dst, template void SparseVanka::vmult(Vector &dst, const Vector &src) const; +template void +SparseVanka::Tvmult(Vector &dst, + const Vector &src) const; +template void +SparseVanka::Tvmult(Vector &dst, + const Vector &src) const; template class SparseBlockVanka; @@ -38,5 +44,11 @@ SparseBlockVanka::vmult(Vector &dst, template void SparseBlockVanka::vmult(Vector &dst, const Vector &src) const; +template void +SparseBlockVanka::Tvmult(Vector &dst, + const Vector &src) const; +template void +SparseBlockVanka::Tvmult(Vector &dst, + const Vector &src) const; DEAL_II_NAMESPACE_CLOSE -- 2.39.5