*
* 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 <int dim, typename Number>
*/
MGTransferBlockMatrixFree (const MGConstrainedDoFs &mg_constrained_dofs);
+ /**
+ * Same as above for the case that each block has its own DoFHandler.
+ */
+ MGTransferBlockMatrixFree (const std::vector<MGConstrainedDoFs> &mg_constrained_dofs);
+
/**
* Destructor.
*/
*/
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<MGConstrainedDoFs> &mg_constrained_dofs);
+
/**
* Reset the object to the state it had right after the default constructor.
*/
*/
void build (const DoFHandler<dim,dim> &mg_dof);
+ /**
+ * Same as above for the case that each block has its own DoFHandler.
+ */
+ void build (const std::vector<const DoFHandler<dim,dim>*> &mg_dof);
+
/**
* Prolongate a vector from level <tt>to_level-1</tt> to level
* <tt>to_level</tt> using the embedding matrices of the underlying finite
MGLevelObject<LinearAlgebra::distributed::BlockVector<Number>> &dst,
const LinearAlgebra::distributed::BlockVector<Number2> &src) const;
+ /**
+ * Same as above for the case that each block has its own DoFHandler.
+ */
+ template <typename Number2, int spacedim>
+ void
+ copy_to_mg (const std::vector<const DoFHandler<dim,spacedim>*> &mg_dof,
+ MGLevelObject<LinearAlgebra::distributed::BlockVector<Number>> &dst,
+ const LinearAlgebra::distributed::BlockVector<Number2> &src) const;
+
/**
* Transfer from multi-level block-vector to normal vector.
*/
LinearAlgebra::distributed::BlockVector<Number2> &dst,
const MGLevelObject<LinearAlgebra::distributed::BlockVector<Number>> &src) const;
+ /**
+ * Same as above for the case that each block has its own DoFHandler.
+ */
+ template <typename Number2, int spacedim>
+ void
+ copy_from_mg (const std::vector<const DoFHandler<dim,spacedim>*> &mg_dof,
+ LinearAlgebra::distributed::BlockVector<Number2> &dst,
+ const MGLevelObject<LinearAlgebra::distributed::BlockVector<Number>> &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<MGTransferMatrixFree<dim,Number>> 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<dim,Number> matrix_free_transfer;
+ const bool same_for_all;
};
copy_to_mg (const DoFHandler<dim,spacedim> &mg_dof,
MGLevelObject<LinearAlgebra::distributed::BlockVector<Number>> &dst,
const LinearAlgebra::distributed::BlockVector<Number2> &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<const DoFHandler<dim, spacedim>*> mg_dofs (src.n_blocks(),&mg_dof);
+
+ copy_to_mg(mg_dofs, dst, src);
+}
+
+template <int dim, typename Number>
+template <typename Number2, int spacedim>
+void
+MGTransferBlockMatrixFree<dim,Number>::
+copy_to_mg (const std::vector<const DoFHandler<dim,spacedim>*> &mg_dof,
+ MGLevelObject<LinearAlgebra::distributed::BlockVector<Number>> &dst,
+ const LinearAlgebra::distributed::BlockVector<Number2> &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();
{
const parallel::Triangulation<dim,spacedim> *tria =
(dynamic_cast<const parallel::Triangulation<dim,spacedim>*>
- (&mg_dof.get_triangulation()));
+ (&(mg_dof[0]->get_triangulation())));
+ for (unsigned int i=1; i<n_blocks; ++i)
+ AssertThrow((dynamic_cast<const parallel::Triangulation<dim,spacedim>*>
+ (&(mg_dof[0]->get_triangulation())) == tria),
+ ExcMessage("The DoFHandler use different Triangulations!"));
for (unsigned int level = min_level; level <= max_level; ++level)
{
for (unsigned int b = 0; b < n_blocks; ++b)
{
LinearAlgebra::distributed::Vector<Number> &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;
}
{
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];
copy_from_mg (const DoFHandler<dim,spacedim> &mg_dof,
LinearAlgebra::distributed::BlockVector<Number2> &dst,
const MGLevelObject<LinearAlgebra::distributed::BlockVector<Number>> &src) const
+{
+ AssertDimension(matrix_free_transfer_vector.size(),1);
+ const std::vector<const DoFHandler<dim, spacedim>*> mg_dofs (dst.n_blocks(), &mg_dof);
+
+ copy_from_mg(mg_dofs, dst, src);
+}
+
+template <int dim, typename Number>
+template <typename Number2, int spacedim>
+void
+MGTransferBlockMatrixFree<dim,Number>::
+copy_from_mg (const std::vector<const DoFHandler<dim,spacedim>*> &mg_dof,
+ LinearAlgebra::distributed::BlockVector<Number2> &dst,
+ const MGLevelObject<LinearAlgebra::distributed::BlockVector<Number>> &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();
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);
}
}
template <int dim, typename Number>
MGTransferBlockMatrixFree<dim,Number>::MGTransferBlockMatrixFree (const MGConstrainedDoFs &mg_c)
- :
- matrix_free_transfer(mg_c)
+ : same_for_all (true)
+{
+ matrix_free_transfer_vector.emplace_back(mg_c);
+}
+
+
+
+template <int dim, typename Number>
+MGTransferBlockMatrixFree<dim,Number>::MGTransferBlockMatrixFree
+(const std::vector<MGConstrainedDoFs> &mg_c)
+ : same_for_all (false)
{
+ for (unsigned int i=0; i<mg_c.size(); ++i)
+ matrix_free_transfer_vector.emplace_back(mg_c[i]);
}
void MGTransferBlockMatrixFree<dim,Number>::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 <int dim, typename Number>
+void MGTransferBlockMatrixFree<dim,Number>::initialize_constraints
+(const std::vector<MGConstrainedDoFs> &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<mg_c.size(); ++i)
+ matrix_free_transfer_vector[i].initialize_constraints(mg_c[i]);
}
template <int dim, typename Number>
void MGTransferBlockMatrixFree<dim,Number>::clear ()
{
- matrix_free_transfer.clear();
+ matrix_free_transfer_vector.clear();
}
void MGTransferBlockMatrixFree<dim,Number>::build
(const DoFHandler<dim,dim> &mg_dof)
{
- matrix_free_transfer.build(mg_dof);
+ AssertDimension(matrix_free_transfer_vector.size(),1);
+ matrix_free_transfer_vector[0].build(mg_dof);
+}
+
+
+
+template <int dim, typename Number>
+void MGTransferBlockMatrixFree<dim,Number>::build
+(const std::vector<const DoFHandler<dim,dim>*> &mg_dof)
+{
+ AssertDimension (matrix_free_transfer_vector.size(),mg_dof.size());
+ for (unsigned int i=0; i<mg_dof.size(); ++i)
+ matrix_free_transfer_vector[i].build(*mg_dof[i]);
}
LinearAlgebra::distributed::BlockVector<Number> &dst,
const LinearAlgebra::distributed::BlockVector<Number> &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));
+ }
}
LinearAlgebra::distributed::BlockVector<Number> &dst,
const LinearAlgebra::distributed::BlockVector<Number> &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));
+ }
}
std::size_t
MGTransferBlockMatrixFree<dim,Number>::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;
}