]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Allow to use MGTransferBlockMatrixFree with a separate DoFHandler for each block
authorDaniel Arndt <daniel.arndt@iwr.uni-heidelberg.de>
Sun, 1 Oct 2017 20:53:10 +0000 (22:53 +0200)
committerDaniel Arndt <daniel.arndt@iwr.uni-heidelberg.de>
Mon, 2 Oct 2017 08:17:34 +0000 (10:17 +0200)
include/deal.II/multigrid/mg_transfer_matrix_free.h
source/multigrid/mg_transfer_matrix_free.cc

index ddb94109c9cf6c74468990ef770351f48abccf41..84936b99fbc036d68027835ae2e7ca03fd0bd2c4 100644 (file)
@@ -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 <int dim, typename Number>
@@ -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<MGConstrainedDoFs> &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<MGConstrainedDoFs> &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<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
@@ -343,6 +359,15 @@ public:
               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.
    */
@@ -352,17 +377,38 @@ public:
                 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;
 };
 
 
@@ -379,8 +425,31 @@ MGTransferBlockMatrixFree<dim,Number>::
 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();
 
@@ -390,7 +459,11 @@ copy_to_mg (const DoFHandler<dim,spacedim> &mg_dof,
   {
     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)
       {
@@ -399,10 +472,10 @@ copy_to_mg (const DoFHandler<dim,spacedim> &mg_dof,
         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;
               }
@@ -421,8 +494,8 @@ copy_to_mg (const DoFHandler<dim,spacedim> &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<dim,Number>::
 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();
 
@@ -454,8 +546,8 @@ copy_from_mg (const DoFHandler<dim,spacedim>                        &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);
     }
 }
 
index 6ec3952ea540637efcc2f71334bfd42acf019e42..ca1a66644fd84de698f973222009b18cf9d85c42 100644 (file)
@@ -550,9 +550,20 @@ MGTransferMatrixFree<dim,Number>::memory_consumption() const
 
 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]);
 }
 
 
@@ -561,7 +572,30 @@ template <int dim, typename Number>
 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]);
 }
 
 
@@ -569,7 +603,7 @@ void MGTransferBlockMatrixFree<dim,Number>::initialize_constraints
 template <int dim, typename Number>
 void MGTransferBlockMatrixFree<dim,Number>::clear ()
 {
-  matrix_free_transfer.clear();
+  matrix_free_transfer_vector.clear();
 }
 
 
@@ -578,7 +612,19 @@ template <int dim, typename Number>
 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]);
 }
 
 
@@ -589,10 +635,17 @@ void MGTransferBlockMatrixFree<dim,Number>
               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));
+    }
 }
 
 
@@ -603,10 +656,17 @@ void MGTransferBlockMatrixFree<dim,Number>
                     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));
+    }
 }
 
 
@@ -615,7 +675,10 @@ template <int dim, typename Number>
 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;
 }
 
 

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.