]> https://gitweb.dealii.org/ - dealii.git/commitdiff
MGTwoLevelTransfer: improve setup of indices
authorPeter Munch <peterrmuench@gmail.com>
Thu, 7 Dec 2023 15:07:58 +0000 (16:07 +0100)
committerPeter Munch <peterrmuench@gmail.com>
Sun, 10 Dec 2023 15:15:19 +0000 (16:15 +0100)
include/deal.II/matrix_free/constraint_info.h
include/deal.II/multigrid/mg_transfer_global_coarsening.templates.h

index d3c7d7727c6bab61895bac625c1de72573143050..184587c5a66142e9dcb36ba64f915b45d8228d39 100644 (file)
@@ -85,15 +85,29 @@ namespace internal
     class ConstraintInfo
     {
     public:
+      /**
+       * Constructor.
+       */
+      ConstraintInfo();
+
+      /**
+       * Set locally owned indices to accelerate the conversion of global to
+       * local indices.
+       */
+      void
+      set_locally_owned_indices(const IndexSet &locally_owned_indices);
+
       /**
        * Version 1: indices are extracted from DoFCellAccessor and
        * constraints are resolved with the help of AffineConstraints.
        */
+      template <typename T = unsigned int>
       void
       reinit(const DoFHandler<dim> &dof_handler,
              const unsigned int     n_cells,
              const bool             use_fast_hanging_node_algorithm = true);
 
+      template <typename T = unsigned int>
       void
       read_dof_indices(
         const unsigned int                                    cell_no,
@@ -106,18 +120,25 @@ namespace internal
       /**
        * Version 2: no constraints, indices are user-provided.
        */
+      template <typename T = unsigned int>
       void
       reinit(const unsigned int n_cells);
 
+      template <typename T = unsigned int>
       void
       read_dof_indices(
         const unsigned int                                        cell_no,
         const std::vector<types::global_dof_index>               &dof_indices,
         const std::shared_ptr<const Utilities::MPI::Partitioner> &partitioner);
 
+      template <typename T = unsigned int>
       void
       finalize();
 
+      template <typename T = unsigned int>
+      std::shared_ptr<const Utilities::MPI::Partitioner>
+      finalize(const MPI_Comm comm);
+
       template <typename T, typename VectorType>
       void
       read_write_operation(const T           &operation,
@@ -146,12 +167,19 @@ namespace internal
       ConstraintValues<double>               constraint_values;
       std::vector<std::vector<unsigned int>> dof_indices_per_cell;
       std::vector<std::vector<unsigned int>> plain_dof_indices_per_cell;
+      std::vector<std::vector<types::global_dof_index>>
+        dof_indices_per_cell_long;
+      std::vector<std::vector<types::global_dof_index>>
+        plain_dof_indices_per_cell_long;
       std::vector<std::vector<std::pair<unsigned short, unsigned short>>>
         constraint_indicator_per_cell;
 
       std::unique_ptr<HangingNodes<dim>>     hanging_nodes;
       std::vector<std::vector<unsigned int>> lexicographic_numbering;
 
+      IndexSet locally_owned_indices;
+      std::pair<types::global_dof_index, types::global_dof_index> local_range;
+
     public:
       // for read_write_operation()
       std::vector<unsigned int> dof_indices;
@@ -176,6 +204,35 @@ namespace internal
 
       inline const typename Number::value_type *
       constraint_pool_end(const unsigned int row) const;
+
+      // TODO
+      std::vector<types::global_dof_index> local_dof_indices;
+      std::vector<types::global_dof_index> local_dof_indices_lex;
+      std::vector<ConstraintKinds>         mask;
+
+      template <typename T>
+      std::vector<std::vector<T>> &
+      get_dof_indices_per_cell()
+      {
+        if constexpr (std::is_same_v<T, unsigned int>)
+          return dof_indices_per_cell;
+        else if constexpr (std::is_same_v<T, types::global_dof_index>)
+          return dof_indices_per_cell_long;
+        else
+          return {};
+      }
+
+      template <typename T>
+      std::vector<std::vector<T>> &
+      get_plain_dof_indices_per_cell()
+      {
+        if constexpr (std::is_same_v<T, unsigned int>)
+          return plain_dof_indices_per_cell;
+        else if constexpr (std::is_same_v<T, types::global_dof_index>)
+          return plain_dof_indices_per_cell_long;
+        else
+          return {};
+      }
     };
 
 
@@ -249,14 +306,39 @@ namespace internal
 
 
     template <int dim, typename Number>
+    ConstraintInfo<dim, Number>::ConstraintInfo()
+      : local_range({0, 0})
+    {}
+
+
+
+    template <int dim, typename Number>
+    void
+    ConstraintInfo<dim, Number>::set_locally_owned_indices(
+      const IndexSet &locally_owned_indices)
+    {
+      this->locally_owned_indices = locally_owned_indices;
+
+      if (locally_owned_indices.is_empty())
+        local_range = {0, 0};
+      else
+        local_range = {locally_owned_indices.nth_index_in_set(0),
+                       locally_owned_indices.nth_index_in_set(0) +
+                         locally_owned_indices.n_elements()};
+    }
+
+
+
+    template <int dim, typename Number>
+    template <typename T>
     inline void
     ConstraintInfo<dim, Number>::reinit(
       const DoFHandler<dim> &dof_handler,
       const unsigned int     n_cells,
       const bool             use_fast_hanging_node_algorithm)
     {
-      this->dof_indices_per_cell.resize(n_cells);
-      this->plain_dof_indices_per_cell.resize(n_cells);
+      this->template get_dof_indices_per_cell<T>().resize(n_cells);
+      this->template get_plain_dof_indices_per_cell<T>().resize(n_cells);
       this->constraint_indicator_per_cell.resize(n_cells);
 
       // note: has_hanging_nodes() is a global operatrion
@@ -299,17 +381,19 @@ namespace internal
 
 
     template <int dim, typename Number>
+    template <typename T>
     inline void
     ConstraintInfo<dim, Number>::reinit(const unsigned int n_cells)
     {
-      this->dof_indices_per_cell.resize(n_cells);
-      this->plain_dof_indices_per_cell.resize(0);
+      this->template get_dof_indices_per_cell<T>().resize(n_cells);
+      this->template get_plain_dof_indices_per_cell<T>().resize(0);
       this->constraint_indicator_per_cell.resize(n_cells);
     }
 
 
 
     template <int dim, typename Number>
+    template <typename T>
     inline void
     ConstraintInfo<dim, Number>::read_dof_indices(
       const unsigned int                                            cell_no,
@@ -318,10 +402,8 @@ namespace internal
       const dealii::AffineConstraints<typename Number::value_type> &constraints,
       const std::shared_ptr<const Utilities::MPI::Partitioner>     &partitioner)
     {
-      std::vector<types::global_dof_index> local_dof_indices(
-        cell->get_fe().n_dofs_per_cell());
-      std::vector<types::global_dof_index> local_dof_indices_lex(
-        cell->get_fe().n_dofs_per_cell());
+      local_dof_indices.resize(cell->get_fe().n_dofs_per_cell());
+      local_dof_indices_lex.resize(cell->get_fe().n_dofs_per_cell());
 
       if (mg_level == numbers::invalid_unsigned_int)
         cell->get_dof_indices(local_dof_indices);
@@ -345,24 +427,34 @@ namespace internal
       std::pair<unsigned short, unsigned short> constraint_iterator(0, 0);
 
       AssertIndexRange(cell_no, this->constraint_indicator_per_cell.size());
-      AssertIndexRange(cell_no, this->dof_indices_per_cell.size());
-      AssertIndexRange(cell_no, this->plain_dof_indices_per_cell.size());
-      AssertIndexRange(cell_no, this->plain_dof_indices_per_cell.size());
+      AssertIndexRange(cell_no,
+                       this->template get_dof_indices_per_cell<T>().size());
+      AssertIndexRange(
+        cell_no, this->template get_plain_dof_indices_per_cell<T>().size());
 
       auto &constraint_indicator = this->constraint_indicator_per_cell[cell_no];
-      auto &dof_indices          = this->dof_indices_per_cell[cell_no];
-      auto &plain_dof_indices    = this->plain_dof_indices_per_cell[cell_no];
+      auto &dof_indices = this->template get_dof_indices_per_cell<T>()[cell_no];
+      auto &plain_dof_indices =
+        this->template get_plain_dof_indices_per_cell<T>()[cell_no];
 
       AssertDimension(constraint_indicator_per_cell[cell_no].size(), 0);
-      AssertDimension(dof_indices_per_cell[cell_no].size(), 0);
-      AssertDimension(plain_dof_indices_per_cell[cell_no].size(), 0);
+      AssertDimension(
+        this->template get_dof_indices_per_cell<T>()[cell_no].size(), 0);
+      AssertDimension(
+        this->template get_plain_dof_indices_per_cell<T>()[cell_no].size(), 0);
 
       const auto global_to_local =
-        [&](const types::global_dof_index global_index) -> unsigned int {
+        [&](const types::global_dof_index global_index) -> T {
         if (partitioner)
           return partitioner->global_to_local(global_index);
         else
-          return global_index;
+          {
+            if (local_range.first <= global_index &&
+                global_index < local_range.second)
+              return global_index - local_range.first;
+            else
+              return global_index + (local_range.second - local_range.first);
+          }
       };
 
       // plain indices
@@ -375,7 +467,9 @@ namespace internal
           AssertIndexRange(cell_no, this->hanging_node_constraint_masks.size());
           AssertIndexRange(cell_no, this->active_fe_indices.size());
 
-          std::vector<ConstraintKinds> mask(cell->get_fe().n_components());
+          mask.assign(
+            cell->get_fe().n_components(),
+            internal::MatrixFreeFunctions::ConstraintKinds::unconstrained);
           hanging_nodes->setup_constraints(
             cell, {}, lexicographic_numbering, local_dof_indices_lex, mask);
 
@@ -437,6 +531,7 @@ namespace internal
 
 
     template <int dim, typename Number>
+    template <typename T>
     inline void
     ConstraintInfo<dim, Number>::read_dof_indices(
       const unsigned int                          cell_no,
@@ -444,17 +539,23 @@ namespace internal
       const std::shared_ptr<const Utilities::MPI::Partitioner> &partitioner)
     {
       const auto global_to_local =
-        [&](const types::global_dof_index global_index) -> unsigned int {
+        [&](const types::global_dof_index global_index) -> T {
         if (partitioner)
           return partitioner->global_to_local(global_index);
         else
-          return global_index;
+          {
+            if (local_range.first <= global_index &&
+                global_index < local_range.second)
+              return global_index - local_range.first;
+            else
+              return global_index + (local_range.second - local_range.first);
+          }
       };
 
       std::pair<unsigned short, unsigned short> constraint_iterator(0, 0);
 
       auto &constraint_indicator = this->constraint_indicator_per_cell[cell_no];
-      auto &dof_indices          = this->dof_indices_per_cell[cell_no];
+      auto &dof_indices = this->template get_dof_indices_per_cell<T>()[cell_no];
 
       for (const auto current_dof : local_dof_indices_lex)
         {
@@ -488,9 +589,15 @@ namespace internal
 
 
     template <int dim, typename Number>
+    template <typename T>
     inline void
     ConstraintInfo<dim, Number>::finalize()
     {
+      Assert((local_range ==
+              std::pair<types::global_dof_index, types::global_dof_index>{0,
+                                                                          0}),
+             ExcNotImplemented());
+
       this->dof_indices          = {};
       this->plain_dof_indices    = {};
       this->constraint_indicator = {};
@@ -498,17 +605,20 @@ namespace internal
       this->row_starts = {};
       this->row_starts.emplace_back(0, 0);
 
-      if (plain_dof_indices_per_cell.empty() == false)
+      if (this->template get_plain_dof_indices_per_cell<T>().empty() == false)
         {
           this->row_starts_plain_indices = {};
           this->row_starts_plain_indices.emplace_back(0);
         }
 
-      for (unsigned int i = 0; i < dof_indices_per_cell.size(); ++i)
+      for (unsigned int i = 0;
+           i < this->template get_dof_indices_per_cell<T>().size();
+           ++i)
         {
-          this->dof_indices.insert(this->dof_indices.end(),
-                                   dof_indices_per_cell[i].begin(),
-                                   dof_indices_per_cell[i].end());
+          this->dof_indices.insert(
+            this->dof_indices.end(),
+            this->template get_dof_indices_per_cell<T>()[i].begin(),
+            this->template get_dof_indices_per_cell<T>()[i].end());
           this->constraint_indicator.insert(
             this->constraint_indicator.end(),
             constraint_indicator_per_cell[i].begin(),
@@ -517,12 +627,13 @@ namespace internal
           this->row_starts.emplace_back(this->dof_indices.size(),
                                         this->constraint_indicator.size());
 
-          if (plain_dof_indices_per_cell.empty() == false)
+          if (this->template get_plain_dof_indices_per_cell<T>().empty() ==
+              false)
             {
               this->plain_dof_indices.insert(
                 this->plain_dof_indices.end(),
-                plain_dof_indices_per_cell[i].begin(),
-                plain_dof_indices_per_cell[i].end());
+                this->template get_plain_dof_indices_per_cell<T>()[i].begin(),
+                this->template get_plain_dof_indices_per_cell<T>()[i].end());
 
               this->row_starts_plain_indices.emplace_back(
                 this->plain_dof_indices.size());
@@ -556,7 +667,7 @@ namespace internal
 
       AssertDimension(constraint_pool_data.size(), length);
 
-      dof_indices_per_cell.clear();
+      this->template get_dof_indices_per_cell<T>().clear();
       constraint_indicator_per_cell.clear();
 
       if (hanging_nodes &&
@@ -569,6 +680,150 @@ namespace internal
     }
 
 
+    template <int dim, typename Number>
+    template <typename T>
+    inline std::shared_ptr<const Utilities::MPI::Partitioner>
+    ConstraintInfo<dim, Number>::finalize(const MPI_Comm comm)
+    {
+      this->dof_indices.clear();
+      this->plain_dof_indices.clear();
+      this->constraint_indicator.clear();
+
+      this->row_starts.clear();
+      this->row_starts.reserve(
+        this->template get_dof_indices_per_cell<T>().size());
+      this->row_starts.emplace_back(0, 0);
+
+      if (this->template get_plain_dof_indices_per_cell<T>().empty() == false)
+        {
+          this->row_starts_plain_indices.clear();
+          this->row_starts_plain_indices.reserve(
+            this->template get_dof_indices_per_cell<T>().size());
+          this->row_starts_plain_indices.emplace_back(0);
+        }
+
+      std::vector<types::global_dof_index>  ghost_dofs;
+      std::pair<unsigned int, unsigned int> counts = {0, 0};
+
+      for (unsigned int i = 0;
+           i < this->template get_dof_indices_per_cell<T>().size();
+           ++i)
+        {
+          counts.first +=
+            this->template get_dof_indices_per_cell<T>()[i].size();
+
+          for (const auto &j : this->template get_dof_indices_per_cell<T>()[i])
+            if (j >= (local_range.second - local_range.first))
+              ghost_dofs.push_back(j -
+                                   (local_range.second - local_range.first));
+
+          if (this->template get_plain_dof_indices_per_cell<T>().empty() ==
+              false)
+            {
+              counts.second +=
+                this->template get_plain_dof_indices_per_cell<T>()[i].size();
+
+              for (const auto &j :
+                   this->template get_plain_dof_indices_per_cell<T>()[i])
+                if (j >= (local_range.second - local_range.first))
+                  ghost_dofs.push_back(
+                    j - (local_range.second - local_range.first));
+            }
+        }
+
+      std::sort(ghost_dofs.begin(), ghost_dofs.end());
+      ghost_dofs.erase(std::unique(ghost_dofs.begin(), ghost_dofs.end()),
+                       ghost_dofs.end());
+
+      IndexSet locally_relevant_dofs(locally_owned_indices.size());
+      locally_relevant_dofs.add_indices(ghost_dofs.begin(), ghost_dofs.end());
+
+      const auto partitioner =
+        std::make_shared<Utilities::MPI::Partitioner>(locally_owned_indices,
+                                                      locally_relevant_dofs,
+                                                      comm);
+
+      this->dof_indices.reserve(counts.first);
+      this->plain_dof_indices.reserve(counts.second);
+
+      for (unsigned int i = 0;
+           i < this->template get_dof_indices_per_cell<T>().size();
+           ++i)
+        {
+          for (const auto &j : this->template get_dof_indices_per_cell<T>()[i])
+            if (j < (local_range.second - local_range.first))
+              this->dof_indices.push_back(j);
+            else
+              this->dof_indices.push_back(partitioner->global_to_local(
+                j - (local_range.second - local_range.first)));
+
+          this->constraint_indicator.insert(
+            this->constraint_indicator.end(),
+            constraint_indicator_per_cell[i].begin(),
+            constraint_indicator_per_cell[i].end());
+
+          this->row_starts.emplace_back(this->dof_indices.size(),
+                                        this->constraint_indicator.size());
+
+          if (this->template get_plain_dof_indices_per_cell<T>().empty() ==
+              false)
+            {
+              for (const auto &j :
+                   this->template get_plain_dof_indices_per_cell<T>()[i])
+                if (j < (local_range.second - local_range.first))
+                  this->plain_dof_indices.push_back(j);
+                else
+                  this->plain_dof_indices.push_back(
+                    partitioner->global_to_local(
+                      j - (local_range.second - local_range.first)));
+
+              this->row_starts_plain_indices.emplace_back(
+                this->plain_dof_indices.size());
+            }
+        }
+
+      std::vector<const std::vector<double> *> constraints(
+        constraint_values.constraints.size());
+      unsigned int length = 0;
+      for (const auto &it : constraint_values.constraints)
+        {
+          AssertIndexRange(it.second, constraints.size());
+          constraints[it.second] = &it.first;
+          length += it.first.size();
+        }
+
+      constraint_pool_data.clear();
+      constraint_pool_data.reserve(length);
+      constraint_pool_row_index.reserve(constraint_values.constraints.size() +
+                                        1);
+      constraint_pool_row_index.resize(1, 0);
+
+      for (const auto &constraint : constraints)
+        {
+          Assert(constraint != nullptr, ExcInternalError());
+          constraint_pool_data.insert(constraint_pool_data.end(),
+                                      constraint->begin(),
+                                      constraint->end());
+          constraint_pool_row_index.push_back(constraint_pool_data.size());
+        }
+
+      AssertDimension(constraint_pool_data.size(), length);
+
+      this->template get_dof_indices_per_cell<T>().clear();
+      constraint_indicator_per_cell.clear();
+
+      if (hanging_nodes &&
+          std::all_of(hanging_node_constraint_masks.begin(),
+                      hanging_node_constraint_masks.end(),
+                      [](const auto i) {
+                        return i == unconstrained_compressed_constraint_kind;
+                      }))
+        hanging_node_constraint_masks.clear();
+
+      return partitioner;
+    }
+
+
 
     template <int dim, typename Number>
     template <typename T, typename VectorType>
@@ -741,6 +996,9 @@ namespace internal
       size += MemoryConsumption::memory_consumption(constraint_values);
       size += MemoryConsumption::memory_consumption(dof_indices_per_cell);
       size += MemoryConsumption::memory_consumption(plain_dof_indices_per_cell);
+      size += MemoryConsumption::memory_consumption(dof_indices_per_cell_long);
+      size +=
+        MemoryConsumption::memory_consumption(plain_dof_indices_per_cell_long);
       size +=
         MemoryConsumption::memory_consumption(constraint_indicator_per_cell);
 
index f4af764f2e394ff0e007c833a886cf18df9c6753..776cc6659319e51f96c1a33ffcff36807df7e92b 100644 (file)
@@ -575,18 +575,6 @@ namespace internal
     virtual FineDoFHandlerViewCell
     get_cell_view(const typename DoFHandler<dim>::cell_iterator &cell,
                   const unsigned int                             c) const = 0;
-
-    /**
-     * Return locally owned DoFs.
-     */
-    virtual const IndexSet &
-    locally_owned_dofs() const = 0;
-
-    /**
-     * Return ghost DoFs.
-     */
-    virtual const IndexSet &
-    locally_active_dofs() const = 0;
   };
 
 
@@ -599,22 +587,7 @@ namespace internal
                                const unsigned int     mg_level_fine)
       : dof_handler_fine(dof_handler_fine)
       , mg_level_fine(mg_level_fine)
-    {
-      if (this->mg_level_fine == numbers::invalid_unsigned_int)
-        {
-          is_locally_owned_dofs = dof_handler_fine.locally_owned_dofs();
-          is_locally_active_dofs =
-            DoFTools::extract_locally_active_dofs(dof_handler_fine);
-        }
-      else
-        {
-          is_locally_owned_dofs =
-            dof_handler_fine.locally_owned_mg_dofs(mg_level_fine);
-          is_locally_active_dofs =
-            DoFTools::extract_locally_active_level_dofs(dof_handler_fine,
-                                                        mg_level_fine);
-        }
-    }
+    {}
 
     virtual ~IdentityFineDoFHandlerView() = default;
 
@@ -663,26 +636,9 @@ namespace internal
         });
     }
 
-    const IndexSet &
-    locally_owned_dofs() const override
-    {
-      return is_locally_owned_dofs;
-    }
-
-    /**
-     * Return ghost DoFs.
-     */
-    const IndexSet &
-    locally_active_dofs() const override
-    {
-      return is_locally_active_dofs;
-    }
-
   private:
     const DoFHandler<dim> &dof_handler_fine;
     const unsigned int     mg_level_fine;
-    IndexSet               is_locally_owned_dofs;
-    IndexSet               is_locally_active_dofs;
   };
 
 
@@ -691,110 +647,11 @@ namespace internal
   class FirstChildPolicyFineDoFHandlerView : public FineDoFHandlerViewBase<dim>
   {
   public:
-    FirstChildPolicyFineDoFHandlerView(
-      const DoFHandler<dim> &dof_handler_fine,
-      const DoFHandler<dim> &dof_handler_coarse,
-      const unsigned int     mg_level_fine)
+    FirstChildPolicyFineDoFHandlerView(const DoFHandler<dim> &dof_handler_fine,
+                                       const unsigned int     mg_level_fine)
       : dof_handler_fine(dof_handler_fine)
       , mg_level_fine(mg_level_fine)
-    {
-      std::vector<types::global_dof_index> locally_active_non_local_indices;
-
-      if (this->mg_level_fine == numbers::invalid_unsigned_int)
-        {
-          is_locally_owned_dofs = dof_handler_fine.locally_owned_dofs();
-
-          std::vector<types::global_dof_index> dof_indices_cell;
-
-          loop_over_active_or_level_cells(
-            dof_handler_coarse,
-            numbers::invalid_unsigned_int,
-            [&](const auto &cell_coarse) {
-              // create fine cell in two steps, since the coarse cell and
-              // the fine cell are associated to different Trinagulation
-              // objects
-              const auto cell_id = cell_coarse->id();
-              const auto cell_fine_raw =
-                dof_handler_fine.get_triangulation().create_cell_iterator(
-                  cell_id);
-
-              if (cell_fine_raw->has_children() == false)
-                {
-                  // cell has no children on fine mesh
-
-                  // convert CellAccessor to DoFCellAccessor
-                  const auto cell_fine =
-                    cell_fine_raw->as_dof_handler_iterator(dof_handler_fine);
-
-                  dof_indices_cell.resize(
-                    cell_fine->get_fe().n_dofs_per_cell());
-                  cell_fine->get_dof_indices(dof_indices_cell);
-                  locally_active_non_local_indices.insert(
-                    locally_active_non_local_indices.end(),
-                    dof_indices_cell.begin(),
-                    dof_indices_cell.end());
-                }
-              else
-                {
-                  // cell has children on fine mesh: loop over all children
-                  for (const auto &child_raw : cell_fine_raw->child_iterators())
-                    {
-                      // convert CellAccessor of child to DoFCellAccessor
-                      const auto child =
-                        child_raw->as_dof_handler_iterator(dof_handler_fine);
-
-                      dof_indices_cell.resize(
-                        child->get_fe().n_dofs_per_cell());
-                      child->get_dof_indices(dof_indices_cell);
-
-                      for (const auto i : dof_indices_cell)
-                        if (is_locally_owned_dofs.is_element(i) == false)
-                          locally_active_non_local_indices.push_back(i);
-                    }
-                }
-            });
-        }
-      else
-        {
-          is_locally_owned_dofs =
-            dof_handler_fine.locally_owned_mg_dofs(mg_level_fine);
-
-          Assert(mg_level_fine > 0, ExcInternalError());
-
-          std::vector<types::global_dof_index> dof_indices_cell;
-
-          loop_over_active_or_level_cells(
-            dof_handler_fine, mg_level_fine - 1, [&](const auto &cell) {
-              if (cell->has_children())
-                {
-                  for (const auto &child : cell->child_iterators())
-                    {
-                      dof_indices_cell.resize(
-                        child->get_fe().n_dofs_per_cell());
-                      child->get_mg_dof_indices(dof_indices_cell);
-
-                      for (const auto i : dof_indices_cell)
-                        if (is_locally_owned_dofs.is_element(i) == false)
-                          locally_active_non_local_indices.push_back(i);
-                    }
-                }
-            });
-        }
-
-      is_locally_active_dofs.set_size(is_locally_owned_dofs.size());
-
-      is_locally_active_dofs.add_indices(is_locally_owned_dofs);
-
-      std::sort(locally_active_non_local_indices.begin(),
-                locally_active_non_local_indices.end());
-      locally_active_non_local_indices.erase(
-        std::unique(locally_active_non_local_indices.begin(),
-                    locally_active_non_local_indices.end()),
-        locally_active_non_local_indices.end());
-      is_locally_active_dofs.add_indices(
-        locally_active_non_local_indices.begin(),
-        locally_active_non_local_indices.end());
-    }
+    {}
 
     virtual ~FirstChildPolicyFineDoFHandlerView() = default;
 
@@ -893,26 +750,9 @@ namespace internal
         });
     }
 
-    const IndexSet &
-    locally_owned_dofs() const override
-    {
-      return is_locally_owned_dofs;
-    }
-
-    /**
-     * Return ghost DoFs.
-     */
-    const IndexSet &
-    locally_active_dofs() const override
-    {
-      return is_locally_active_dofs;
-    }
-
   private:
     const DoFHandler<dim> &dof_handler_fine;
     const unsigned int     mg_level_fine;
-    IndexSet               is_locally_owned_dofs;
-    IndexSet               is_locally_active_dofs;
   };
 
 
@@ -1008,13 +848,7 @@ namespace internal
       this->is_dst_remote        = is_dst_remote;
       this->is_src_locally_owned = is_src_locally_owned;
 
-      this->is_extended_locally_owned =
-        mg_level_fine == numbers::invalid_unsigned_int ?
-          dof_handler_fine.locally_owned_dofs() :
-          dof_handler_fine.locally_owned_mg_dofs(mg_level_fine);
-
       const auto targets_with_indexset = process.get_requesters();
-      std::vector<types::global_dof_index> ghost_indices;
 
 #ifndef DEAL_II_WITH_MPI
       Assert(targets_with_indexset.empty(), ExcInternalError());
@@ -1072,31 +906,6 @@ namespace internal
       }
 
       // process local cells
-      {
-        std::vector<types::global_dof_index> indices;
-
-        for (const auto id : is_dst_locally_owned)
-          {
-            const auto cell_id = cell_id_translator.to_cell_id(id);
-
-            typename DoFHandler<dim>::cell_iterator cell_(
-              *dof_handler_fine.get_triangulation().create_cell_iterator(
-                cell_id),
-              &dof_handler_fine);
-
-            indices.resize(cell_->get_fe().n_dofs_per_cell());
-
-            if (mg_level_fine == numbers::invalid_unsigned_int)
-              cell_->get_dof_indices(indices);
-            else
-              cell_->get_mg_dof_indices(indices);
-
-            for (const auto i : indices)
-              if (!is_extended_locally_owned.is_element(i))
-                ghost_indices.push_back(i);
-          }
-      }
-
       {
         std::map<unsigned int, std::vector<types::global_dof_index>>
           rank_to_ids;
@@ -1139,19 +948,6 @@ namespace internal
               MPI_STATUS_IGNORE);
             AssertThrowMPI(ierr_3);
 
-            for (unsigned int i = 0; i < buffer.size();)
-              {
-                const unsigned int dofs_per_cell =
-                  dof_handler_fine.get_fe(buffer[i]).n_dofs_per_cell();
-                ++i;
-                for (unsigned int j = 0; j < dofs_per_cell; ++j, ++i)
-                  {
-                    AssertIndexRange(i, buffer.size());
-                    if (!is_extended_locally_owned.is_element(buffer[i]))
-                      ghost_indices.push_back(buffer[i]);
-                  }
-              }
-
             const unsigned int rank = status.MPI_SOURCE;
 
             const auto ids = rank_to_ids[rank];
@@ -1180,16 +976,6 @@ namespace internal
           }
       }
 #endif
-
-      std::sort(ghost_indices.begin(), ghost_indices.end());
-
-      this->is_extended_ghosts =
-        IndexSet(mg_level_fine == numbers::invalid_unsigned_int ?
-                   dof_handler_fine.n_dofs() :
-                   dof_handler_fine.n_dofs(mg_level_fine));
-      this->is_extended_ghosts.add_indices(ghost_indices.begin(),
-                                           ghost_indices.end());
-      this->is_extended_ghosts.subtract_set(this->is_extended_locally_owned);
     }
 
     FineDoFHandlerViewCell
@@ -1318,18 +1104,6 @@ namespace internal
         });
     }
 
-    const IndexSet &
-    locally_owned_dofs() const override
-    {
-      return is_extended_locally_owned;
-    }
-
-    const IndexSet &
-    locally_active_dofs() const override
-    {
-      return is_extended_ghosts;
-    }
-
   private:
     const DoFHandler<dim> &dof_handler_fine;
     const DoFHandler<dim> &dof_handler_coarse;
@@ -1344,9 +1118,6 @@ namespace internal
   private:
     IndexSet is_src_locally_owned;
 
-    IndexSet is_extended_locally_owned;
-    IndexSet is_extended_ghosts;
-
     std::map<types::global_cell_index,
              std::pair<unsigned int, std::vector<types::global_dof_index>>>
       map;
@@ -1789,7 +1560,7 @@ namespace internal
                                                        mg_level_coarse))
         dof_handler_fine_view =
           std::make_unique<FirstChildPolicyFineDoFHandlerView<dim>>(
-            dof_handler_fine, dof_handler_coarse, mg_level_fine);
+            dof_handler_fine, mg_level_fine);
       else
         dof_handler_fine_view =
           std::make_unique<GlobalCoarseningFineDoFHandlerView<dim>>(
@@ -1857,28 +1628,6 @@ namespace internal
 
       const auto reference_cell = dof_handler_fine.get_fe(0).reference_cell();
 
-      // create partitioners and vectors for internal purposes
-      {
-        // ... for fine mesh
-        {
-          transfer.partitioner_fine =
-            std::make_shared<Utilities::MPI::Partitioner>(
-              dof_handler_fine_view->locally_owned_dofs(),
-              dof_handler_fine_view->locally_active_dofs(),
-              dof_handler_fine.get_communicator());
-          transfer.vec_fine.reinit(transfer.partitioner_fine);
-        }
-
-        // ... coarse mesh (needed since user vector might be const)
-        {
-          transfer.partitioner_coarse =
-            create_coarse_partitioner(dof_handler_coarse,
-                                      constraints_coarse,
-                                      mg_level_coarse);
-          transfer.vec_coarse.reinit(transfer.partitioner_coarse);
-        }
-      }
-
       // helper function: to process the fine level cells; function @p fu_non_refined is
       // performed on cells that are not refined and @fu_refined is performed on
       // children of cells that are refined
@@ -2032,28 +1781,38 @@ namespace internal
         unsigned int cell_no_0 = 0;
         unsigned int cell_no_1 = transfer.schemes[0].n_coarse_cells;
 
-        transfer.constraint_info_coarse.reinit(
-          dof_handler_coarse,
-          transfer.schemes[0].n_coarse_cells +
-            transfer.schemes[1].n_coarse_cells,
-          constraints_coarse.n_constraints() > 0 &&
-            use_fast_hanging_node_algorithm(dof_handler_coarse,
-                                            mg_level_coarse));
-
-        transfer.constraint_info_fine.reinit(
+        transfer.constraint_info_coarse
+          .template reinit<types::global_dof_index>(
+            dof_handler_coarse,
+            transfer.schemes[0].n_coarse_cells +
+              transfer.schemes[1].n_coarse_cells,
+            constraints_coarse.n_constraints() > 0 &&
+              use_fast_hanging_node_algorithm(dof_handler_coarse,
+                                              mg_level_coarse));
+        transfer.constraint_info_coarse.set_locally_owned_indices(
+          (mg_level_coarse == numbers::invalid_unsigned_int) ?
+            dof_handler_coarse.locally_owned_dofs() :
+            dof_handler_coarse.locally_owned_mg_dofs(mg_level_coarse));
+
+        transfer.constraint_info_fine.template reinit<types::global_dof_index>(
           transfer.schemes[0].n_coarse_cells +
           transfer.schemes[1].n_coarse_cells);
+        transfer.constraint_info_fine.set_locally_owned_indices(
+          (mg_level_fine == numbers::invalid_unsigned_int) ?
+            dof_handler_fine.locally_owned_dofs() :
+            dof_handler_fine.locally_owned_mg_dofs(mg_level_fine));
 
         process_cells(
           [&](const auto &cell_coarse, const auto &cell_fine) {
             // parent
             {
-              transfer.constraint_info_coarse.read_dof_indices(
-                cell_no_0,
-                mg_level_coarse,
-                cell_coarse,
-                constraints_coarse,
-                transfer.partitioner_coarse);
+              transfer.constraint_info_coarse
+                .template read_dof_indices<types::global_dof_index>(
+                  cell_no_0,
+                  mg_level_coarse,
+                  cell_coarse,
+                  constraints_coarse,
+                  {});
             }
 
             // child
@@ -2065,8 +1824,9 @@ namespace internal
                 level_dof_indices_fine_0[i] =
                   local_dof_indices[lexicographic_numbering_fine[i]];
 
-              transfer.constraint_info_fine.read_dof_indices(
-                cell_no_0, level_dof_indices_fine_0, transfer.partitioner_fine);
+              transfer.constraint_info_fine
+                .template read_dof_indices<types::global_dof_index>(
+                  cell_no_0, level_dof_indices_fine_0, {});
             }
 
             // move pointers
@@ -2078,12 +1838,13 @@ namespace internal
             // parent (only once at the beginning)
             if (c == 0)
               {
-                transfer.constraint_info_coarse.read_dof_indices(
-                  cell_no_1,
-                  mg_level_coarse,
-                  cell_coarse,
-                  constraints_coarse,
-                  transfer.partitioner_coarse);
+                transfer.constraint_info_coarse
+                  .template read_dof_indices<types::global_dof_index>(
+                    cell_no_1,
+                    mg_level_coarse,
+                    cell_coarse,
+                    constraints_coarse,
+                    {});
 
                 level_dof_indices_fine_1.assign(level_dof_indices_fine_1.size(),
                                                 numbers::invalid_dof_index);
@@ -2114,17 +1875,27 @@ namespace internal
             // move pointers (only once at the end)
             if (c + 1 == GeometryInfo<dim>::max_children_per_cell)
               {
-                transfer.constraint_info_fine.read_dof_indices(
-                  cell_no_1,
-                  level_dof_indices_fine_1,
-                  transfer.partitioner_fine);
+                transfer.constraint_info_fine
+                  .template read_dof_indices<types::global_dof_index>(
+                    cell_no_1, level_dof_indices_fine_1, {});
 
                 cell_no_1++;
               }
           });
+      }
+
+      {
+        transfer.partitioner_coarse =
+          transfer.constraint_info_coarse
+            .template finalize<types::global_dof_index>(
+              dof_handler_coarse.get_communicator());
+        transfer.vec_coarse.reinit(transfer.partitioner_coarse);
 
-        transfer.constraint_info_coarse.finalize();
-        transfer.constraint_info_fine.finalize();
+        transfer.partitioner_fine =
+          transfer.constraint_info_fine
+            .template finalize<types::global_dof_index>(
+              dof_handler_fine.get_communicator());
+        transfer.vec_fine.reinit(transfer.partitioner_fine);
       }
 
 
@@ -2391,22 +2162,6 @@ namespace internal
             dof_handler_fine.get_fe(fe_index_pair.first.second).degree;
         }
 
-      {
-        transfer.partitioner_fine =
-          std::make_shared<Utilities::MPI::Partitioner>(
-            dof_handler_fine_view->locally_owned_dofs(),
-            dof_handler_fine_view->locally_active_dofs(),
-            dof_handler_fine.get_communicator());
-        transfer.vec_fine.reinit(transfer.partitioner_fine);
-      }
-      {
-        transfer.partitioner_coarse =
-          create_coarse_partitioner(dof_handler_coarse,
-                                    constraints_coarse,
-                                    mg_level_coarse);
-        transfer.vec_coarse.reinit(transfer.partitioner_coarse);
-      }
-
       std::vector<unsigned int> n_dof_indices_fine(fe_index_pairs.size() + 1);
       std::vector<unsigned int> n_dof_indices_coarse(fe_index_pairs.size() + 1);
       std::vector<unsigned int> cell_no(fe_index_pairs.size() + 1, 0);
@@ -2513,8 +2268,16 @@ namespace internal
           constraints_coarse.n_constraints() > 0 &&
             use_fast_hanging_node_algorithm(dof_handler_coarse,
                                             mg_level_coarse));
+        transfer.constraint_info_coarse.set_locally_owned_indices(
+          (mg_level_coarse == numbers::invalid_unsigned_int) ?
+            dof_handler_coarse.locally_owned_dofs() :
+            dof_handler_coarse.locally_owned_mg_dofs(mg_level_coarse));
 
         transfer.constraint_info_fine.reinit(cell_no.back());
+        transfer.constraint_info_fine.set_locally_owned_indices(
+          (mg_level_fine == numbers::invalid_unsigned_int) ?
+            dof_handler_fine.locally_owned_dofs() :
+            dof_handler_fine.locally_owned_mg_dofs(mg_level_fine));
 
         process_cells([&](const auto &cell_coarse, const auto &cell_fine) {
           const auto fe_pair_no =
@@ -2523,12 +2286,13 @@ namespace internal
 
           // parent
           {
-            transfer.constraint_info_coarse.read_dof_indices(
-              cell_no[fe_pair_no],
-              mg_level_coarse,
-              cell_coarse,
-              constraints_coarse,
-              transfer.partitioner_coarse);
+            transfer.constraint_info_coarse
+              .template read_dof_indices<types::global_dof_index>(
+                cell_no[fe_pair_no],
+                mg_level_coarse,
+                cell_coarse,
+                constraints_coarse,
+                {});
           }
 
           // child
@@ -2541,10 +2305,11 @@ namespace internal
               local_dof_indices_fine_lex[fe_pair_no][i] = local_dof_indices_fine
                 [fe_pair_no][lexicographic_numbering_fine[fe_pair_no][i]];
 
-            transfer.constraint_info_fine.read_dof_indices(
-              cell_no[fe_pair_no],
-              local_dof_indices_fine_lex[fe_pair_no],
-              transfer.partitioner_fine);
+            transfer.constraint_info_fine
+              .template read_dof_indices<types::global_dof_index>(
+                cell_no[fe_pair_no],
+                local_dof_indices_fine_lex[fe_pair_no],
+                {});
           }
 
           // move pointers
@@ -2552,9 +2317,20 @@ namespace internal
             cell_no[fe_pair_no]++;
           }
         });
+      }
 
-        transfer.constraint_info_coarse.finalize();
-        transfer.constraint_info_fine.finalize();
+      {
+        transfer.partitioner_coarse =
+          transfer.constraint_info_coarse
+            .template finalize<types::global_dof_index>(
+              dof_handler_coarse.get_communicator());
+        transfer.vec_coarse.reinit(transfer.partitioner_coarse);
+
+        transfer.partitioner_fine =
+          transfer.constraint_info_fine
+            .template finalize<types::global_dof_index>(
+              dof_handler_fine.get_communicator());
+        transfer.vec_fine.reinit(transfer.partitioner_fine);
       }
 
       // ------------------------- prolongation matrix -------------------------

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.