]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Make 'CellWeights' available for any 'p::Triangulation'. 7305/head
authorMarc Fehling <marc.fehling@gmx.net>
Mon, 8 Oct 2018 22:29:13 +0000 (16:29 -0600)
committerMarc Fehling <marc.fehling@gmx.net>
Wed, 10 Oct 2018 00:42:56 +0000 (18:42 -0600)
include/deal.II/distributed/cell_weights.h
source/distributed/cell_weights.cc
source/distributed/cell_weights.inst.in
tests/mpi/hp_cell_weights_01.cc
tests/mpi/hp_cell_weights_02.cc
tests/mpi/hp_cell_weights_03.cc [new file with mode: 0644]
tests/mpi/hp_cell_weights_03.with_mpi=true.with_metis=true.mpirun=2.output [new file with mode: 0644]

index 904ffd3672e9446decae452608bc2954dcbcf8d1..f41765bfa99d36f2c8f236807ae21abadc95da94 100644 (file)
@@ -18,7 +18,7 @@
 
 #include <deal.II/base/config.h>
 
-#include <deal.II/distributed/tria.h>
+#include <deal.II/distributed/tria_base.h>
 
 #include <deal.II/hp/dof_handler.h>
 
@@ -27,166 +27,161 @@ DEAL_II_NAMESPACE_OPEN
 
 namespace parallel
 {
-  namespace distributed
+  /**
+   * Anytime a parallel::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::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::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:
     /**
-     * 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.
+     * Constructor.
      *
-     * 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.
+     * @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.
      *
-     * An object of this class needs to exist for every DoFHandler associated
-     * with the Triangulation we work on to achieve satisfying work balancing
-     * results.
+     * @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.
      *
-     * 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.
+     * We store both to make sure to always work on the correct combination of
+     * both.
+     */
+    SmartPointer<const parallel::Triangulation<dim, spacedim>, CellWeights>
+      triangulation;
+
+    /**
+     * Function that will determine each cell's weight.
      *
-     * @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.
+     * Can be set using the register_constant_weighting(),
+     * register_ndofs_weighting(), register_ndofs_squared_weighting(), and
+     * register_custom_weighting() member functions.
      *
-     * @ingroup distributed
-     * @author Marc Fehling, 2018
+     * 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.
      */
-    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
+    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 parallel
 
 
index 3c29bf7156680380178c9365b63f9d797449ac91..c031c4654303f459927a7f7253a877ab7ad04324 100644 (file)
@@ -14,8 +14,6 @@
 // ---------------------------------------------------------------------
 
 
-#include <deal.II/base/config.h>
-
 #include <deal.II/distributed/cell_weights.h>
 
 #include <deal.II/dofs/dof_accessor.h>
@@ -25,159 +23,152 @@ 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::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::CellWeights requires a parallel::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)
   {
-    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_in_collection(
-                             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
+    // 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 Triangulation<dim, spacedim>::CELL_PERSIST:
+        case Triangulation<dim, spacedim>::CELL_REFINE:
+        case Triangulation<dim, spacedim>::CELL_INVALID:
+          fe_index = cell->active_fe_index();
+          break;
+
+        case 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_in_collection(
+                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 parallel
 
 
index b630293009f5e53462c916be100fd4462bb30514..4d0ffae2be391825e57574cdbb2ade0c494f0485 100644 (file)
@@ -19,11 +19,8 @@ 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>;
+      template class CellWeights<deal_II_dimension, deal_II_space_dimension>;
 #endif
-      \}
     \}
   }
index e99980c08a99625cd0a97760cc2e0f59fccffc4e..838e48cc85f51ae5bdd225ebaef0cfaaec734ce0 100644 (file)
@@ -26,6 +26,8 @@
 // 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 works on a parallel::distributed::Triangulation.
 
 
 #include <deal.II/distributed/active_fe_indices_transfer.h>
@@ -67,7 +69,7 @@ test()
   parallel::distributed::ActiveFEIndicesTransfer<dim> feidx_transfer(dh);
   feidx_transfer.prepare_for_transfer();
 
-  parallel::distributed::CellWeights<dim> cell_weights(dh);
+  parallel::CellWeights<dim> cell_weights(dh);
   cell_weights.register_ndofs_weighting(100000);
 
 
index 81e00a44ea9a5c5d1ed1b364de062aa3734e2cfb..52b8a6608007a6be3b9cf3acbc6db80623478e9d 100644 (file)
@@ -27,6 +27,8 @@
 // We employ a large proportionality factor on our weighting function
 // to neglect the standard weight of '1000' per cell.
 //
+// This test works on a parallel::distributed::Triangulation.
+//
 // 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
@@ -73,7 +75,7 @@ test()
   parallel::distributed::ActiveFEIndicesTransfer<dim> feidx_transfer(dh);
   feidx_transfer.prepare_for_transfer();
 
-  parallel::distributed::CellWeights<dim> cell_weights(dh);
+  parallel::CellWeights<dim> cell_weights(dh);
   cell_weights.register_ndofs_weighting(100000);
 
 
diff --git a/tests/mpi/hp_cell_weights_03.cc b/tests/mpi/hp_cell_weights_03.cc
new file mode 100644 (file)
index 0000000..c2eccc2
--- /dev/null
@@ -0,0 +1,123 @@
+// ---------------------------------------------------------------------
+//
+// 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.
+//
+// This test works on a parallel::shared::Triangulation with METIS
+// as a partitioner. Cell weighting with ZOLTAN was not available
+// during the time this test was written.
+
+
+#include <deal.II/distributed/cell_weights.h>
+#include <deal.II/distributed/shared_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::shared::Triangulation<dim> tria(
+    MPI_COMM_WORLD,
+    ::Triangulation<dim>::none,
+    false,
+    parallel::shared::Triangulation<dim>::Settings::partition_metis);
+
+  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::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;
+  }
+
+
+  // we didn't mark any cells, but we want to repartition our domain
+  tria.execute_coarsening_and_refinement();
+  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_03.with_mpi=true.with_metis=true.mpirun=2.output b/tests/mpi/hp_cell_weights_03.with_mpi=true.with_metis=true.mpirun=2.output
new file mode 100644 (file)
index 0000000..d1d118c
--- /dev/null
@@ -0,0 +1,23 @@
+
+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::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::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::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::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.