]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Refactored parallel::CellWeights. 8687/head
authorMarc Fehling <marc.fehling@gmx.net>
Tue, 3 Sep 2019 10:25:02 +0000 (12:25 +0200)
committerMarc Fehling <marc.fehling@gmx.net>
Wed, 30 Oct 2019 09:47:54 +0000 (10:47 +0100)
15 files changed:
doc/news/changes/incompatibilities/20190903Fehling [new file with mode: 0644]
include/deal.II/distributed/cell_weights.h
source/distributed/cell_weights.cc
source/hp/dof_handler.cc
tests/mpi/hp_active_fe_indices_transfer_03.cc
tests/mpi/hp_cell_weights_01.cc
tests/mpi/hp_cell_weights_01.with_p4est=true.mpirun=2.output
tests/mpi/hp_cell_weights_02.cc
tests/mpi/hp_cell_weights_02.with_p4est=true.mpirun=3.output
tests/mpi/hp_cell_weights_03.cc
tests/mpi/hp_cell_weights_03.with_mpi=true.with_metis=true.mpirun=2.output
tests/mpi/hp_cell_weights_03.with_mpi=true.with_metis=true.mpirun=2.output.1
tests/mpi/hp_cell_weights_04.cc
tests/mpi/hp_cell_weights_04.with_mpi=true.with_metis=true.mpirun=2.output
tests/mpi/hp_cell_weights_04.with_mpi=true.with_metis=true.mpirun=2.output.1

diff --git a/doc/news/changes/incompatibilities/20190903Fehling b/doc/news/changes/incompatibilities/20190903Fehling
new file mode 100644 (file)
index 0000000..ef66391
--- /dev/null
@@ -0,0 +1,6 @@
+Changed: The interface to the parallel::CellWeights class was generalized
+by introducing static member functions. This way, users are capable of
+connecting and disconnecting callback functions manually if desired without
+instantiating an object of this class. The legacy interface has been deprecated.
+<br>
+(Marc Fehling, 2019/09/03)
index 05eaae20d7aef7f005c9430bd9f7833b4396fc46..4ccd92b51c43ebaaa7dc586f3f39715e7f6a454b 100644 (file)
@@ -1,6 +1,6 @@
 // ---------------------------------------------------------------------
 //
-// Copyright (C) 2018 by the deal.II authors
+// Copyright (C) 2018 - 2019 by the deal.II authors
 //
 // This file is part of the deal.II library.
 //
@@ -38,80 +38,180 @@ namespace parallel
    * This class allows computing these weights for load balancing by
    * consulting the FiniteElement that is associated with each cell of
    * a hp::DoFHandler. 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::TriangulationBase via callback.
+   * algorithms provided by this class or provide a custom one.
    *
-   * 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.
+   * This class offers two different ways of connecting the chosen weighting
+   * function to the corresponding signal of the linked
+   * parallel::TriangulationBase. The recommended way involves creating an
+   * object of this class which will automatically take care of registering the
+   * weighting function upon creation and de-registering it once destroyed. An
+   * object of this class needs to exist for every DoFHandler associated with
+   * the Triangulation we work on to achieve satisfying work balancing results.
+   * The connected weighting function may be changed anytime using the
+   * CellWeights::reinit() function. The following code snippet demonstrates how
+   * to achieve each cell being weighted by its current number of degrees of
+   * freedom. We chose a factor of `1000` that corresponds to the initial weight
+   * each cell is assigned to upon creation.
    * @code
-   * parallel::CellWeights<dim, spacedim> cell_weights(hp_dof_handler);
-   * cell_weights.register_ndofs_weighting();
+   * parallel::CellWeights<dim, spacedim> cell_weights(
+   *   hp_dof_handler,
+   *   parallel::CellWeights<dim, spacedim>::ndofs_weighting({1000, 1}));
    * @endcode
-   * The weighting function can be changed anytime. Even more ambitious
-   * approaches are possible by submitting customized functions, e.g.
+   *
+   * On the other hand, you are also able to take care of handling the signal
+   * connection manually by using the static member function of this class. In
+   * this case, an analogous code example looks as follows.
    * @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);
-   * });
+   * boost::signals2::connection connection =
+   *   hp_dof_handler.get_triangulation().signals.cell_weight.connect(
+   *     parallel::CellWeights<dim, spacedim>::make_weighting_callback(
+   *       hp_dof_handler,
+   *       parallel::CellWeights<dim, spacedim>::ndofs_weighting(
+   *         {1000, 1}));
    * @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 See Triangulation::Signals::cell_weight for more information on
+   * weighting and load balancing.
    *
    * @note Be aware that this class connects the weight function to the
-   * Triangulation during its construction. If the Triangulation
+   * Triangulation during this class's constructor. 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.
+   * latter via hp::DoFHandler::initialize(), an assertion will be triggered in
+   * the weight_callback() function. Use CellWeights::reinit() to deregister the
+   * weighting function on the old Triangulation and connect it to the new one.
    *
    * @ingroup distributed
-   * @author Marc Fehling, 2018
+   * @author Marc Fehling, 2018, 2019
    */
   template <int dim, int spacedim = dim>
   class CellWeights
   {
   public:
+    /**
+     * An alias that defines the characteristics of a function that can be used
+     * for weighting cells during load balancing.
+     *
+     * Such weighting functions take as arguments an iterator to a cell and the
+     * future finite element that will be assigned to it after repartitioning.
+     * They return an unsigned integer which is interpreted as the cell's
+     * weight or, in other words, the additional computational load associated
+     * with it.
+     */
+    using WeightingFunction = std::function<unsigned int(
+      const typename hp::DoFHandler<dim, spacedim>::cell_iterator &,
+      const FiniteElement<dim, spacedim> &)>;
+
     /**
      * Constructor.
      *
      * @param[in] dof_handler The hp::DoFHandler which will be used to
      *    determine each cell's finite element.
+     * @param[in] weighting_function The function that determines each
+     *    cell's weight during load balancing.
      */
-    CellWeights(const dealii::hp::DoFHandler<dim, spacedim> &dof_handler);
+    CellWeights(const dealii::hp::DoFHandler<dim, spacedim> &dof_handler,
+                const WeightingFunction &weighting_function);
 
     /**
      * Destructor.
+     *
+     * Disconnects the function previously connected to the weighting signal.
      */
     ~CellWeights();
 
     /**
-     * Choose a constant weight @p factor on each cell.
+     * Connect a different @p weighting_function to the Triangulation
+     * associated with the @p dof_handler.
+     *
+     * Disconnects the function previously connected to the weighting signal.
      */
     void
+    reinit(const hp::DoFHandler<dim, spacedim> &dof_handler,
+           const WeightingFunction &            weighting_function);
+
+    /**
+     * Converts a @p weighting_function to a different type that qualifies as
+     * a callback function, which can be connected to a weighting signal of a
+     * Triangulation.
+     *
+     * This function does <b>not</b> connect the converted function to the
+     * Triangulation associated with the @p dof_handler.
+     */
+    static std::function<unsigned int(
+      const typename Triangulation<dim, spacedim>::cell_iterator &cell,
+      const typename Triangulation<dim, spacedim>::CellStatus     status)>
+    make_weighting_callback(const hp::DoFHandler<dim, spacedim> &dof_handler,
+                            const WeightingFunction &weighting_function);
+
+    /**
+     * @name Selection of weighting functions
+     * @{
+     */
+
+    /**
+     * Choose a constant weight @p factor on each cell.
+     */
+    static WeightingFunction
+    constant_weighting(const unsigned int factor = 1000);
+
+    /**
+     * The pair of floating point numbers $(a,b)$ provided via
+     * @p coefficients determines the weight $w_K$ of each cell $K$ with
+     * $n_K$ degrees of freedom in the following way: \f[ w_K =
+     * a \, n_K^b \f]
+     *
+     * The right hand side will be rounded to the nearest integer since cell
+     * weights are required to be integers.
+     */
+    static WeightingFunction
+    ndofs_weighting(const std::pair<float, float> &coefficients);
+
+    /**
+     * The container @p coefficients provides pairs of floating point numbers
+     * $(a_i, b_i)$ that determine the weight $w_K$ of each cell
+     * $K$ with $n_K$ degrees of freedom in the following way: \f[ w_K =
+     * \sum_i a_i \, n_K^{b_i} \f]
+     *
+     * The right hand side will be rounded to the nearest integer since cell
+     * weights are required to be integers.
+     */
+    static WeightingFunction
+    ndofs_weighting(const std::vector<std::pair<float, float>> &coefficients);
+
+    /**
+     * @}
+     */
+
+    /**
+     * @name Deprecated functions
+     * @{
+     */
+
+    /**
+     * Constructor.
+     *
+     * @param[in] dof_handler The hp::DoFHandler which will be used to
+     *    determine each cell's finite element.
+     */
+    DEAL_II_DEPRECATED
+    CellWeights(const dealii::hp::DoFHandler<dim, spacedim> &dof_handler);
+
+    /**
+     * @copydoc CellWeights::constant_weighting()
+     */
+    DEAL_II_DEPRECATED 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.
+     * @copydoc CellWeights::ndofs_weighting()
      */
-    void
+    DEAL_II_DEPRECATED 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.
+     * @copydoc CellWeights::ndofs_squared_weighting()
      */
-    void
+    DEAL_II_DEPRECATED void
     register_ndofs_squared_weighting(const unsigned int factor = 1000);
 
     /**
@@ -126,17 +226,27 @@ namespace parallel
      *    get the right active FiniteElement on each cell in case that we
      *    coarsen the Triangulation.
      */
-    void
+    DEAL_II_DEPRECATED void
     register_custom_weighting(
       const std::function<unsigned int(
         const FiniteElement<dim, spacedim> &,
         const typename hp::DoFHandler<dim, spacedim>::cell_iterator &)>
         custom_function);
 
+    /**
+     * @}
+     */
+
   private:
+    /**
+     * @name Deprecated members
+     * @{
+     */
+
     /**
      * Pointer to the degree of freedom handler.
      */
+    DEAL_II_DEPRECATED
     SmartPointer<const dealii::hp::DoFHandler<dim, spacedim>, CellWeights>
       dof_handler;
 
@@ -147,40 +257,33 @@ namespace parallel
      * We store both to make sure to always work on the correct combination of
      * both.
      */
+    DEAL_II_DEPRECATED
     SmartPointer<const parallel::TriangulationBase<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.
+     * A connection to the corresponding cell_weight signal of the Triangulation
+     * which is attached to the DoFHandler.
      */
-    boost::signals2::connection tria_listener;
+    boost::signals2::connection connection;
 
     /**
-     * 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.
+     * A callback function that will be connected to the cell_weight signal of
+     * the @p triangulation, to which the @p dof_handler is attached. Ultimately
+     * returns the weight for each cell, determined by the @p weighting_function
+     * provided as a parameter.
      */
-    unsigned int
-    weight_callback(
+    static unsigned int
+    weighting_callback(
       const typename Triangulation<dim, spacedim>::cell_iterator &cell,
-      const typename Triangulation<dim, spacedim>::CellStatus     status);
+      const typename Triangulation<dim, spacedim>::CellStatus     status,
+      const hp::DoFHandler<dim, spacedim> &                       dof_handler,
+      const parallel::TriangulationBase<dim, spacedim> &          triangulation,
+      const WeightingFunction &weighting_function);
   };
 } // namespace parallel
 
index 0e92fa111c99d03a226e3faf3d04341fffc234b4..2eb34fcaec9f1ab82fd87903903ad77a4c7bcdff 100644 (file)
 
 DEAL_II_NAMESPACE_OPEN
 
+
 namespace parallel
 {
   template <int dim, int spacedim>
   CellWeights<dim, spacedim>::CellWeights(
-    const hp::DoFHandler<dim, spacedim> &dof_handler)
-    : dof_handler(&dof_handler, typeid(*this).name())
+    const hp::DoFHandler<dim, spacedim> &dof_handler,
+    const typename CellWeights<dim, spacedim>::WeightingFunction
+      &weighting_function)
   {
-    triangulation = (dynamic_cast<parallel::TriangulationBase<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(
-          [this](
-            const typename Triangulation<dim, spacedim>::cell_iterator &cell_,
-            const typename Triangulation<dim, spacedim>::CellStatus status) {
-            return this->weight_callback(cell_, status);
-          });
-      }
-    else
-      Assert(
-        triangulation != nullptr,
-        ExcMessage(
-          "parallel::CellWeights requires a parallel::TriangulationBase object."));
+    reinit(dof_handler, weighting_function);
   }
 
 
+
   template <int dim, int spacedim>
   CellWeights<dim, spacedim>::~CellWeights()
   {
-    tria_listener.disconnect();
+    connection.disconnect();
   }
 
 
 
   template <int dim, int spacedim>
   void
-  CellWeights<dim, spacedim>::register_constant_weighting(
-    const unsigned int factor)
+  CellWeights<dim, spacedim>::reinit(
+    const hp::DoFHandler<dim, spacedim> &dof_handler,
+    const typename CellWeights<dim, spacedim>::WeightingFunction
+      &weighting_function)
   {
-    weighting_function =
-      [factor](const FiniteElement<dim, spacedim> &,
-               const typename hp::DoFHandler<dim, spacedim>::cell_iterator &)
-      -> unsigned int { return factor; };
+    connection.disconnect();
+    connection = dof_handler.get_triangulation().signals.cell_weight.connect(
+      make_weighting_callback(dof_handler, weighting_function));
   }
 
 
+
+  // ---------- static member functions ----------
+
+  // ---------- selection of weighting functions ----------
+
   template <int dim, int spacedim>
-  void
-  CellWeights<dim, spacedim>::register_ndofs_weighting(
-    const unsigned int factor)
+  typename CellWeights<dim, spacedim>::WeightingFunction
+  CellWeights<dim, spacedim>::constant_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; };
+    return
+      [factor](const typename hp::DoFHandler<dim, spacedim>::cell_iterator &,
+               const FiniteElement<dim, spacedim> &) -> unsigned int {
+        return factor;
+      };
   }
 
 
+
   template <int dim, int spacedim>
-  void
-  CellWeights<dim, spacedim>::register_ndofs_squared_weighting(
-    const unsigned int factor)
+  typename CellWeights<dim, spacedim>::WeightingFunction
+  CellWeights<dim, spacedim>::ndofs_weighting(
+    const std::pair<float, float> &coefficients)
   {
-    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;
+    return [&coefficients](
+             const typename hp::DoFHandler<dim, spacedim>::cell_iterator &,
+             const FiniteElement<dim, spacedim> &future_fe) -> unsigned int {
+      const float result =
+        std::trunc(coefficients.first *
+                   std::pow(future_fe.dofs_per_cell, coefficients.second));
+
+      Assert(result >= 0 && result <= std::numeric_limits<unsigned int>::max(),
+             ExcMessage(
+               "Cannot cast determined weight for this cell to unsigned int!"));
+
+      return static_cast<unsigned int>(result);
     };
   }
 
 
+
   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)
+  typename CellWeights<dim, spacedim>::WeightingFunction
+  CellWeights<dim, spacedim>::ndofs_weighting(
+    const std::vector<std::pair<float, float>> &coefficients)
+  {
+    return [&coefficients](
+             const typename hp::DoFHandler<dim, spacedim>::cell_iterator &,
+             const FiniteElement<dim, spacedim> &future_fe) -> unsigned int {
+      float result = 0;
+      for (const auto &pair : coefficients)
+        result += pair.first * std::pow(future_fe.dofs_per_cell, pair.second);
+      result = std::trunc(result);
+
+      Assert(result >= 0 && result <= std::numeric_limits<unsigned int>::max(),
+             ExcMessage(
+               "Cannot cast determined weight for this cell to unsigned int!"));
+
+      return static_cast<unsigned int>(result);
+    };
+  }
+
+
+
+  // ---------- handling callback functions ----------
+
+  template <int dim, int spacedim>
+  std::function<unsigned int(
+    const typename Triangulation<dim, spacedim>::cell_iterator &cell,
+    const typename Triangulation<dim, spacedim>::CellStatus     status)>
+  CellWeights<dim, spacedim>::make_weighting_callback(
+    const hp::DoFHandler<dim, spacedim> &dof_handler,
+    const typename CellWeights<dim, spacedim>::WeightingFunction
+      &weighting_function)
   {
-    weighting_function = custom_function;
+    const parallel::TriangulationBase<dim, spacedim> *tria =
+      dynamic_cast<const parallel::TriangulationBase<dim, spacedim> *>(
+        &(dof_handler.get_triangulation()));
+
+    Assert(
+      tria != nullptr,
+      ExcMessage(
+        "parallel::CellWeights requires a parallel::TriangulationBase object."));
+
+    return [&dof_handler, tria, &weighting_function](
+             const typename Triangulation<dim, spacedim>::cell_iterator &cell,
+             const typename Triangulation<dim, spacedim>::CellStatus     status)
+             -> unsigned int {
+      return CellWeights<dim, spacedim>::weighting_callback(cell,
+                                                            status,
+                                                            std::cref(
+                                                              dof_handler),
+                                                            std::cref(*tria),
+                                                            weighting_function);
+    };
   }
 
 
 
   template <int dim, int spacedim>
   unsigned int
-  CellWeights<dim, spacedim>::weight_callback(
+  CellWeights<dim, spacedim>::weighting_callback(
     const typename Triangulation<dim, spacedim>::cell_iterator &cell_,
-    const typename Triangulation<dim, spacedim>::CellStatus     status)
+    const typename Triangulation<dim, spacedim>::CellStatus     status,
+    const hp::DoFHandler<dim, spacedim> &                       dof_handler,
+    const parallel::TriangulationBase<dim, spacedim> &          triangulation,
+    const typename CellWeights<dim, spacedim>::WeightingFunction
+      &weighting_function)
   {
     // 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!"));
+    AssertThrow(&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);
+      *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.
@@ -157,7 +202,7 @@ namespace parallel
             Assert(!fe_indices_children.empty(), ExcInternalError());
 
             fe_index =
-              dof_handler->get_fe_collection().find_dominated_fe_extended(
+              dof_handler.get_fe_collection().find_dominated_fe_extended(
                 fe_indices_children, /*codim=*/0);
 
             Assert(fe_index != numbers::invalid_unsigned_int,
@@ -172,7 +217,128 @@ namespace parallel
       }
 
     // Return the cell weight determined by the function of choice.
-    return weighting_function(dof_handler->get_fe(fe_index), cell);
+    return weighting_function(cell, dof_handler.get_fe(fe_index));
+  }
+
+
+
+  // ---------- deprecated functions ----------
+
+  template <int dim, int spacedim>
+  CellWeights<dim, spacedim>::CellWeights(
+    const hp::DoFHandler<dim, spacedim> &dof_handler)
+    : dof_handler(&dof_handler)
+    , triangulation(
+        dynamic_cast<const parallel::TriangulationBase<dim, spacedim> *>(
+          &(dof_handler.get_triangulation())))
+  {
+    Assert(
+      triangulation != nullptr,
+      ExcMessage(
+        "parallel::CellWeights requires a parallel::TriangulationBase object."));
+
+    register_constant_weighting();
+  }
+
+
+
+  template <int dim, int spacedim>
+  void
+  CellWeights<dim, spacedim>::register_constant_weighting(
+    const unsigned int factor)
+  {
+    connection.disconnect();
+
+    connection = triangulation->signals.cell_weight.connect(
+      [this,
+       factor](const typename Triangulation<dim, spacedim>::cell_iterator &cell,
+               const typename Triangulation<dim, spacedim>::CellStatus status)
+        -> unsigned int {
+        return this->weighting_callback(cell,
+                                        status,
+                                        std::cref(*(this->dof_handler)),
+                                        std::cref(*(this->triangulation)),
+                                        this->constant_weighting(factor));
+      });
+  }
+
+
+
+  template <int dim, int spacedim>
+  void
+  CellWeights<dim, spacedim>::register_ndofs_weighting(
+    const unsigned int factor)
+  {
+    connection.disconnect();
+
+    connection = triangulation->signals.cell_weight.connect(
+      [this,
+       factor](const typename Triangulation<dim, spacedim>::cell_iterator &cell,
+               const typename Triangulation<dim, spacedim>::CellStatus status)
+        -> unsigned int {
+        return this->weighting_callback(cell,
+                                        status,
+                                        std::cref(*(this->dof_handler)),
+                                        std::cref(*(this->triangulation)),
+                                        this->ndofs_weighting({factor, 1}));
+      });
+  }
+
+
+
+  template <int dim, int spacedim>
+  void
+  CellWeights<dim, spacedim>::register_ndofs_squared_weighting(
+    const unsigned int factor)
+  {
+    connection.disconnect();
+
+    connection = triangulation->signals.cell_weight.connect(
+      [this,
+       factor](const typename Triangulation<dim, spacedim>::cell_iterator &cell,
+               const typename Triangulation<dim, spacedim>::CellStatus status)
+        -> unsigned int {
+        return this->weighting_callback(cell,
+                                        status,
+                                        std::cref(*(this->dof_handler)),
+                                        std::cref(*(this->triangulation)),
+                                        this->ndofs_weighting({factor, 2}));
+      });
+  }
+
+
+
+  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)
+  {
+    connection.disconnect();
+
+    const std::function<unsigned int(
+      const typename hp::DoFHandler<dim, spacedim>::cell_iterator &,
+      const FiniteElement<dim, spacedim> &)>
+      converted_function =
+        [&custom_function](
+          const typename hp::DoFHandler<dim, spacedim>::cell_iterator &cell,
+          const FiniteElement<dim, spacedim> &future_fe) -> unsigned int {
+      return custom_function(future_fe, cell);
+    };
+
+    connection = triangulation->signals.cell_weight.connect(
+      [this, converted_function](
+        const typename Triangulation<dim, spacedim>::cell_iterator &cell,
+        const typename Triangulation<dim, spacedim>::CellStatus     status)
+        -> unsigned int {
+        return this->weighting_callback(cell,
+                                        status,
+                                        std::cref(*(this->dof_handler)),
+                                        std::cref(*(this->triangulation)),
+                                        converted_function);
+      });
   }
 } // namespace parallel
 
index c5e1678176b17c59454d4b0a0cd1ab417468e70f..023d9315552bdfcece034a1773656671b438ca98 100644 (file)
@@ -1587,6 +1587,8 @@ namespace hp
     const Triangulation<dim, spacedim> &   tria,
     const hp::FECollection<dim, spacedim> &fe)
   {
+    clear();
+
     if (this->tria != &tria)
       {
         for (auto &connection : tria_listeners)
index c85962f28da343d053be126916bca08a725cbfbb..3dc5dec6b2edbd685e340b82816b9971f9c155e4 100644 (file)
@@ -65,8 +65,8 @@ test()
       }
 
   // ----- transfer -----
-  parallel::CellWeights<dim> cell_weights(dh);
-  cell_weights.register_ndofs_weighting(100000);
+  const parallel::CellWeights<dim> cell_weights(
+    dh, parallel::CellWeights<dim>::ndofs_weighting({100000, 1}));
 
   tria.repartition();
 
index 1bead3b9ab1d287bf261b5026a80d74ddbeb25d2..8e924f566fc1a3dc78cefe0dc470a773afaa037d 100644 (file)
@@ -75,8 +75,8 @@ test()
   }
 
 
-  parallel::CellWeights<dim> cell_weights(dh);
-  cell_weights.register_ndofs_weighting(100000);
+  const parallel::CellWeights<dim> cell_weights(
+    dh, parallel::CellWeights<dim>::ndofs_weighting({100000, 1}));
 
   tria.repartition();
 
@@ -91,6 +91,23 @@ test()
     deallog << "  Cumulative dofs per cell: " << dof_counter << std::endl;
   }
 
+#ifdef DEBUG
+  parallel::distributed::Triangulation<dim> other_tria(MPI_COMM_WORLD);
+  GridGenerator::hyper_cube(other_tria);
+  other_tria.refine_global(2);
+
+  dh.initialize(other_tria, fe_collection);
+
+  try
+    {
+      tria.repartition();
+    }
+  catch (ExcMessage &)
+    {
+      deallog << "Triangulation changed" << std::endl;
+    }
+#endif
+
   // make sure no processor is hanging
   MPI_Barrier(MPI_COMM_WORLD);
 
index e04bfd18cd3e481431e008775663fbef9089bcc9..5f8c76b45b5ddadf6fc089b50e21935ab22f5e13 100644 (file)
@@ -3,21 +3,25 @@ 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::Triangulation changed
 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::Triangulation changed
 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::Triangulation changed
 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::Triangulation changed
 DEAL:1:3d::OK
 
index 662cf1de04f011fde040d4975c66ddf217bb3f00..5767e9b7f537cd8d628072b361766f5bc744ef0d 100644 (file)
@@ -81,8 +81,8 @@ test()
   }
 
 
-  parallel::CellWeights<dim> cell_weights(dh);
-  cell_weights.register_ndofs_weighting(100000);
+  const parallel::CellWeights<dim> cell_weights(
+    dh, parallel::CellWeights<dim>::ndofs_weighting({100000, 1}));
 
   tria.repartition();
 
@@ -97,6 +97,23 @@ test()
     deallog << "  Cumulative dofs per cell: " << dof_counter << std::endl;
   }
 
+#ifdef DEBUG
+  parallel::distributed::Triangulation<dim> other_tria(MPI_COMM_WORLD);
+  GridGenerator::hyper_cube(other_tria);
+  other_tria.refine_global(3);
+
+  dh.initialize(other_tria, fe_collection);
+
+  try
+    {
+      tria.repartition();
+    }
+  catch (ExcMessage &)
+    {
+      deallog << "Triangulation changed" << std::endl;
+    }
+#endif
+
   // make sure no processor is hanging
   MPI_Barrier(MPI_COMM_WORLD);
 
index bf02d552a17c9f429751ab42598803a37d138c27..c45b65557593891cbd8ff2f84a1c6c1fab39891b 100644 (file)
@@ -3,22 +3,26 @@ 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::Triangulation changed
 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::Triangulation changed
 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::Triangulation changed
 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::Triangulation changed
 DEAL:1:3d::OK
 
 
@@ -26,10 +30,12 @@ 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::Triangulation changed
 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::Triangulation changed
 DEAL:2:3d::OK
 
index 1017803178df24a464c724418c5f1751f6ae9ba5..4cdc621b937d1c61b7a12d8ff62b92a1063151d3 100644 (file)
@@ -81,8 +81,8 @@ test()
   }
 
 
-  parallel::CellWeights<dim> cell_weights(dh);
-  cell_weights.register_ndofs_weighting(100000);
+  const parallel::CellWeights<dim> cell_weights(
+    dh, parallel::CellWeights<dim>::ndofs_weighting({100000, 1}));
 
   // we didn't mark any cells, but we want to repartition our domain
   tria.execute_coarsening_and_refinement();
@@ -98,6 +98,27 @@ test()
     deallog << "  Cumulative dofs per cell: " << dof_counter << std::endl;
   }
 
+#ifdef DEBUG
+  parallel::shared::Triangulation<dim> other_tria(
+    MPI_COMM_WORLD,
+    ::Triangulation<dim>::none,
+    false,
+    parallel::shared::Triangulation<dim>::Settings::partition_metis);
+  GridGenerator::hyper_cube(other_tria);
+  other_tria.refine_global(3);
+
+  dh.initialize(other_tria, fe_collection);
+
+  try
+    {
+      tria.execute_coarsening_and_refinement();
+    }
+  catch (ExcMessage &)
+    {
+      deallog << "Triangulation changed" << std::endl;
+    }
+#endif
+
   // make sure no processor is hanging
   MPI_Barrier(MPI_COMM_WORLD);
 
index d1d118c4663310ad76118ea90ead0e7c640a3ed3..7f1e8bb661022609e492d5947b0a7eb7e78f3796 100644 (file)
@@ -3,21 +3,25 @@ DEAL:0:2d::Number of cells before repartitioning: 8
 DEAL:0:2d::  Cumulative dofs per cell: 32
 DEAL:0:2d::Number of cells after repartitioning: 13
 DEAL:0:2d::  Cumulative dofs per cell: 52
+DEAL:0:2d::Triangulation changed
 DEAL:0:2d::OK
 DEAL:0:3d::Number of cells before repartitioning: 33
 DEAL:0:3d::  Cumulative dofs per cell: 264
 DEAL:0:3d::Number of cells after repartitioning: 47
 DEAL:0:3d::  Cumulative dofs per cell: 376
+DEAL:0:3d::Triangulation changed
 DEAL:0:3d::OK
 
 DEAL:1:2d::Number of cells before repartitioning: 8
 DEAL:1:2d::  Cumulative dofs per cell: 64
 DEAL:1:2d::Number of cells after repartitioning: 3
 DEAL:1:2d::  Cumulative dofs per cell: 44
+DEAL:1:2d::Triangulation changed
 DEAL:1:2d::OK
 DEAL:1:3d::Number of cells before repartitioning: 31
 DEAL:1:3d::  Cumulative dofs per cell: 456
 DEAL:1:3d::Number of cells after repartitioning: 17
 DEAL:1:3d::  Cumulative dofs per cell: 344
+DEAL:1:3d::Triangulation changed
 DEAL:1:3d::OK
 
index f42d31ee46c3516f9ae1985eac2a4aa22ba5f8ab..0f26ba9e6e637bbac0debbeb85d33ff9d5e8fcc7 100644 (file)
@@ -3,21 +3,25 @@ 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: 13
 DEAL:0:2d::  Cumulative dofs per cell: 52
+DEAL:0:2d::Triangulation changed
 DEAL:0:2d::OK
 DEAL:0:3d::Number of cells before repartitioning: 32
 DEAL:0:3d::  Cumulative dofs per cell: 256
 DEAL:0:3d::Number of cells after repartitioning: 47
 DEAL:0:3d::  Cumulative dofs per cell: 376
+DEAL:0:3d::Triangulation changed
 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: 3
 DEAL:1:2d::  Cumulative dofs per cell: 44
+DEAL:1:2d::Triangulation changed
 DEAL:1:2d::OK
 DEAL:1:3d::Number of cells before repartitioning: 32
 DEAL:1:3d::  Cumulative dofs per cell: 464
 DEAL:1:3d::Number of cells after repartitioning: 17
 DEAL:1:3d::  Cumulative dofs per cell: 344
+DEAL:1:3d::Triangulation changed
 DEAL:1:3d::OK
 
index 3b0461bbea1112d5cf0cf15913dcfd56c77c5a4c..7e1818bb259475ba10fb3f7ceba5a2808407a5b9 100644 (file)
@@ -83,8 +83,8 @@ test()
   }
 
 
-  parallel::CellWeights<dim> cell_weights(dh);
-  cell_weights.register_ndofs_weighting(100000);
+  const parallel::CellWeights<dim> cell_weights(
+    dh, parallel::CellWeights<dim>::ndofs_weighting({100000, 1}));
 
   // we didn't mark any cells, but we want to repartition our domain
   tria.execute_coarsening_and_refinement();
@@ -100,6 +100,27 @@ test()
     deallog << "  Cumulative dofs per cell: " << dof_counter << std::endl;
   }
 
+#ifdef DEBUG
+  parallel::shared::Triangulation<dim> other_tria(
+    MPI_COMM_WORLD,
+    ::Triangulation<dim>::none,
+    false,
+    parallel::shared::Triangulation<dim>::Settings::partition_metis);
+  GridGenerator::hyper_cube(other_tria);
+  other_tria.refine_global(3);
+
+  dh.initialize(other_tria, fe_collection);
+
+  try
+    {
+      tria.execute_coarsening_and_refinement();
+    }
+  catch (ExcMessage &)
+    {
+      deallog << "Triangulation changed" << std::endl;
+    }
+#endif
+
   // make sure no processor is hanging
   MPI_Barrier(MPI_COMM_WORLD);
 
index d1d118c4663310ad76118ea90ead0e7c640a3ed3..7f1e8bb661022609e492d5947b0a7eb7e78f3796 100644 (file)
@@ -3,21 +3,25 @@ DEAL:0:2d::Number of cells before repartitioning: 8
 DEAL:0:2d::  Cumulative dofs per cell: 32
 DEAL:0:2d::Number of cells after repartitioning: 13
 DEAL:0:2d::  Cumulative dofs per cell: 52
+DEAL:0:2d::Triangulation changed
 DEAL:0:2d::OK
 DEAL:0:3d::Number of cells before repartitioning: 33
 DEAL:0:3d::  Cumulative dofs per cell: 264
 DEAL:0:3d::Number of cells after repartitioning: 47
 DEAL:0:3d::  Cumulative dofs per cell: 376
+DEAL:0:3d::Triangulation changed
 DEAL:0:3d::OK
 
 DEAL:1:2d::Number of cells before repartitioning: 8
 DEAL:1:2d::  Cumulative dofs per cell: 64
 DEAL:1:2d::Number of cells after repartitioning: 3
 DEAL:1:2d::  Cumulative dofs per cell: 44
+DEAL:1:2d::Triangulation changed
 DEAL:1:2d::OK
 DEAL:1:3d::Number of cells before repartitioning: 31
 DEAL:1:3d::  Cumulative dofs per cell: 456
 DEAL:1:3d::Number of cells after repartitioning: 17
 DEAL:1:3d::  Cumulative dofs per cell: 344
+DEAL:1:3d::Triangulation changed
 DEAL:1:3d::OK
 
index f42d31ee46c3516f9ae1985eac2a4aa22ba5f8ab..0f26ba9e6e637bbac0debbeb85d33ff9d5e8fcc7 100644 (file)
@@ -3,21 +3,25 @@ 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: 13
 DEAL:0:2d::  Cumulative dofs per cell: 52
+DEAL:0:2d::Triangulation changed
 DEAL:0:2d::OK
 DEAL:0:3d::Number of cells before repartitioning: 32
 DEAL:0:3d::  Cumulative dofs per cell: 256
 DEAL:0:3d::Number of cells after repartitioning: 47
 DEAL:0:3d::  Cumulative dofs per cell: 376
+DEAL:0:3d::Triangulation changed
 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: 3
 DEAL:1:2d::  Cumulative dofs per cell: 44
+DEAL:1:2d::Triangulation changed
 DEAL:1:2d::OK
 DEAL:1:3d::Number of cells before repartitioning: 32
 DEAL:1:3d::  Cumulative dofs per cell: 464
 DEAL:1:3d::Number of cells after repartitioning: 17
 DEAL:1:3d::  Cumulative dofs per cell: 344
+DEAL:1:3d::Triangulation changed
 DEAL:1: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.