From 1735f0496073e71f6d60630841a9be81521274b2 Mon Sep 17 00:00:00 2001 From: Peter Munch Date: Sun, 10 Apr 2022 16:03:26 +0200 Subject: [PATCH] Introduce MGTransferBlockMatrixFreeBase --- .../multigrid/mg_transfer_matrix_free.h | 184 ++++++++++-------- source/multigrid/mg_transfer_matrix_free.cc | 76 ++++---- .../multigrid/mg_transfer_matrix_free.inst.in | 4 + 3 files changed, 142 insertions(+), 122 deletions(-) diff --git a/include/deal.II/multigrid/mg_transfer_matrix_free.h b/include/deal.II/multigrid/mg_transfer_matrix_free.h index 3e558e0954..84ea6ba8b0 100644 --- a/include/deal.II/multigrid/mg_transfer_matrix_free.h +++ b/include/deal.II/multigrid/mg_transfer_matrix_free.h @@ -310,78 +310,17 @@ private: }; + /** - * Implementation of the MGTransferBase interface for which the transfer - * operations is implemented in a matrix-free way based on the interpolation - * matrices of the underlying finite element. This requires considerably less - * memory than MGTransferPrebuilt and can also be considerably faster than - * that variant. - * - * This class works with LinearAlgebra::distributed::BlockVector and - * performs exactly the same transfer operations for each block as - * MGTransferMatrixFree. - * Both the cases that the same DoFHandler is used for all the blocks - * and that each block uses its own DoFHandler are supported. + * Base class of MGTransferBlockMatrixFree. While MGTransferBlockMatrixFree + * constains all the setup routines of the transfer operators for the blocks, + * this class simply applies them, e.g., for restricting and prolongating. */ -template -class MGTransferBlockMatrixFree +template +class MGTransferBlockMatrixFreeBase : public MGTransferBase> { public: - /** - * Constructor without constraint matrices. Use this constructor only with - * discontinuous finite elements or with no local refinement. - */ - MGTransferBlockMatrixFree() = default; - - /** - * Constructor with constraints. Equivalent to the default constructor - * followed by initialize_constraints(). - */ - MGTransferBlockMatrixFree(const MGConstrainedDoFs &mg_constrained_dofs); - - /** - * Same as above for the case that each block has its own DoFHandler. - */ - MGTransferBlockMatrixFree( - const std::vector &mg_constrained_dofs); - - /** - * Destructor. - */ - virtual ~MGTransferBlockMatrixFree() override = default; - - /** - * Initialize the constraints to be used in build(). - */ - void - initialize_constraints(const MGConstrainedDoFs &mg_constrained_dofs); - - /** - * Same as above for the case that each block has its own DoFHandler. - */ - void - initialize_constraints( - const std::vector &mg_constrained_dofs); - - /** - * Reset the object to the state it had right after the default constructor. - */ - void - clear(); - - /** - * Actually build the information for the prolongation for each level. - */ - void - build(const DoFHandler &dof_handler); - - /** - * Same as above for the case that each block has its own DoFHandler. - */ - void - build(const std::vector *> &dof_handler); - /** * Prolongate a vector from level to_level-1 to level * to_level using the embedding matrices of the underlying finite @@ -481,29 +420,106 @@ public: const MGLevelObject> &src) const; - /** - * Memory used by this object. - */ - std::size_t - memory_consumption() const; - /** * This class can both be used with a single DoFHandler * or a separate DoFHandler for each block. */ static const bool supports_dof_handler_vector = true; -private: +protected: /** * Non-block matrix-free versions of transfer operation. */ - std::vector> matrix_free_transfer_vector; + std::vector matrix_free_transfer_vector; /** * A flag to indicate whether the same DoFHandler is used for all * the components or if each block has its own DoFHandler. */ - const bool same_for_all; + bool same_for_all; +}; + + + +/** + * Implementation of the MGTransferBase interface for which the transfer + * operations is implemented in a matrix-free way based on the interpolation + * matrices of the underlying finite element. This requires considerably less + * memory than MGTransferPrebuilt and can also be considerably faster than + * that variant. + * + * This class works with LinearAlgebra::distributed::BlockVector and + * performs exactly the same transfer operations for each block as + * MGTransferMatrixFree. + * Both the cases that the same DoFHandler is used for all the blocks + * and that each block uses its own DoFHandler are supported. + */ +template +class MGTransferBlockMatrixFree + : public MGTransferBlockMatrixFreeBase> +{ +public: + /** + * Constructor without constraint matrices. Use this constructor only with + * discontinuous finite elements or with no local refinement. + */ + MGTransferBlockMatrixFree() = default; + + /** + * Constructor with constraints. Equivalent to the default constructor + * followed by initialize_constraints(). + */ + MGTransferBlockMatrixFree(const MGConstrainedDoFs &mg_constrained_dofs); + + /** + * Same as above for the case that each block has its own DoFHandler. + */ + MGTransferBlockMatrixFree( + const std::vector &mg_constrained_dofs); + + /** + * Destructor. + */ + virtual ~MGTransferBlockMatrixFree() override = default; + + /** + * Initialize the constraints to be used in build(). + */ + void + initialize_constraints(const MGConstrainedDoFs &mg_constrained_dofs); + + /** + * Same as above for the case that each block has its own DoFHandler. + */ + void + initialize_constraints( + const std::vector &mg_constrained_dofs); + + /** + * Reset the object to the state it had right after the default constructor. + */ + void + clear(); + + /** + * Actually build the information for the prolongation for each level. + */ + void + build(const DoFHandler &dof_handler); + + /** + * Same as above for the case that each block has its own DoFHandler. + */ + void + build(const std::vector *> &dof_handler); + + /** + * Memory used by this object. + */ + std::size_t + memory_consumption() const; }; @@ -607,10 +623,10 @@ MGTransferMatrixFree::interpolate_to_mg( -template +template template void -MGTransferBlockMatrixFree::copy_to_mg( +MGTransferBlockMatrixFreeBase::copy_to_mg( const DoFHandler & dof_handler, MGLevelObject> &dst, const LinearAlgebra::distributed::BlockVector & src) const @@ -629,10 +645,10 @@ MGTransferBlockMatrixFree::copy_to_mg( -template +template template void -MGTransferBlockMatrixFree::copy_to_mg( +MGTransferBlockMatrixFreeBase::copy_to_mg( const std::vector *> & dof_handler, MGLevelObject> &dst, const LinearAlgebra::distributed::BlockVector & src) const @@ -720,10 +736,10 @@ MGTransferBlockMatrixFree::copy_to_mg( } } -template +template template void -MGTransferBlockMatrixFree::copy_from_mg( +MGTransferBlockMatrixFreeBase::copy_from_mg( const DoFHandler & dof_handler, LinearAlgebra::distributed::BlockVector &dst, const MGLevelObject> &src) @@ -736,10 +752,10 @@ MGTransferBlockMatrixFree::copy_from_mg( copy_from_mg(mg_dofs, dst, src); } -template +template template void -MGTransferBlockMatrixFree::copy_from_mg( +MGTransferBlockMatrixFreeBase::copy_from_mg( const std::vector *> &dof_handler, LinearAlgebra::distributed::BlockVector & dst, const MGLevelObject> &src) diff --git a/source/multigrid/mg_transfer_matrix_free.cc b/source/multigrid/mg_transfer_matrix_free.cc index 6cb3169617..19fcb4295f 100644 --- a/source/multigrid/mg_transfer_matrix_free.cc +++ b/source/multigrid/mg_transfer_matrix_free.cc @@ -707,9 +707,9 @@ MGTransferMatrixFree::memory_consumption() const template MGTransferBlockMatrixFree::MGTransferBlockMatrixFree( const MGConstrainedDoFs &mg_c) - : same_for_all(true) { - matrix_free_transfer_vector.emplace_back(mg_c); + this->same_for_all = true; + this->matrix_free_transfer_vector.emplace_back(mg_c); } @@ -717,10 +717,10 @@ MGTransferBlockMatrixFree::MGTransferBlockMatrixFree( template MGTransferBlockMatrixFree::MGTransferBlockMatrixFree( const std::vector &mg_c) - : same_for_all(false) { + this->same_for_all = false; for (const auto &constrained_block_dofs : mg_c) - matrix_free_transfer_vector.emplace_back(constrained_block_dofs); + this->matrix_free_transfer_vector.emplace_back(constrained_block_dofs); } @@ -730,13 +730,13 @@ void MGTransferBlockMatrixFree::initialize_constraints( const MGConstrainedDoFs &mg_c) { - Assert(same_for_all, + Assert(this->same_for_all, ExcMessage("This object was initialized with support for usage with " "one DoFHandler for each block, but this method assumes " "that the same DoFHandler is used for all the blocks!")); - AssertDimension(matrix_free_transfer_vector.size(), 1); + AssertDimension(this->matrix_free_transfer_vector.size(), 1); - matrix_free_transfer_vector[0].initialize_constraints(mg_c); + this->matrix_free_transfer_vector[0].initialize_constraints(mg_c); } @@ -746,24 +746,26 @@ void MGTransferBlockMatrixFree::initialize_constraints( const std::vector &mg_c) { - Assert(!same_for_all, + Assert(!this->same_for_all, ExcMessage("This object was initialized with support for using " "the same DoFHandler for all the blocks, but this " "method assumes that there is a separate DoFHandler " "for each block!")); - AssertDimension(matrix_free_transfer_vector.size(), mg_c.size()); + AssertDimension(this->matrix_free_transfer_vector.size(), mg_c.size()); for (unsigned int i = 0; i < mg_c.size(); ++i) - matrix_free_transfer_vector[i].initialize_constraints(mg_c[i]); + this->matrix_free_transfer_vector[i].initialize_constraints(mg_c[i]); } template void -MGTransferBlockMatrixFree::clear() +MGTransferBlockMatrixFree::build( + const DoFHandler &dof_handler) { - matrix_free_transfer_vector.clear(); + AssertDimension(this->matrix_free_transfer_vector.size(), 1); + this->matrix_free_transfer_vector[0].build(dof_handler); } @@ -771,29 +773,39 @@ MGTransferBlockMatrixFree::clear() template void MGTransferBlockMatrixFree::build( - const DoFHandler &dof_handler) + const std::vector *> &dof_handler) { - AssertDimension(matrix_free_transfer_vector.size(), 1); - matrix_free_transfer_vector[0].build(dof_handler); + AssertDimension(this->matrix_free_transfer_vector.size(), dof_handler.size()); + for (unsigned int i = 0; i < dof_handler.size(); ++i) + this->matrix_free_transfer_vector[i].build(*dof_handler[i]); } template -void -MGTransferBlockMatrixFree::build( - const std::vector *> &dof_handler) +std::size_t +MGTransferBlockMatrixFree::memory_consumption() const { - AssertDimension(matrix_free_transfer_vector.size(), dof_handler.size()); - for (unsigned int i = 0; i < dof_handler.size(); ++i) - matrix_free_transfer_vector[i].build(*dof_handler[i]); + std::size_t total_memory_consumption = 0; + for (const auto &el : this->matrix_free_transfer_vector) + total_memory_consumption += el.memory_consumption(); + return total_memory_consumption; } template void -MGTransferBlockMatrixFree::prolongate( +MGTransferBlockMatrixFree::clear() +{ + this->matrix_free_transfer_vector.clear(); +} + + + +template +void +MGTransferBlockMatrixFreeBase::prolongate( const unsigned int to_level, LinearAlgebra::distributed::BlockVector & dst, const LinearAlgebra::distributed::BlockVector &src) const @@ -815,9 +827,9 @@ MGTransferBlockMatrixFree::prolongate( -template +template void -MGTransferBlockMatrixFree::prolongate_and_add( +MGTransferBlockMatrixFreeBase::prolongate_and_add( const unsigned int to_level, LinearAlgebra::distributed::BlockVector & dst, const LinearAlgebra::distributed::BlockVector &src) const @@ -839,9 +851,9 @@ MGTransferBlockMatrixFree::prolongate_and_add( -template +template void -MGTransferBlockMatrixFree::restrict_and_add( +MGTransferBlockMatrixFreeBase::restrict_and_add( const unsigned int from_level, LinearAlgebra::distributed::BlockVector & dst, const LinearAlgebra::distributed::BlockVector &src) const @@ -863,18 +875,6 @@ MGTransferBlockMatrixFree::restrict_and_add( -template -std::size_t -MGTransferBlockMatrixFree::memory_consumption() const -{ - std::size_t total_memory_consumption = 0; - for (const auto &el : matrix_free_transfer_vector) - total_memory_consumption += el.memory_consumption(); - return total_memory_consumption; -} - - - // explicit instantiation #include "mg_transfer_matrix_free.inst" diff --git a/source/multigrid/mg_transfer_matrix_free.inst.in b/source/multigrid/mg_transfer_matrix_free.inst.in index 2967b4c71a..868b82d799 100644 --- a/source/multigrid/mg_transfer_matrix_free.inst.in +++ b/source/multigrid/mg_transfer_matrix_free.inst.in @@ -18,5 +18,9 @@ for (deal_II_dimension : DIMENSIONS; S1 : REAL_SCALARS) { template class MGTransferMatrixFree; + template class MGTransferBlockMatrixFreeBase< + deal_II_dimension, + S1, + MGTransferMatrixFree>; template class MGTransferBlockMatrixFree; } -- 2.39.5