From: Daniel Arndt Date: Sun, 1 Oct 2017 20:53:10 +0000 (+0200) Subject: Allow to use MGTransferBlockMatrixFree with a separate DoFHandler for each block X-Git-Tag: v9.0.0-rc1~995^2~2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=646d734fb0f0f864c35c3c1551172533d236996d;p=dealii.git Allow to use MGTransferBlockMatrixFree with a separate DoFHandler for each block --- diff --git a/include/deal.II/multigrid/mg_transfer_matrix_free.h b/include/deal.II/multigrid/mg_transfer_matrix_free.h index ddb94109c9..84936b99fb 100644 --- a/include/deal.II/multigrid/mg_transfer_matrix_free.h +++ b/include/deal.II/multigrid/mg_transfer_matrix_free.h @@ -245,10 +245,11 @@ private: * * This class works with LinearAlgebra::distributed::BlockVector and * performs exactly the same transfer operations for each block as - * MGTransferMatrixFree. This implies that each block should cover the - * same index space of DoFs and have the same partitioning. + * MGTransferMatrixFree. + * Both the cases that the same DoFHandler is used for all the blocks + * and that each block uses its own DoFHandler are supported. * - * @author Denis Davydov + * @author Denis Davydov, Daniel Arndt * @date 2017 */ template @@ -267,6 +268,11 @@ public: */ 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. */ @@ -277,6 +283,11 @@ public: */ 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. */ @@ -287,6 +298,11 @@ public: */ void build (const DoFHandler &mg_dof); + /** + * Same as above for the case that each block has its own DoFHandler. + */ + void build (const std::vector*> &mg_dof); + /** * Prolongate a vector from level to_level-1 to level * to_level using the embedding matrices of the underlying finite @@ -343,6 +359,15 @@ public: MGLevelObject> &dst, const LinearAlgebra::distributed::BlockVector &src) const; + /** + * Same as above for the case that each block has its own DoFHandler. + */ + template + void + copy_to_mg (const std::vector*> &mg_dof, + MGLevelObject> &dst, + const LinearAlgebra::distributed::BlockVector &src) const; + /** * Transfer from multi-level block-vector to normal vector. */ @@ -352,17 +377,38 @@ public: LinearAlgebra::distributed::BlockVector &dst, const MGLevelObject> &src) const; + /** + * Same as above for the case that each block has its own DoFHandler. + */ + template + void + copy_from_mg (const std::vector*> &mg_dof, + LinearAlgebra::distributed::BlockVector &dst, + 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: /** - * A non-block matrix-free version of transfer operation. + * Non-block matrix-free versions of transfer operation. + */ + 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. */ - MGTransferMatrixFree matrix_free_transfer; + const bool same_for_all; }; @@ -379,8 +425,31 @@ MGTransferBlockMatrixFree:: copy_to_mg (const DoFHandler &mg_dof, MGLevelObject> &dst, const LinearAlgebra::distributed::BlockVector &src) const +{ + AssertDimension(matrix_free_transfer_vector.size(),1); + Assert(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!")); + const std::vector*> mg_dofs (src.n_blocks(),&mg_dof); + + copy_to_mg(mg_dofs, dst, src); +} + +template +template +void +MGTransferBlockMatrixFree:: +copy_to_mg (const std::vector*> &mg_dof, + MGLevelObject> &dst, + const LinearAlgebra::distributed::BlockVector &src) const { const unsigned int n_blocks = src.n_blocks(); + AssertDimension(mg_dof.size(), n_blocks); + + if (n_blocks==0) + return; + const unsigned int min_level = dst.min_level(); const unsigned int max_level = dst.max_level(); @@ -390,7 +459,11 @@ copy_to_mg (const DoFHandler &mg_dof, { const parallel::Triangulation *tria = (dynamic_cast*> - (&mg_dof.get_triangulation())); + (&(mg_dof[0]->get_triangulation()))); + for (unsigned int i=1; i*> + (&(mg_dof[0]->get_triangulation())) == tria), + ExcMessage("The DoFHandler use different Triangulations!")); for (unsigned int level = min_level; level <= max_level; ++level) { @@ -399,10 +472,10 @@ copy_to_mg (const DoFHandler &mg_dof, for (unsigned int b = 0; b < n_blocks; ++b) { LinearAlgebra::distributed::Vector &v = dst[level].block(b); - if (v.size() != mg_dof.locally_owned_mg_dofs(level).size() || - v.local_size() != mg_dof.locally_owned_mg_dofs(level).n_elements()) + if (v.size() != mg_dof[b]->locally_owned_mg_dofs(level).size() || + v.local_size() != mg_dof[b]->locally_owned_mg_dofs(level).n_elements()) { - v.reinit(mg_dof.locally_owned_mg_dofs(level), tria != nullptr ? + v.reinit(mg_dof[b]->locally_owned_mg_dofs(level), tria != nullptr ? tria->get_communicator() : MPI_COMM_SELF); collect_size = true; } @@ -421,8 +494,8 @@ copy_to_mg (const DoFHandler &mg_dof, { for (unsigned int l = min_level; l <= max_level; ++l) dst_non_block[l].reinit(dst[l].block(b)); - - matrix_free_transfer.copy_to_mg(mg_dof, dst_non_block, src.block(b)); + const unsigned int data_block = same_for_all ? 0 : b; + matrix_free_transfer_vector[data_block].copy_to_mg(*mg_dof[b], dst_non_block, src.block(b)); for (unsigned int l = min_level; l <= max_level; ++l) dst[l].block(b) = dst_non_block[l]; @@ -436,8 +509,27 @@ MGTransferBlockMatrixFree:: copy_from_mg (const DoFHandler &mg_dof, LinearAlgebra::distributed::BlockVector &dst, const MGLevelObject> &src) const +{ + AssertDimension(matrix_free_transfer_vector.size(),1); + const std::vector*> mg_dofs (dst.n_blocks(), &mg_dof); + + copy_from_mg(mg_dofs, dst, src); +} + +template +template +void +MGTransferBlockMatrixFree:: +copy_from_mg (const std::vector*> &mg_dof, + LinearAlgebra::distributed::BlockVector &dst, + const MGLevelObject> &src) const { const unsigned int n_blocks = dst.n_blocks(); + AssertDimension(mg_dof.size(), n_blocks); + + if (n_blocks==0) + return; + const unsigned int min_level = src.min_level(); const unsigned int max_level = src.max_level(); @@ -454,8 +546,8 @@ copy_from_mg (const DoFHandler &mg_dof, src_non_block[l].reinit(src[l].block(b)); src_non_block[l] = src[l].block(b); } - - matrix_free_transfer.copy_from_mg(mg_dof, dst.block(b), src_non_block); + const unsigned int data_block = same_for_all ? 0 : b; + matrix_free_transfer_vector[data_block].copy_from_mg(*mg_dof[b], dst.block(b), src_non_block); } } diff --git a/source/multigrid/mg_transfer_matrix_free.cc b/source/multigrid/mg_transfer_matrix_free.cc index 6ec3952ea5..ca1a66644f 100644 --- a/source/multigrid/mg_transfer_matrix_free.cc +++ b/source/multigrid/mg_transfer_matrix_free.cc @@ -550,9 +550,20 @@ MGTransferMatrixFree::memory_consumption() const template MGTransferBlockMatrixFree::MGTransferBlockMatrixFree (const MGConstrainedDoFs &mg_c) - : - matrix_free_transfer(mg_c) + : same_for_all (true) +{ + matrix_free_transfer_vector.emplace_back(mg_c); +} + + + +template +MGTransferBlockMatrixFree::MGTransferBlockMatrixFree +(const std::vector &mg_c) + : same_for_all (false) { + for (unsigned int i=0; i void MGTransferBlockMatrixFree::initialize_constraints (const MGConstrainedDoFs &mg_c) { - matrix_free_transfer.initialize_constraints(mg_c); + Assert(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); + + matrix_free_transfer_vector[0].initialize_constraints(mg_c); +} + + + +template +void MGTransferBlockMatrixFree::initialize_constraints +(const std::vector &mg_c) +{ + Assert(!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()); + + for (unsigned int i = 0; i::initialize_constraints template void MGTransferBlockMatrixFree::clear () { - matrix_free_transfer.clear(); + matrix_free_transfer_vector.clear(); } @@ -578,7 +612,19 @@ template void MGTransferBlockMatrixFree::build (const DoFHandler &mg_dof) { - matrix_free_transfer.build(mg_dof); + AssertDimension(matrix_free_transfer_vector.size(),1); + matrix_free_transfer_vector[0].build(mg_dof); +} + + + +template +void MGTransferBlockMatrixFree::build +(const std::vector*> &mg_dof) +{ + AssertDimension (matrix_free_transfer_vector.size(),mg_dof.size()); + for (unsigned int i=0; i LinearAlgebra::distributed::BlockVector &dst, const LinearAlgebra::distributed::BlockVector &src) const { - AssertDimension(dst.n_blocks(), src.n_blocks()); - const unsigned int n_blocks = dst.n_blocks(); + const unsigned int n_blocks = src.n_blocks(); + AssertDimension(dst.n_blocks(), n_blocks); + + if (!same_for_all) + AssertDimension (matrix_free_transfer_vector.size(), n_blocks); + for (unsigned int b = 0; b < n_blocks; ++b) - matrix_free_transfer.prolongate(to_level, dst.block(b), src.block(b)); + { + const unsigned int data_block = same_for_all ? 0 : b; + matrix_free_transfer_vector[data_block].prolongate(to_level, dst.block(b), src.block(b)); + } } @@ -603,10 +656,17 @@ void MGTransferBlockMatrixFree LinearAlgebra::distributed::BlockVector &dst, const LinearAlgebra::distributed::BlockVector &src) const { - AssertDimension(dst.n_blocks(), src.n_blocks()); - const unsigned int n_blocks = dst.n_blocks(); + const unsigned int n_blocks = src.n_blocks(); + AssertDimension(dst.n_blocks(), n_blocks); + + if (!same_for_all) + AssertDimension (matrix_free_transfer_vector.size(), n_blocks); + for (unsigned int b = 0; b < n_blocks; ++b) - matrix_free_transfer.restrict_and_add(from_level, dst.block(b), src.block(b)); + { + const unsigned int data_block = same_for_all ? 0 : b; + matrix_free_transfer_vector[data_block].restrict_and_add(from_level, dst.block(b), src.block(b)); + } } @@ -615,7 +675,10 @@ template std::size_t MGTransferBlockMatrixFree::memory_consumption() const { - return matrix_free_transfer.memory_consumption(); + 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; }