* <tt>void copy_from_mg(VectorType&)</tt> to store the result of the v-cycle
* in @p dst.
*
- * @author Guido Kanschat, 1999, 2000, 2001, 2002
+ * If VectorType is in fact a block vector and the TRANSFER object supports
+ * use of a separate DoFHandler for each block, this class also allows
+ * to be initialized with a separate DoFHandler for each block.
+ *
+ * @author Guido Kanschat, Daniel Arndt, 1999, 2000, 2001, 2002, 2017
*/
template <int dim, typename VectorType, class TRANSFER>
class PreconditionMG : public Subscriptor
Multigrid<VectorType> &mg,
const TRANSFER &transfer);
+ /**
+ * Same as above in case every component of a block vector
+ * uses its own DoFHandler.
+ */
+ PreconditionMG(const std::vector<const DoFHandler<dim>*> &dof_handler,
+ Multigrid<VectorType> &mg,
+ const TRANSFER &transfer);
+
/**
* Dummy function needed by other classes.
*/
const OtherVectorType &src) const;
/**
- * Return the partitioning of the range space of this preconditioner, i.e., the partitioning of the vectors that are result from matrix-vector products.
+ * Return the partitioning of the range space of this preconditioner, i.e.,
+ * the partitioning of the vectors that are result from matrix-vector products.
+ * By default, the respective information for the first DoFHandler object
+ * are returned.
*/
- IndexSet locally_owned_range_indices() const;
+ IndexSet locally_owned_range_indices(const unsigned int block=0) const;
/**
- * Return the partitioning of the domain space of this preconditioner, i.e., the partitioning of the vectors this matrix has to be multiplied with.
+ * Return the partitioning of the domain space of this preconditioner, i.e.,
+ * the partitioning of the vectors this matrix has to be multiplied with.
+ * By default, the respective information for the first DoFHandler object
+ * are returned.
*/
- IndexSet locally_owned_domain_indices() const;
+ IndexSet locally_owned_domain_indices(const unsigned int block=0) const;
/**
* Return the MPI communicator object in use with this preconditioner.
/**
* Associated @p DoFHandler.
*/
- SmartPointer<const DoFHandler<dim>,PreconditionMG<dim,VectorType,TRANSFER> > dof_handler;
+ std::vector<SmartPointer<const DoFHandler<dim>,PreconditionMG<dim,VectorType,TRANSFER> > > dof_handler_vector;
+
+ /**
+ * Storage for the pointers to the DoFHandler objects
+ * without SmartPointer wrapper.
+ */
+ std::vector<const DoFHandler<dim>*> dof_handler_vector_raw;
/**
* The multigrid object.
* Object for grid tranfer.
*/
SmartPointer<const TRANSFER,PreconditionMG<dim,VectorType,TRANSFER> > transfer;
+
+ /**
+ * Flag to indicate if the object is initialized with a single DoFHandler
+ * or with one for each block.
+ */
+ const bool uses_dof_handler_vector;
};
/*@}*/
/* --------------------------- inline functions --------------------- */
+namespace internal
+{
+ namespace PreconditionMG
+ {
+ template <int dim, typename VectorType, class TRANSFER, typename OtherVectorType>
+ typename std::enable_if<TRANSFER::supports_dof_handler_vector>::type
+ vmult(const std::vector<const dealii::DoFHandler<dim>*> &dof_handler_vector,
+ dealii::Multigrid<VectorType> &multigrid,
+ const TRANSFER &transfer,
+ OtherVectorType &dst,
+ const OtherVectorType &src,
+ const bool uses_dof_handler_vector, int)
+ {
+ if (uses_dof_handler_vector)
+ transfer.copy_to_mg(dof_handler_vector,
+ multigrid.defect,
+ src);
+ else
+ transfer.copy_to_mg(*dof_handler_vector[0],
+ multigrid.defect,
+ src);
+
+ multigrid.cycle();
+ if (uses_dof_handler_vector)
+ transfer.copy_from_mg(dof_handler_vector,
+ dst,
+ multigrid.solution);
+ else
+ transfer.copy_from_mg(*dof_handler_vector[0],
+ dst,
+ multigrid.solution);
+ }
+
+ template <int dim, typename VectorType, class TRANSFER, typename OtherVectorType>
+ void
+ vmult(const std::vector<const dealii::DoFHandler<dim>*> &dof_handler_vector,
+ dealii::Multigrid<VectorType> &multigrid,
+ const TRANSFER &transfer,
+ OtherVectorType &dst,
+ const OtherVectorType &src,
+ const bool uses_dof_handler_vector,...)
+ {
+ (void) uses_dof_handler_vector;
+ Assert (!uses_dof_handler_vector, ExcInternalError());
+ transfer.copy_to_mg(*dof_handler_vector[0],
+ multigrid.defect,
+ src);
+ multigrid.cycle();
+ transfer.copy_from_mg(*dof_handler_vector[0],
+ dst,
+ multigrid.solution);
+ }
+
+ template <int dim, typename VectorType, class TRANSFER, typename OtherVectorType>
+ typename std::enable_if<TRANSFER::supports_dof_handler_vector>::type
+ vmult_add(const std::vector<const dealii::DoFHandler<dim>*> &dof_handler_vector,
+ dealii::Multigrid<VectorType> &multigrid,
+ const TRANSFER &transfer,
+ OtherVectorType &dst,
+ const OtherVectorType &src,
+ const bool uses_dof_handler_vector, int)
+ {
+ if (uses_dof_handler_vector)
+ transfer.copy_to_mg(dof_handler_vector,
+ multigrid.defect,
+ src);
+ else
+ transfer.copy_to_mg(*dof_handler_vector[0],
+ multigrid.defect,
+ src);
+
+ multigrid.cycle();
+ if (uses_dof_handler_vector)
+ transfer.copy_from_mg_add(dof_handler_vector,
+ dst,
+ multigrid.solution);
+ else
+ transfer.copy_from_mg_add(*dof_handler_vector[0],
+ dst,
+ multigrid.solution);
+ }
+
+ template <int dim, typename VectorType, class TRANSFER, typename OtherVectorType>
+ void
+ vmult_add(const std::vector<const dealii::DoFHandler<dim>*> &dof_handler_vector,
+ dealii::Multigrid<VectorType> &multigrid,
+ const TRANSFER &transfer,
+ OtherVectorType &dst,
+ const OtherVectorType &src,
+ const bool uses_dof_handler_vector,...)
+ {
+ (void) uses_dof_handler_vector;
+ Assert (!uses_dof_handler_vector, ExcInternalError());
+ transfer.copy_to_mg(*dof_handler_vector[0],
+ multigrid.defect,
+ src);
+ multigrid.cycle();
+ transfer.copy_from_mg_add(*dof_handler_vector[0],
+ dst,
+ multigrid.solution);
+ }
+ }
+}
+
template <int dim, typename VectorType, class TRANSFER>
PreconditionMG<dim, VectorType, TRANSFER>
::PreconditionMG(const DoFHandler<dim> &dof_handler,
Multigrid<VectorType> &mg,
const TRANSFER &transfer)
:
- dof_handler(&dof_handler),
+ dof_handler_vector(1,&dof_handler),
+ dof_handler_vector_raw(1,&dof_handler),
multigrid(&mg),
- transfer(&transfer)
+ transfer(&transfer),
+ uses_dof_handler_vector(false)
{}
+template <int dim, typename VectorType, class TRANSFER>
+PreconditionMG<dim, VectorType, TRANSFER>
+::PreconditionMG(const std::vector<const DoFHandler<dim>*> &dof_handler,
+ Multigrid<VectorType> &mg,
+ const TRANSFER &transfer)
+ :
+ dof_handler_vector(dof_handler.size()),
+ dof_handler_vector_raw(dof_handler.size()),
+ multigrid(&mg),
+ transfer(&transfer),
+ uses_dof_handler_vector(true)
+{
+ for (unsigned int i = 0; i< dof_handler.size() ; ++i)
+ {
+ dof_handler_vector[i] = dof_handler[i];
+ dof_handler_vector_raw[i] = dof_handler[i];
+ }
+}
+
template <int dim, typename VectorType, class TRANSFER>
inline bool
PreconditionMG<dim, VectorType, TRANSFER>::empty () const
(OtherVectorType &dst,
const OtherVectorType &src) const
{
- transfer->copy_to_mg(*dof_handler,
- multigrid->defect,
- src);
- multigrid->cycle();
-
- transfer->copy_from_mg(*dof_handler,
- dst,
- multigrid->solution);
+ internal::PreconditionMG::vmult(dof_handler_vector_raw,*multigrid,*transfer,
+ dst,src,uses_dof_handler_vector,0);
}
template <int dim, typename VectorType, class TRANSFER>
IndexSet
-PreconditionMG<dim, VectorType, TRANSFER>::locally_owned_range_indices() const
+PreconditionMG<dim, VectorType, TRANSFER>::locally_owned_range_indices(const unsigned int block) const
{
- return dof_handler->locally_owned_dofs();
+ AssertIndexRange(block,dof_handler_vector.size());
+ return dof_handler_vector[block]->locally_owned_dofs();
}
template <int dim, typename VectorType, class TRANSFER>
IndexSet
-PreconditionMG<dim, VectorType, TRANSFER>::locally_owned_domain_indices() const
+PreconditionMG<dim, VectorType, TRANSFER>::locally_owned_domain_indices(const unsigned int block) const
{
- return dof_handler->locally_owned_dofs();
+ AssertIndexRange(block,dof_handler_vector.size());
+ return dof_handler_vector[block]->locally_owned_dofs();
}
{
// currently parallel GMG works with distributed Triangulation only,
// so it should be a safe bet to use it to query MPI communicator:
- const Triangulation<dim> &tria = dof_handler->get_triangulation();
+ const Triangulation<dim> &tria = dof_handler_vector[0]->get_triangulation();
const parallel::distributed::Triangulation<dim> *ptria = dynamic_cast<const parallel::distributed::Triangulation<dim> *>(&tria);
Assert (ptria != nullptr, ExcInternalError());
return ptria->get_communicator ();
(OtherVectorType &dst,
const OtherVectorType &src) const
{
- transfer->copy_to_mg(*dof_handler,
- multigrid->defect,
- src);
- multigrid->cycle();
- transfer->copy_from_mg_add(*dof_handler,
- dst,
- multigrid->solution);
+ internal::PreconditionMG::vmult_add(dof_handler_vector_raw,*multigrid,*transfer,
+ dst,src,uses_dof_handler_vector,0);
}