]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Introduced CellWeights class for load balancing by the number of DoFs. 7264/head
authorMarc Fehling <marc.fehling@gmx.net>
Sun, 30 Sep 2018 05:03:35 +0000 (23:03 -0600)
committerMarc Fehling <marc.fehling@gmx.net>
Thu, 4 Oct 2018 22:37:25 +0000 (16:37 -0600)
doc/news/changes/minor/20181001MarcFehling [new file with mode: 0644]
include/deal.II/distributed/cell_weights.h [new file with mode: 0644]
source/distributed/CMakeLists.txt
source/distributed/cell_weights.cc [new file with mode: 0644]
source/distributed/cell_weights.inst.in [new file with mode: 0644]
tests/mpi/hp_cell_weights_01.cc [new file with mode: 0644]
tests/mpi/hp_cell_weights_01.with_p4est=true.mpirun=2.output [new file with mode: 0644]
tests/mpi/hp_cell_weights_02.cc [new file with mode: 0644]
tests/mpi/hp_cell_weights_02.with_p4est=true.mpirun=3.output [new file with mode: 0644]

diff --git a/doc/news/changes/minor/20181001MarcFehling b/doc/news/changes/minor/20181001MarcFehling
new file mode 100644 (file)
index 0000000..3431476
--- /dev/null
@@ -0,0 +1,6 @@
+New: Class parallel::distributed::CellWeights allows to set the
+weight of a cell depending on the number of degrees of freedom
+associated with it. This guarantees an improved load balancing
+in case a hp::DoFHandler is used on a p::d::Triangulation.
+<br>
+(Marc Fehling, 2018/10/01)
diff --git a/include/deal.II/distributed/cell_weights.h b/include/deal.II/distributed/cell_weights.h
new file mode 100644 (file)
index 0000000..904ffd3
--- /dev/null
@@ -0,0 +1,195 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2018 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+#ifndef dealii_distributed_cell_weights_h
+#define dealii_distributed_cell_weights_h
+
+#include <deal.II/base/config.h>
+
+#include <deal.II/distributed/tria.h>
+
+#include <deal.II/hp/dof_handler.h>
+
+
+DEAL_II_NAMESPACE_OPEN
+
+namespace parallel
+{
+  namespace distributed
+  {
+    /**
+     * Anytime a parallel::distributed::Triangulation is repartitioned, either
+     * upon request or by refinement/coarsening, cells will be distributed
+     * amongst all subdomains to achieve an equally balanced workload. If the
+     * workload per cell varies, which is in general the case for hp::DoFHandler
+     * objects, we can take those into account by introducing individual weights
+     * for different cells.
+     *
+     * This class allows to consult the FiniteElement, that it associated with
+     * each cell by the hp::DoFHandler, to determine the weight of the cell for
+     * load balancing. One can choose from predefined weighting algorithms
+     * provided by this class or provide a custom one. The chosen weighting
+     * function will be connected to the corresponding signal of the linked
+     * parallel::distributed::Triangulation via callback.
+     *
+     * An object of this class needs to exist for every DoFHandler associated
+     * with the Triangulation we work on to achieve satisfying work balancing
+     * results.
+     *
+     * A small code snippet follows explaining how to achieve each cell
+     * being weighted by its current number of degrees of freedom.
+     * @code
+     * parallel::distributed::CellWeights<dim, spacedim>
+     *   cell_weights(hp_dof_handler);
+     * cell_weights.register_ndofs_weighting();
+     * @endcode
+     * The weighting function can be changed anytime. Even more ambitious
+     * approaches are possible by submitting customized functions, e.g.
+     * @code
+     * cell_weights.register_custom_weighting(
+     *   [](const FiniteElement<dim, spacedim> &active_fe,
+     *      const typename hp::DoFHandler<dim, spacedim>::cell_iterator &)
+     *   -> unsigned int {
+     *   return 1000 * std::pow(active_fe.dofs_per_cell, 1.6);
+     * });
+     * @endcode
+     * The returned value has to be an unsigned integer and is thus limited by
+     * some large number. It is interpreted as the additional computational load
+     * of each cell. See Triangulation::Signals::cell_weight for a discussion on
+     * this topic.
+     *
+     * @note Be aware that this class connects the weight function to the
+     * Triangulation during its construction. If the Triangulation
+     * associated with the DoFHandler changes during the lifetime of the
+     * latter, an assertion will be triggered in the weight_callback() function.
+     * It is recommended to create a separate object in this case and to destroy
+     * the previous one.
+     *
+     * @ingroup distributed
+     * @author Marc Fehling, 2018
+     */
+    template <int dim, int spacedim = dim>
+    class CellWeights
+    {
+    public:
+      /**
+       * Constructor.
+       *
+       * @param[in] dof_handler The hp::DoFHandler which will be used to
+       *    determine each cell's finite element.
+       */
+      CellWeights(const dealii::hp::DoFHandler<dim, spacedim> &dof_handler);
+
+      /**
+       * Destructor.
+       */
+      ~CellWeights();
+
+      /**
+       * Choose a constant weight @p factor on each cell.
+       */
+      void
+      register_constant_weighting(const unsigned int factor = 1000);
+
+      /**
+       * Choose a weight for each cell that is proportional to its number of
+       * degrees of freedom with a factor @p factor.
+       */
+      void
+      register_ndofs_weighting(const unsigned int factor = 1000);
+
+      /**
+       * Choose a weight for each cell that is proportional to its number of
+       * degrees of freedom <i>squared</i> with a factor @p factor.
+       */
+      void
+      register_ndofs_squared_weighting(const unsigned int factor = 1000);
+
+      /**
+       * Register a custom weight for each cell by providing a function as a
+       * parameter.
+       *
+       * @param[in] custom_function The custom weighting function returning
+       *    the weight of each cell as an unsigned integer. It is required
+       *    to have two arguments, namely the FiniteElement that will be
+       *    active on the particular cell, and the cell itself of type
+       *    hp::DoFHandler::cell_iterator. We require both to make sure to
+       *    get the right active FiniteElement on each cell in case that we
+       *    coarsen the Triangulation.
+       */
+      void
+      register_custom_weighting(
+        const std::function<unsigned int(
+          const FiniteElement<dim, spacedim> &,
+          const typename hp::DoFHandler<dim, spacedim>::cell_iterator &)>
+          custom_function);
+
+    private:
+      /**
+       * Pointer to the degree of freedom handler.
+       */
+      SmartPointer<const dealii::hp::DoFHandler<dim, spacedim>, CellWeights>
+        dof_handler;
+
+      /**
+       * Pointer to the Triangulation associated with the degree of freedom
+       * handler.
+       *
+       * We store both to make sure to always work on the correct combination of
+       * both.
+       */
+      SmartPointer<const parallel::distributed::Triangulation<dim, spacedim>,
+                   CellWeights>
+        triangulation;
+
+      /**
+       * Function that will determine each cell's weight.
+       *
+       * Can be set using the register_constant_weighting(),
+       * register_ndofs_weighting(), register_ndofs_squared_weighting(), and
+       * register_custom_weighting() member functions.
+       *
+       * The function requires the active FiniteElement object on each cell
+       * as an argument, as well as the cell itself of type
+       * hp::DoFHandler::cell_iterator.
+       */
+      std::function<unsigned int(
+        const FiniteElement<dim, spacedim> &,
+        const typename hp::DoFHandler<dim, spacedim>::cell_iterator &)>
+        weighting_function;
+
+      /**
+       * A connection to the Triangulation of the DoFHandler.
+       */
+      boost::signals2::connection tria_listener;
+
+      /**
+       * A callback function that will be attached to the cell_weight signal of
+       * the Triangulation, that is a member of the DoFHandler. Ultimately
+       * returns the weight for each cell, determined by the weighting_function
+       * member.
+       */
+      unsigned int
+      weight_callback(
+        const typename Triangulation<dim, spacedim>::cell_iterator &cell,
+        const typename Triangulation<dim, spacedim>::CellStatus     status);
+    };
+  } // namespace distributed
+} // namespace parallel
+
+
+DEAL_II_NAMESPACE_CLOSE
+
+#endif
index fa97e30c5e94cc94efe4c5329578717adb0ee49b..c2707d59871747fc5b8a2dfc03d459bd25eac3f4 100644 (file)
@@ -17,6 +17,7 @@ INCLUDE_DIRECTORIES(BEFORE ${CMAKE_CURRENT_BINARY_DIR})
 
 SET(_unity_include_src
   grid_refinement.cc
+  cell_weights.cc
   active_fe_indices_transfer.cc
   solution_transfer.cc
   tria.cc
@@ -39,6 +40,7 @@ SETUP_SOURCE_LIST("${_unity_include_src}"
 
 SET(_inst
   grid_refinement.inst.in
+  cell_weights.inst.in
   active_fe_indices_transfer.inst.in
   solution_transfer.inst.in
   tria.inst.in
diff --git a/source/distributed/cell_weights.cc b/source/distributed/cell_weights.cc
new file mode 100644 (file)
index 0000000..750375b
--- /dev/null
@@ -0,0 +1,186 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2018 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+
+#include <deal.II/base/config.h>
+
+#include <deal.II/distributed/cell_weights.h>
+
+#include <deal.II/dofs/dof_accessor.h>
+
+
+DEAL_II_NAMESPACE_OPEN
+
+namespace parallel
+{
+  namespace distributed
+  {
+    template <int dim, int spacedim>
+    CellWeights<dim, spacedim>::CellWeights(
+      const hp::DoFHandler<dim, spacedim> &dof_handler)
+      : dof_handler(&dof_handler, typeid(*this).name())
+    {
+      triangulation =
+        (dynamic_cast<parallel::distributed::Triangulation<dim, spacedim> *>(
+          const_cast<dealii::Triangulation<dim, spacedim> *>(
+            &(this->dof_handler->get_triangulation()))));
+
+      if (triangulation != nullptr)
+        {
+          // Choose default weighting.
+          register_constant_weighting();
+
+          // Add callback function to the cell_weight signal and store its
+          // connection.
+          tria_listener = triangulation->signals.cell_weight.connect(
+            std::bind(&CellWeights<dim, spacedim>::weight_callback,
+                      std::ref(*this),
+                      std::placeholders::_1,
+                      std::placeholders::_2));
+        }
+      else
+        Assert(
+          triangulation != nullptr,
+          ExcMessage(
+            "parallel::distributed::CellWeights requires a parallel::distributed::Triangulation object."));
+    }
+
+
+    template <int dim, int spacedim>
+    CellWeights<dim, spacedim>::~CellWeights()
+    {
+      tria_listener.disconnect();
+    }
+
+
+
+    template <int dim, int spacedim>
+    void
+    CellWeights<dim, spacedim>::register_constant_weighting(
+      const unsigned int factor)
+    {
+      weighting_function =
+        [factor](const FiniteElement<dim, spacedim> &,
+                 const typename hp::DoFHandler<dim, spacedim>::cell_iterator &)
+        -> unsigned int { return factor; };
+    }
+
+
+    template <int dim, int spacedim>
+    void
+    CellWeights<dim, spacedim>::register_ndofs_weighting(
+      const unsigned int factor)
+    {
+      weighting_function =
+        [factor](const FiniteElement<dim, spacedim> &active_fe,
+                 const typename hp::DoFHandler<dim, spacedim>::cell_iterator &)
+        -> unsigned int { return factor * active_fe.dofs_per_cell; };
+    }
+
+
+    template <int dim, int spacedim>
+    void
+    CellWeights<dim, spacedim>::register_ndofs_squared_weighting(
+      const unsigned int factor)
+    {
+      weighting_function =
+        [factor](const FiniteElement<dim, spacedim> &active_fe,
+                 const typename hp::DoFHandler<dim, spacedim>::cell_iterator &)
+        -> unsigned int {
+        return factor * active_fe.dofs_per_cell * active_fe.dofs_per_cell;
+      };
+    }
+
+
+    template <int dim, int spacedim>
+    void
+    CellWeights<dim, spacedim>::register_custom_weighting(
+      const std::function<unsigned int(
+        const FiniteElement<dim, spacedim> &,
+        const typename hp::DoFHandler<dim, spacedim>::cell_iterator &)>
+        custom_function)
+    {
+      weighting_function = custom_function;
+    }
+
+
+
+    template <int dim, int spacedim>
+    unsigned int
+    CellWeights<dim, spacedim>::weight_callback(
+      const typename Triangulation<dim, spacedim>::cell_iterator &cell_,
+      const typename Triangulation<dim, spacedim>::CellStatus     status)
+    {
+      // Check if we are still working with the correct combination of
+      // Triangulation and DoFHandler.
+      Assert(&(*triangulation) == &(dof_handler->get_triangulation()),
+             ExcMessage(
+               "Triangulation associated with the DoFHandler has changed!"));
+
+      // Convert cell type from Triangulation to DoFHandler to be able
+      // to access the information about the degrees of freedom.
+      const typename hp::DoFHandler<dim, spacedim>::cell_iterator cell(
+        *cell_, dof_handler);
+
+      // Determine which FiniteElement object will be present on this cell after
+      // refinement and will thus specify the number of degrees of freedom.
+      unsigned int fe_index = numbers::invalid_unsigned_int;
+      switch (status)
+        {
+          case parallel::distributed::Triangulation<dim,
+                                                    spacedim>::CELL_PERSIST:
+          case parallel::distributed::Triangulation<dim, spacedim>::CELL_REFINE:
+          case parallel::distributed::Triangulation<dim,
+                                                    spacedim>::CELL_INVALID:
+            fe_index = cell->active_fe_index();
+            break;
+
+          case parallel::distributed::Triangulation<dim,
+                                                    spacedim>::CELL_COARSEN:
+            {
+              std::set<unsigned int> fe_indices_children;
+              for (unsigned int child_index = 0;
+                   child_index < GeometryInfo<dim>::max_children_per_cell;
+                   ++child_index)
+                fe_indices_children.insert(
+                  cell->child(child_index)->active_fe_index());
+
+              fe_index = dof_handler->get_fe().find_least_face_dominating_fe(
+                fe_indices_children);
+
+              Assert(fe_index != numbers::invalid_unsigned_int,
+                     ExcMessage(
+                       "No FiniteElement has been found in your FECollection "
+                       "that dominates all children of a cell you are trying "
+                       "to coarsen!"));
+            }
+            break;
+
+          default:
+            Assert(false, ExcInternalError());
+            break;
+        }
+
+      // Return the cell weight determined by the function of choice.
+      return weighting_function(dof_handler->get_fe(fe_index), cell);
+    }
+  } // namespace distributed
+} // namespace parallel
+
+
+// explicit instantiations
+#include "cell_weights.inst"
+
+DEAL_II_NAMESPACE_CLOSE
diff --git a/source/distributed/cell_weights.inst.in b/source/distributed/cell_weights.inst.in
new file mode 100644 (file)
index 0000000..b630293
--- /dev/null
@@ -0,0 +1,29 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2018 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+
+
+for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS)
+  {
+    namespace parallel
+    \{
+      namespace distributed
+      \{
+#if deal_II_dimension <= deal_II_space_dimension
+        template class CellWeights<deal_II_dimension, deal_II_space_dimension>;
+#endif
+      \}
+    \}
+  }
diff --git a/tests/mpi/hp_cell_weights_01.cc b/tests/mpi/hp_cell_weights_01.cc
new file mode 100644 (file)
index 0000000..e99980c
--- /dev/null
@@ -0,0 +1,119 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2018 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+
+
+// Cell Weights Test
+// -----------------
+// Create a 4x4(x4) grid, on which all cells are associated with a Q1
+// element besides the very first one, which has a Q5 element.
+// We choose a cell weighting algorithm based on the number of degrees
+// of freedom and check if load is balanced as expected after
+// repartitioning the triangulation. The expected accumulated weight on
+// each processor should correlate to the sum of all degrees of
+// freedom on all cells of the corresponding subdomain.
+// We employ a large proportionality factor on our weighting function
+// to neglect the standard weight of '1000' per cell.
+
+
+#include <deal.II/distributed/active_fe_indices_transfer.h>
+#include <deal.II/distributed/cell_weights.h>
+#include <deal.II/distributed/tria.h>
+
+#include <deal.II/fe/fe_q.h>
+
+#include <deal.II/grid/grid_generator.h>
+
+#include <deal.II/hp/dof_handler.h>
+
+#include "../tests.h"
+
+
+template <int dim>
+void
+test()
+{
+  parallel::distributed::Triangulation<dim> tria(MPI_COMM_WORLD);
+
+  GridGenerator::hyper_cube(tria);
+  tria.refine_global(2);
+
+  // Apply ndof cell weights.
+  hp::FECollection<dim> fe_collection;
+  fe_collection.push_back(FE_Q<dim>(1));
+  fe_collection.push_back(FE_Q<dim>(5));
+
+  hp::DoFHandler<dim> dh(tria);
+  // default: active_fe_index = 0
+  for (auto &cell : dh.active_cell_iterators())
+    if (cell->is_locally_owned())
+      if (cell->id().to_string() == "0_2:00")
+        cell->set_active_fe_index(1);
+  dh.distribute_dofs(fe_collection);
+
+
+  parallel::distributed::ActiveFEIndicesTransfer<dim> feidx_transfer(dh);
+  feidx_transfer.prepare_for_transfer();
+
+  parallel::distributed::CellWeights<dim> cell_weights(dh);
+  cell_weights.register_ndofs_weighting(100000);
+
+
+  deallog << "Number of cells before repartitioning: "
+          << tria.n_locally_owned_active_cells() << std::endl;
+  {
+    unsigned int dof_counter = 0;
+    for (auto &cell : dh.active_cell_iterators())
+      if (cell->is_locally_owned())
+        dof_counter += cell->get_fe().dofs_per_cell;
+    deallog << "  Cumulative dofs per cell: " << dof_counter << std::endl;
+  }
+
+
+  tria.repartition();
+  feidx_transfer.unpack();
+  dh.distribute_dofs(fe_collection);
+
+
+  deallog << "Number of cells after repartitioning: "
+          << tria.n_locally_owned_active_cells() << std::endl;
+  {
+    unsigned int dof_counter = 0;
+    for (auto &cell : dh.active_cell_iterators())
+      if (cell->is_locally_owned())
+        dof_counter += cell->get_fe().dofs_per_cell;
+    deallog << "  Cumulative dofs per cell: " << dof_counter << std::endl;
+  }
+
+  // make sure no processor is hanging
+  MPI_Barrier(MPI_COMM_WORLD);
+
+  deallog << "OK" << std::endl;
+}
+
+
+int
+main(int argc, char *argv[])
+{
+  Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1);
+  MPILogInitAll                    log;
+
+  deallog.push("2d");
+  test<2>();
+  deallog.pop();
+  deallog.push("3d");
+  test<3>();
+  deallog.pop();
+}
diff --git a/tests/mpi/hp_cell_weights_01.with_p4est=true.mpirun=2.output b/tests/mpi/hp_cell_weights_01.with_p4est=true.mpirun=2.output
new file mode 100644 (file)
index 0000000..e04bfd1
--- /dev/null
@@ -0,0 +1,23 @@
+
+DEAL:0:2d::Number of cells before repartitioning: 8
+DEAL:0:2d::  Cumulative dofs per cell: 64
+DEAL:0:2d::Number of cells after repartitioning: 4
+DEAL:0:2d::  Cumulative dofs per cell: 48
+DEAL:0:2d::OK
+DEAL:0:3d::Number of cells before repartitioning: 32
+DEAL:0:3d::  Cumulative dofs per cell: 464
+DEAL:0:3d::Number of cells after repartitioning: 24
+DEAL:0:3d::  Cumulative dofs per cell: 400
+DEAL:0:3d::OK
+
+DEAL:1:2d::Number of cells before repartitioning: 8
+DEAL:1:2d::  Cumulative dofs per cell: 32
+DEAL:1:2d::Number of cells after repartitioning: 12
+DEAL:1:2d::  Cumulative dofs per cell: 48
+DEAL:1:2d::OK
+DEAL:1:3d::Number of cells before repartitioning: 32
+DEAL:1:3d::  Cumulative dofs per cell: 256
+DEAL:1:3d::Number of cells after repartitioning: 40
+DEAL:1:3d::  Cumulative dofs per cell: 320
+DEAL:1:3d::OK
+
diff --git a/tests/mpi/hp_cell_weights_02.cc b/tests/mpi/hp_cell_weights_02.cc
new file mode 100644 (file)
index 0000000..81e00a4
--- /dev/null
@@ -0,0 +1,125 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2018 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+
+
+// Cell Weights Test
+// -----------------
+// Create a larger 8x8(x8) grid, on which all cells are associated with
+// a Q1 element, besides the very first one which has a Q7 element.
+// We choose a cell weighting algorithm based on the number of degrees
+// of freedom and check if load is balanced as expected after
+// repartitioning the triangulation. The expected accumulated weight on
+// each processor should correlate to the sum of all degrees of
+// freedom on all cells of the corresponding subdomain.
+// We employ a large proportionality factor on our weighting function
+// to neglect the standard weight of '1000' per cell.
+//
+// This test runs on a larger domain with a Lagrangian element of higher
+// order, compared to the previous one. If we would have used a Q7
+// element on the smaller grid, load balancing would fail in such a way
+// that the second (last) processor owns the whole domain -- p4est wants
+// to 'cut' its tree on a parent branch that does not exist in this case.
+
+
+#include <deal.II/distributed/active_fe_indices_transfer.h>
+#include <deal.II/distributed/cell_weights.h>
+#include <deal.II/distributed/tria.h>
+
+#include <deal.II/fe/fe_q.h>
+
+#include <deal.II/grid/grid_generator.h>
+
+#include <deal.II/hp/dof_handler.h>
+
+#include "../tests.h"
+
+
+template <int dim>
+void
+test()
+{
+  parallel::distributed::Triangulation<dim> tria(MPI_COMM_WORLD);
+
+  GridGenerator::hyper_cube(tria);
+  tria.refine_global(3);
+
+  // Apply ndof cell weights.
+  hp::FECollection<dim> fe_collection;
+  fe_collection.push_back(FE_Q<dim>(1));
+  fe_collection.push_back(FE_Q<dim>(7));
+
+  hp::DoFHandler<dim> dh(tria);
+  // default: active_fe_index = 0
+  for (auto &cell : dh.active_cell_iterators())
+    if (cell->is_locally_owned())
+      if (cell->id().to_string() == "0_3:000")
+        cell->set_active_fe_index(1);
+  dh.distribute_dofs(fe_collection);
+
+
+  parallel::distributed::ActiveFEIndicesTransfer<dim> feidx_transfer(dh);
+  feidx_transfer.prepare_for_transfer();
+
+  parallel::distributed::CellWeights<dim> cell_weights(dh);
+  cell_weights.register_ndofs_weighting(100000);
+
+
+  deallog << "Number of cells before repartitioning: "
+          << tria.n_locally_owned_active_cells() << std::endl;
+  {
+    unsigned int dof_counter = 0;
+    for (auto &cell : dh.active_cell_iterators())
+      if (cell->is_locally_owned())
+        dof_counter += cell->get_fe().dofs_per_cell;
+    deallog << "  Cumulative dofs per cell: " << dof_counter << std::endl;
+  }
+
+
+  tria.repartition();
+  feidx_transfer.unpack();
+  dh.distribute_dofs(fe_collection);
+
+
+  deallog << "Number of cells after repartitioning: "
+          << tria.n_locally_owned_active_cells() << std::endl;
+  {
+    unsigned int dof_counter = 0;
+    for (auto &cell : dh.active_cell_iterators())
+      if (cell->is_locally_owned())
+        dof_counter += cell->get_fe().dofs_per_cell;
+    deallog << "  Cumulative dofs per cell: " << dof_counter << std::endl;
+  }
+
+  // make sure no processor is hanging
+  MPI_Barrier(MPI_COMM_WORLD);
+
+  deallog << "OK" << std::endl;
+}
+
+
+int
+main(int argc, char *argv[])
+{
+  Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1);
+  MPILogInitAll                    log;
+
+  deallog.push("2d");
+  test<2>();
+  deallog.pop();
+  deallog.push("3d");
+  test<3>();
+  deallog.pop();
+}
diff --git a/tests/mpi/hp_cell_weights_02.with_p4est=true.mpirun=3.output b/tests/mpi/hp_cell_weights_02.with_p4est=true.mpirun=3.output
new file mode 100644 (file)
index 0000000..bf02d55
--- /dev/null
@@ -0,0 +1,35 @@
+
+DEAL:0:2d::Number of cells before repartitioning: 20
+DEAL:0:2d::  Cumulative dofs per cell: 140
+DEAL:0:2d::Number of cells after repartitioning: 12
+DEAL:0:2d::  Cumulative dofs per cell: 108
+DEAL:0:2d::OK
+DEAL:0:3d::Number of cells before repartitioning: 168
+DEAL:0:3d::  Cumulative dofs per cell: 1848
+DEAL:0:3d::Number of cells after repartitioning: 128
+DEAL:0:3d::  Cumulative dofs per cell: 1528
+DEAL:0:3d::OK
+
+DEAL:1:2d::Number of cells before repartitioning: 24
+DEAL:1:2d::  Cumulative dofs per cell: 96
+DEAL:1:2d::Number of cells after repartitioning: 28
+DEAL:1:2d::  Cumulative dofs per cell: 112
+DEAL:1:2d::OK
+DEAL:1:3d::Number of cells before repartitioning: 176
+DEAL:1:3d::  Cumulative dofs per cell: 1408
+DEAL:1:3d::Number of cells after repartitioning: 192
+DEAL:1:3d::  Cumulative dofs per cell: 1536
+DEAL:1:3d::OK
+
+
+DEAL:2:2d::Number of cells before repartitioning: 20
+DEAL:2:2d::  Cumulative dofs per cell: 80
+DEAL:2:2d::Number of cells after repartitioning: 24
+DEAL:2:2d::  Cumulative dofs per cell: 96
+DEAL:2:2d::OK
+DEAL:2:3d::Number of cells before repartitioning: 168
+DEAL:2:3d::  Cumulative dofs per cell: 1344
+DEAL:2:3d::Number of cells after repartitioning: 192
+DEAL:2:3d::  Cumulative dofs per cell: 1536
+DEAL:2:3d::OK
+

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.