]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Reworked CoarseningStrategies and moved it into a separate header.
authorMarc Fehling <marc.fehling@gmx.net>
Wed, 10 Jul 2019 14:53:37 +0000 (16:53 +0200)
committerMarc Fehling <marc.fehling@gmx.net>
Wed, 10 Jul 2019 14:53:43 +0000 (16:53 +0200)
Deprecated the old p::d::CellDataTransfer::CoarseningStrategies struct.

Defined exceptions that recur in code.

Minor improvements to p::d::CellDataTransfer and hp::DoFHandler.

13 files changed:
doc/news/changes/incompatibilities/20190626Fehling [new file with mode: 0644]
include/deal.II/distributed/cell_data_transfer.h
include/deal.II/distributed/cell_data_transfer.templates.h
include/deal.II/grid/tria.h
include/deal.II/hp/dof_handler.h
include/deal.II/hp/fe_collection.h
include/deal.II/numerics/coarsening_strategies.h [new file with mode: 0644]
source/distributed/cell_data_transfer.inst.in
source/distributed/cell_weights.cc
source/distributed/solution_transfer.cc
source/hp/dof_handler.cc
source/numerics/solution_transfer.cc
tests/mpi/cell_data_transfer_03.cc

diff --git a/doc/news/changes/incompatibilities/20190626Fehling b/doc/news/changes/incompatibilities/20190626Fehling
new file mode 100644 (file)
index 0000000..f7f159c
--- /dev/null
@@ -0,0 +1,10 @@
+Changed: The CoarseningStrategies struct has been moved out of the
+parallel::distributed::CellDataTransfer class into a separate header and
+is now treated as a namespace. Its static member functions are now free
+functions and only take a `std::vector` as a parameter that contains all
+the data from the children. Therefore, the `coarsening_strategy`
+parameter for the constructor of the
+parallel::distributed::CellDataTransfer class has to be adjusted
+accordingly as well.
+<br>
+(Marc Fehling, 2019/06/26)
index af04d75becb4661df3cc6e53c2d085d4ef85d80d..1fffffedc0ef4bf51cce569bbcaed1a015f05c83 100644 (file)
@@ -22,6 +22,8 @@
 
 #include <deal.II/distributed/tria.h>
 
+#include <deal.II/numerics/coarsening_strategies.h>
+
 #include <boost/range/iterator_range.hpp>
 
 #include <functional>
@@ -43,11 +45,11 @@ namespace parallel
      * freedom defined on a parallel::distributed::Triangulation.
      *
      * This class has been designed to operate on any kind of datatype that is
-     * serializable. A non-distributed Vector container has to be provided,
-     * which holds the cell-wise data in the same order as active cells are
-     * traversed. In other words, each entry corresponds to the cell with the
-     * same index CellAccessor::active_cell_index(), and the container has to be
-     * of size Triangulation::n_active_cells().
+     * serializable. A non-distributed container (like Vector or `std::vector`)
+     * has to be provided, which holds the cell-wise data in the same order as
+     * active cells are traversed. In other words, each entry corresponds to the
+     * cell with the same index CellAccessor::active_cell_index(), and the
+     * container has to be of size Triangulation::n_active_cells().
      *
      * <h3>Transferring cell-wise data</h3>
      *
@@ -60,7 +62,7 @@ namespace parallel
      *
      * // prepare the CellDataTransfer object for coarsening and refinement
      * // and give the cell data vector that we intend to unpack later,
-     * Vector<double> data_to_transfer (triangulation.n_active_cells());
+     * Vector<double> data_to_transfer(triangulation.n_active_cells());
      * //[fill data_to_transfer with cell-wise values...]
      *
      * parallel::distributed::CellDataTransfer<dim, spacedim, Vector<double>>
@@ -71,7 +73,7 @@ namespace parallel
      * triangulation.execute_coarsening_and_refinement();
      *
      * // unpack transferred data,
-     * Vector<double> transferred_data (triangulation.n_active_cells());
+     * Vector<double> transferred_data(triangulation.n_active_cells());
      * cell_data_trans.unpack(transferred_data);
      *
      * @endcode
@@ -85,7 +87,7 @@ namespace parallel
      * For serialization, the following code snippet saves not only the
      * triangulation itself, but also the cell-wise data attached:
      * @code
-     * Vector<double> data_to_transfer (triangulation.n_active_cells());
+     * Vector<double> data_to_transfer(triangulation.n_active_cells());
      * //[fill data_to_transfer with cell-wise values...]
      *
      * parallel::distributed::CellDataTransfer<dim, spacedim, Vector<double>>
@@ -103,7 +105,7 @@ namespace parallel
      *
      * parallel::distributed::CellDataTransfer<dim, spacedim, Vector<double>>
      *   cell_data_trans(triangulation);
-     * Vector<double> transferred_data (triangulation.n_active_cells());
+     * Vector<double> transferred_data(triangulation.n_active_cells());
      * cell_data_trans.deserialize(transferred_data);
      * @endcode
      *
@@ -125,23 +127,26 @@ namespace parallel
     template <int dim, int spacedim = dim, typename VectorType = Vector<double>>
     class CellDataTransfer
     {
+    private:
+      /**
+       * An alias that defines the data type of provided container template.
+       */
+      using value_type = typename VectorType::value_type;
+
     public:
       /**
-       * When data is transferred during coarsening, it is not trivial to decide
-       * how to handle data of child cells which will be coarsened. Or in other
-       * words, which data should be stored in the corresponding parent cell.
+       * @copydoc dealii::CoarseningStrategies
        *
-       * In this namespace, we offer a few strategies that cope with this
-       * problem. Such strategies can be passed to the
-       * parallel::distributed::CellDataTransfer constructor.
+       * @deprecated Use dealii::CoarseningStrategies instead.
        */
-      struct CoarseningStrategies
+      struct DEAL_II_DEPRECATED CoarseningStrategies
       {
         /**
-         * Evaluate data from @p input_vector on children of @p parent.
-         * Check if data on all children match, and return that value.
+         * @copydoc dealii::CoarseningStrategies::check_equality
+         *
+         * @deprecated Use dealii::CoarseningStrategies::check_equality() instead.
          */
-        static typename VectorType::value_type
+        DEAL_II_DEPRECATED static typename VectorType::value_type
         check_equality(const typename parallel::distributed::
                          Triangulation<dim, spacedim>::cell_iterator &parent,
                        const VectorType &input_vector)
@@ -161,10 +166,11 @@ namespace parallel
         }
 
         /**
-         * Evaluate data from @p input_vector on children of @p parent.
-         * Return sum.
+         * @copydoc dealii::CoarseningStrategies::sum
+         *
+         * @deprecated Use dealii::CoarseningStrategies::sum() instead.
          */
-        static typename VectorType::value_type
+        DEAL_II_DEPRECATED static typename VectorType::value_type
         sum(const typename parallel::distributed::Triangulation<dim, spacedim>::
               cell_iterator & parent,
             const VectorType &input_vector)
@@ -180,10 +186,11 @@ namespace parallel
         }
 
         /**
-         * Evaluate data from @p input_vector on children of @p parent.
-         * Return mean value.
+         * @copydoc dealii::CoarseningStrategies::mean
+         *
+         * @deprecated Use dealii::CoarseningStrategies::mean() instead.
          */
-        static typename VectorType::value_type
+        DEAL_II_DEPRECATED static typename VectorType::value_type
         mean(const typename parallel::distributed::
                Triangulation<dim, spacedim>::cell_iterator &parent,
              const VectorType &                             input_vector)
@@ -205,21 +212,35 @@ namespace parallel
        * @param[in] coarsening_strategy Function deciding which data to store on
        *   a cell whose children will get coarsened into.
        */
+      CellDataTransfer(
+        const parallel::distributed::Triangulation<dim, spacedim>
+          &        triangulation,
+        const bool transfer_variable_size_data = false,
+        const std::function<value_type(
+          const std::vector<value_type> &children_values)> coarsening_strategy =
+          &dealii::CoarseningStrategies::check_equality<value_type>);
+
+      /**
+       * @copydoc CellDataTransfer::CellDataTransfer
+       *
+       * @copydoc Use the above constructor instead.
+       */
+      DEAL_II_DEPRECATED
       CellDataTransfer(const parallel::distributed::Triangulation<dim, spacedim>
                          &        triangulation,
-                       const bool transfer_variable_size_data = false,
-                       const std::function<typename VectorType::value_type(
+                       const bool transfer_variable_size_data,
+                       const std::function<value_type(
                          const typename parallel::distributed::
                            Triangulation<dim, spacedim>::cell_iterator &parent,
-                         const VectorType &input_vector)> coarsening_strategy =
-                         &CoarseningStrategies::check_equality);
+                         const VectorType &input_vector)> coarsening_strategy);
 
       /**
-       * Prepare the current object for coarsening and refinement. It registers
-       * the data transfer of @p in on the underlying triangulation. @p in includes
-       * data to be interpolated onto the new (refined and/or coarsened) grid.
-       * See documentation of this class for more information on how to use this
-       * functionality.
+       * Prepare the current object for coarsening and refinement.
+       *
+       * It registers the data transfer of @p in on the underlying triangulation.
+       * @p in includes data to be interpolated onto the new (refined and/or
+       * coarsened) grid. See documentation of this class for more information
+       * on how to use this functionality.
        *
        * This function can be called only once for the specified container
        * until data transfer has been completed. If multiple vectors shall be
@@ -236,9 +257,10 @@ namespace parallel
         const std::vector<const VectorType *> &all_in);
 
       /**
-       * Prepare the serialization of the given vector. The serialization is
-       * done by Triangulation::save(). See documentation of this class for more
-       * information on how to use this functionality.
+       * Prepare the serialization of the given vector.
+       *
+       * The serialization is done by Triangulation::save(). See documentation
+       * of this class for more information on how to use this functionality.
        *
        * This function can be called only once for the specified container
        * until data transfer has been completed. If multiple vectors shall be
@@ -293,13 +315,11 @@ namespace parallel
       const bool transfer_variable_size_data;
 
       /**
-       * Function deciding which data to store from @p input_vector on a
-       * cell @p parent, whose children will get coarsened into.
+       * Function deciding on how to process data from children to be stored on
+       * the parent cell.
        */
-      const std::function<typename VectorType::value_type(
-        const typename parallel::distributed::Triangulation<dim, spacedim>::
-          cell_iterator & parent,
-        const VectorType &input_vector)>
+      const std::function<value_type(
+        const std::vector<value_type> &children_values)>
         coarsening_strategy;
 
       /**
@@ -353,4 +373,4 @@ namespace parallel
 
 DEAL_II_NAMESPACE_CLOSE
 
-#endif
+#endif /* dealii_distributed_cell_data_transfer_h */
index 5e66849ef5b56a50895d1ed1fd045b5ac2cc17f0..08523d833b3f07872b539912ef0c8dfd810e3c34 100644 (file)
@@ -63,17 +63,60 @@ namespace parallel
     template <int dim, int spacedim, typename VectorType>
     CellDataTransfer<dim, spacedim, VectorType>::CellDataTransfer(
       const parallel::distributed::Triangulation<dim, spacedim> &triangulation,
-      const bool                         transfer_variable_size_data,
-      const std::function<typename VectorType::value_type(
-        const typename parallel::distributed::Triangulation<dim, spacedim>::
-          cell_iterator & parent,
-        const VectorType &input_vector)> coarsening_strategy)
+      const bool transfer_variable_size_data,
+      const std::function<value_type(
+        const std::vector<value_type> &children_values)> coarsening_strategy)
       : triangulation(&triangulation, typeid(*this).name())
       , transfer_variable_size_data(transfer_variable_size_data)
       , coarsening_strategy(coarsening_strategy)
       , handle(numbers::invalid_unsigned_int)
     {}
 
+    template <int dim, int spacedim, typename VectorType>
+    DEAL_II_DEPRECATED
+    CellDataTransfer<dim, spacedim, VectorType>::CellDataTransfer(
+      const parallel::distributed::Triangulation<dim, spacedim> &triangulation,
+      const bool transfer_variable_size_data,
+      const std::function<
+        value_type(const typename parallel::distributed::
+                     Triangulation<dim, spacedim>::cell_iterator &parent,
+                   const VectorType &input_vector)> coarsening_strategy)
+      : triangulation(&triangulation, typeid(*this).name())
+      , transfer_variable_size_data(transfer_variable_size_data)
+      , handle(numbers::invalid_unsigned_int)
+    {
+      value_type (*const *old_strategy)(
+        const typename parallel::distributed::Triangulation<dim, spacedim>::
+          cell_iterator &,
+        const VectorType &) =
+        coarsening_strategy.template target<
+          value_type (*)(const typename parallel::distributed::
+                           Triangulation<dim, spacedim>::cell_iterator &,
+                         const VectorType &)>();
+
+      if (*old_strategy == CoarseningStrategies::check_equality)
+        const_cast<
+          std::function<value_type(const std::vector<value_type> &)> &>(
+          this->coarsening_strategy) =
+          &dealii::CoarseningStrategies::check_equality<value_type>;
+      else if (*old_strategy == CoarseningStrategies::sum)
+        const_cast<
+          std::function<value_type(const std::vector<value_type> &)> &>(
+          this->coarsening_strategy) =
+          &dealii::CoarseningStrategies::sum<value_type>;
+      else if (*old_strategy == CoarseningStrategies::mean)
+        const_cast<
+          std::function<value_type(const std::vector<value_type> &)> &>(
+          this->coarsening_strategy) =
+          &dealii::CoarseningStrategies::mean<value_type>;
+      else
+        Assert(
+          false,
+          ExcMessage(
+            "The constructor using the former function type of the "
+            "coarsening_strategy parameter is no longer supported. Please use "
+            "the latest function type instead"));
+    }
 
 
     // Interface for packing
@@ -226,8 +269,7 @@ namespace parallel
                                                           spacedim>::CellStatus
         status)
     {
-      std::vector<typename VectorType::value_type> cell_data(
-        input_vectors.size());
+      std::vector<value_type> cell_data(input_vectors.size());
 
       // Extract data from input_vectors for this particular cell.
       auto it_input  = input_vectors.cbegin();
@@ -247,10 +289,26 @@ namespace parallel
 
               case parallel::distributed::Triangulation<dim,
                                                         spacedim>::CELL_COARSEN:
-                // Cell is parent whose children will get coarsened to.
-                // Decide data to store on parent by provided strategy.
-                *it_output = coarsening_strategy(cell, *(*it_input));
-                break;
+                {
+                  // Cell is parent whose children will get coarsened to.
+                  // Decide data to store on parent by provided strategy.
+                  std::vector<value_type> children_values(cell->n_children());
+                  for (unsigned int child_index = 0;
+                       child_index < cell->n_children();
+                       ++child_index)
+                    {
+                      const auto child = cell->child(child_index);
+                      Assert(child->active() && child->coarsen_flag_set(),
+                             typename dealii::Triangulation<
+                               dim>::ExcInconsistentCoarseningFlags());
+
+                      children_values[child_index] =
+                        (**it_input)[child->active_cell_index()];
+                    }
+
+                  *it_output = coarsening_strategy(children_values);
+                  break;
+                }
 
               default:
                 Assert(false, ExcInternalError());
@@ -281,21 +339,20 @@ namespace parallel
         &                        data_range,
       std::vector<VectorType *> &all_out)
     {
-      std::vector<typename VectorType::value_type> cell_data;
+      std::vector<value_type> cell_data;
 
       // We have to unpack the corresponding datatype that has been packed
       // beforehand.
       if (all_out.size() == 1)
-        cell_data.push_back(Utilities::unpack<typename VectorType::value_type>(
+        cell_data.push_back(Utilities::unpack<value_type>(
           data_range.begin(),
           data_range.end(),
           /*allow_compression=*/transfer_variable_size_data));
       else
-        cell_data =
-          Utilities::unpack<std::vector<typename VectorType::value_type>>(
-            data_range.begin(),
-            data_range.end(),
-            /*allow_compression=*/transfer_variable_size_data);
+        cell_data = Utilities::unpack<std::vector<value_type>>(
+          data_range.begin(),
+          data_range.end(),
+          /*allow_compression=*/transfer_variable_size_data);
 
       // Check if sizes match.
       Assert(cell_data.size() == all_out.size(), ExcInternalError());
index 88875403bef3e11f6dcba7a732b634fe4d5e14c2..515db001e7fffe43f0f32312523194e29a225317 100644 (file)
@@ -3556,6 +3556,18 @@ public:
                  << "The given boundary_id " << arg1
                  << " is not defined in this Triangulation!");
 
+  /**
+   * Exception
+   *
+   * @ingroup Exceptions
+   */
+  DeclExceptionMsg(
+    ExcInconsistentCoarseningFlags,
+    "A cell is flagged for coarsening, but either not all of its siblings "
+    "are active or flagged for coarsening as well. Please clean up all "
+    "coarsen flags on your triangulation via "
+    "Triangulation::prepare_coarsening_and_refinement() beforehand!");
+
   /*
    * @}
    */
index 3a042711b2dc0e20e4732848aa78883a1306914c..abb18fdbbcaa92e2cb5032a6d3195adb1091100b 100644 (file)
@@ -26,8 +26,6 @@
 #include <deal.II/base/smartpointer.h>
 #include <deal.II/base/template_constraints.h>
 
-#include <deal.II/distributed/cell_data_transfer.templates.h>
-
 #include <deal.II/dofs/deprecated_function_map.h>
 #include <deal.II/dofs/dof_accessor.h>
 #include <deal.II/dofs/dof_iterator_selector.h>
@@ -47,6 +45,15 @@ DEAL_II_NAMESPACE_OPEN
 template <int dim, int spacedim>
 class Triangulation;
 
+namespace parallel
+{
+  namespace distributed
+  {
+    template <int dim, int spacedim, typename VectorType>
+    class CellDataTransfer;
+  }
+} // namespace parallel
+
 namespace internal
 {
   namespace DoFHandlerImplementation
index 392f872db41d1a970fc0545a9d9c331f809f4d8b..a422679dd4974d351445b80644858c7435edcaeb 100644 (file)
@@ -768,12 +768,32 @@ namespace hp
     BlockMask
     block_mask(const ComponentMask &component_mask) const;
 
+    /**
+     * @name Exceptions
+     * @{
+     */
 
     /**
      * Exception
+     *
+     * @ingroup Exceptions
      */
     DeclException0(ExcNoFiniteElements);
 
+    /**
+     * Exception
+     *
+     * @ingroup Exceptions
+     */
+    DeclExceptionMsg(
+      ExcNoDominatedFiniteElementAmongstChildren,
+      "No FiniteElement has been found in your FECollection that is "
+      "dominated by all children of a cell you are trying to coarsen!");
+
+    /**
+     * @}
+     */
+
   private:
     /**
      * Array of pointers to the finite elements stored by this collection.
diff --git a/include/deal.II/numerics/coarsening_strategies.h b/include/deal.II/numerics/coarsening_strategies.h
new file mode 100644 (file)
index 0000000..0a29483
--- /dev/null
@@ -0,0 +1,139 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2019 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_coarsening_strategies_h
+#define dealii_coarsening_strategies_h
+
+
+#include <deal.II/base/config.h>
+
+#include <deal.II/base/exceptions.h>
+
+#include <algorithm>
+#include <numeric>
+#include <typeinfo>
+#include <vector>
+
+DEAL_II_NAMESPACE_OPEN
+
+/**
+ * When data is transferred during coarsening, it is not trivial to decide
+ * how to handle data of child cells which will be coarsened. Or in other
+ * words, which data should be stored in the corresponding parent cell.
+ *
+ * In this namespace, we offer a few strategies that cope with this
+ * problem. Such strategies can be passed to the CellDataTransfer and
+ * parallel::distributed::CellDataTransfer constructors.
+ *
+ * @author Marc Fehling, 2019
+ */
+namespace CoarseningStrategies
+{
+  /**
+   * Check if data on all children match, and return that value.
+   */
+  template <typename value_type>
+  value_type
+  check_equality(const std::vector<value_type> &children_values);
+
+  /**
+   * Return sum of data on all children.
+   */
+  template <typename value_type>
+  value_type
+  sum(const std::vector<value_type> &children_values);
+
+  /**
+   * Return mean value of data on all children.
+   */
+  template <typename value_type>
+  value_type
+  mean(const std::vector<value_type> &children_values);
+
+  /**
+   * Return maximum value of data on all children.
+   */
+  template <typename value_type>
+  value_type
+  max(const std::vector<value_type> &children_values);
+} // namespace CoarseningStrategies
+
+
+
+/* ---------------- template functions ---------------- */
+
+#ifndef DOXYGEN
+
+namespace CoarseningStrategies
+{
+  template <typename value_type>
+  value_type
+  check_equality(const std::vector<value_type> &children_values)
+  {
+    Assert(!children_values.empty(), ExcInternalError());
+
+    const auto first_child = children_values.cbegin();
+    for (auto other_child = first_child + 1;
+         other_child != children_values.cend();
+         ++other_child)
+      Assert(*first_child == *other_child,
+             ExcMessage(
+               "Values on cells that will be coarsened are not equal!"));
+
+    return *first_child;
+  }
+
+
+
+  template <typename value_type>
+  value_type
+  sum(const std::vector<value_type> &children_values)
+  {
+    static_assert(std::is_arithmetic<value_type>::value &&
+                    !std::is_same<value_type, bool>::value,
+                  "The provided value_type may not meet the requirements "
+                  "of this function.");
+
+    Assert(!children_values.empty(), ExcInternalError());
+    return std::accumulate(children_values.cbegin(),
+                           children_values.cend(),
+                           static_cast<value_type>(0));
+  }
+
+
+
+  template <typename value_type>
+  value_type
+  mean(const std::vector<value_type> &children_values)
+  {
+    return sum<value_type>(children_values) / children_values.size();
+  }
+
+
+
+  template <typename value_type>
+  value_type
+  max(const std::vector<value_type> &children_values)
+  {
+    Assert(!children_values.empty(), ExcInternalError());
+    return *std::max_element(children_values.cbegin(), children_values.cend());
+  }
+} // namespace CoarseningStrategies
+
+#endif
+
+DEAL_II_NAMESPACE_CLOSE
+
+#endif /* dealii_coarsening_strategies_h */
index 83def7d8df13003356e9584e09e4443e5a41260d..87274cc2f576b08618c6b7b3c2a2206971207425 100644 (file)
 
 
 
-for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS)
+for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS;
+     SCALAR : REAL_SCALARS)
   {
 #if deal_II_dimension <= deal_II_space_dimension
     template class parallel::distributed::CellDataTransfer<
       deal_II_dimension,
       deal_II_space_dimension,
-      Vector<double>>;
-    template class parallel::distributed::CellDataTransfer<
-      deal_II_dimension,
-      deal_II_space_dimension,
-      Vector<float>>;
+      Vector<SCALAR>>;
 #endif
   }
index 69f2c267beb9b6312cc038b756f48e9e0d9df4ea..70797886f6eec87ba0d6d7c824525080691288bd 100644 (file)
@@ -145,18 +145,23 @@ namespace parallel
             std::set<unsigned int> fe_indices_children;
             for (unsigned int child_index = 0; child_index < cell->n_children();
                  ++child_index)
-              fe_indices_children.insert(
-                cell->child(child_index)->future_fe_index());
+              {
+                const auto child = cell->child(child_index);
+                Assert(child->active() && child->coarsen_flag_set(),
+                       typename dealii::Triangulation<
+                         dim>::ExcInconsistentCoarseningFlags());
+
+                fe_indices_children.insert(child->future_fe_index());
+              }
+            Assert(!fe_indices_children.empty(), ExcInternalError());
 
             fe_index =
               dof_handler->get_fe_collection().find_dominated_fe_extended(
                 fe_indices_children, /*codim=*/0);
 
             Assert(fe_index != numbers::invalid_unsigned_int,
-                   ExcMessage(
-                     "No FiniteElement has been found in your FECollection "
-                     "that is dominated by all children of a cell you are "
-                     "trying to coarsen!"));
+                   typename dealii::hp::FECollection<
+                     dim>::ExcNoDominatedFiniteElementAmongstChildren());
           }
           break;
 
index 8b27eb4eb80cc87e8828ad2d8a3758a83029eb3d..9bf0af39cc425283c4005f5126d957d040dacc30 100644 (file)
@@ -266,19 +266,23 @@ namespace parallel
                   for (unsigned int child_index = 0;
                        child_index < cell->n_children();
                        ++child_index)
-                    fe_indices_children.insert(
-                      cell->child(child_index)->future_fe_index());
+                    {
+                      const auto child = cell->child(child_index);
+                      Assert(child->active() && child->coarsen_flag_set(),
+                             typename dealii::Triangulation<
+                               dim>::ExcInconsistentCoarseningFlags());
+
+                      fe_indices_children.insert(child->future_fe_index());
+                    }
+                  Assert(!fe_indices_children.empty(), ExcInternalError());
 
                   fe_index =
                     dof_handler->get_fe_collection().find_dominated_fe_extended(
                       fe_indices_children, /*codim=*/0);
 
-                  Assert(
-                    fe_index != numbers::invalid_unsigned_int,
-                    ExcMessage(
-                      "No FiniteElement has been found in your FECollection "
-                      "that is dominated by all children of a cell you are "
-                      "trying to coarsen!"));
+                  Assert(fe_index != numbers::invalid_unsigned_int,
+                         typename dealii::hp::FECollection<
+                           dim>::ExcNoDominatedFiniteElementAmongstChildren());
 
                   break;
                 }
index e84bab152da3ea82b5c745e21cf2781d51e361e1..1b67554d37b445545b6f2c90d61934ca709bbb5f 100644 (file)
@@ -18,6 +18,7 @@
 #include <deal.II/base/std_cxx14/memory.h>
 #include <deal.II/base/thread_management.h>
 
+#include <deal.II/distributed/cell_data_transfer.templates.h>
 #include <deal.II/distributed/shared_tria.h>
 #include <deal.II/distributed/tria.h>
 
@@ -36,6 +37,8 @@
 #include <deal.II/hp/dof_handler.h>
 #include <deal.II/hp/dof_level.h>
 
+#include <deal.II/numerics/coarsening_strategies.h>
+
 #include <boost/serialization/array.hpp>
 
 #include <algorithm>
@@ -1120,6 +1123,7 @@ namespace internal
                     // particular cell has a parent at all.
                     Assert(cell->level() > 0, ExcInternalError());
                     const auto &parent = cell->parent();
+
                     // Check if the active_fe_index for the current cell has
                     // been determined already.
                     if (fe_transfer->coarsened_cells_fe_index.find(parent) ==
@@ -1134,11 +1138,14 @@ namespace internal
                              child_index < parent->n_children();
                              ++child_index)
                           {
-                            Assert(parent->child(child_index)->active(),
-                                   ExcInternalError());
+                            const auto sibling = parent->child(child_index);
+                            Assert(sibling->active() &&
+                                     sibling->coarsen_flag_set(),
+                                   typename dealii::Triangulation<
+                                     dim>::ExcInconsistentCoarseningFlags());
 
                             fe_indices_children.insert(
-                              parent->child(child_index)->future_fe_index());
+                              sibling->future_fe_index());
                           }
                         Assert(!fe_indices_children.empty(),
                                ExcInternalError());
@@ -1147,12 +1154,9 @@ namespace internal
                           dof_handler.fe_collection.find_dominated_fe_extended(
                             fe_indices_children, /*codim=*/0);
 
-                        Assert(
-                          fe_index != numbers::invalid_unsigned_int,
-                          ExcMessage(
-                            "No FiniteElement has been found in your FECollection "
-                            "that is dominated by all children of a cell you are "
-                            "trying to coarsen!"));
+                        Assert(fe_index != numbers::invalid_unsigned_int,
+                               typename dealii::hp::FECollection<dim>::
+                                 ExcNoDominatedFiniteElementAmongstChildren());
 
                         fe_transfer->coarsened_cells_fe_index.insert(
                           {parent, fe_index});
@@ -1184,22 +1188,22 @@ namespace internal
           const auto &fe_transfer = dof_handler.active_fe_index_transfer;
 
           // Set active_fe_indices on persisting cells.
-          for (const auto &pair : fe_transfer->persisting_cells_fe_index)
+          for (const auto &persist : fe_transfer->persisting_cells_fe_index)
             {
-              const auto &cell = pair.first;
+              const auto &cell = persist.first;
 
               if (cell->is_locally_owned())
                 {
                   Assert(cell->active(), ExcInternalError());
-                  cell->set_active_fe_index(pair.second);
+                  cell->set_active_fe_index(persist.second);
                 }
             }
 
           // Distribute active_fe_indices from all refined cells on their
           // respective children.
-          for (const auto &pair : fe_transfer->refined_cells_fe_index)
+          for (const auto &refine : fe_transfer->refined_cells_fe_index)
             {
-              const auto &parent = pair.first;
+              const auto &parent = refine.first;
 
               for (unsigned int child_index = 0;
                    child_index < parent->n_children();
@@ -1208,18 +1212,18 @@ namespace internal
                   const auto &child = parent->child(child_index);
                   Assert(child->is_locally_owned() && child->active(),
                          ExcInternalError());
-                  child->set_active_fe_index(pair.second);
+                  child->set_active_fe_index(refine.second);
                 }
             }
 
           // Set active_fe_indices on coarsened cells that have been determined
           // before the actual coarsening happened.
-          for (const auto &pair : fe_transfer->coarsened_cells_fe_index)
+          for (const auto &coarsen : fe_transfer->coarsened_cells_fe_index)
             {
-              const auto &cell = pair.first;
+              const auto &cell = coarsen.first;
               Assert(cell->is_locally_owned() && cell->active(),
                      ExcInternalError());
-              cell->set_active_fe_index(pair.second);
+              cell->set_active_fe_index(coarsen.second);
             }
         }
       };
@@ -2085,10 +2089,7 @@ namespace hp
             CellDataTransfer<dim, spacedim, std::vector<unsigned int>>>(
           *distributed_tria,
           /*transfer_variable_size_data=*/false,
-          &parallel::distributed::CellDataTransfer<
-            dim,
-            spacedim,
-            std::vector<unsigned int>>::CoarseningStrategies::check_equality);
+          &CoarseningStrategies::check_equality<unsigned int>);
 
         active_fe_index_transfer->cell_data_transfer
           ->prepare_for_coarsening_and_refinement(
@@ -2191,10 +2192,7 @@ namespace hp
             CellDataTransfer<dim, spacedim, std::vector<unsigned int>>>(
           *distributed_tria,
           /*transfer_variable_size_data=*/false,
-          &parallel::distributed::CellDataTransfer<
-            dim,
-            spacedim,
-            std::vector<unsigned int>>::CoarseningStrategies::check_equality);
+          &CoarseningStrategies::check_equality<unsigned int>);
 
         // If we work on a p::d::Triangulation, we have to transfer all
         // active fe indices since ownership of cells may change.
@@ -2267,10 +2265,7 @@ namespace hp
             CellDataTransfer<dim, spacedim, std::vector<unsigned int>>>(
           *distributed_tria,
           /*transfer_variable_size_data=*/false,
-          &parallel::distributed::CellDataTransfer<
-            dim,
-            spacedim,
-            std::vector<unsigned int>>::CoarseningStrategies::check_equality);
+          &CoarseningStrategies::check_equality<unsigned int>);
 
         // Unpack active_fe_indices.
         active_fe_index_transfer->active_fe_indices.resize(
index 792647f17cc51b0d81f6c6cb8cd5848c1eb9190b..3456035e602c40f445981cbdf9edd8de97c1df33 100644 (file)
@@ -356,17 +356,6 @@ SolutionTransfer<dim, VectorType, DoFHandlerType>::
       // CASE 2: cell is inactive but will become active
       else if (cell->has_children() && cell->child(0)->coarsen_flag_set())
         {
-          // note that if one child has the
-          // coarsen flag, then all should
-          // have if Tria::prepare_* has
-          // worked correctly
-          for (unsigned int i = 1; i < cell->n_children(); ++i)
-            Assert(cell->child(i)->coarsen_flag_set(),
-                   ExcMessage(
-                     "It looks like you didn't call "
-                     "Triangulation::prepare_coarsening_and_refinement before "
-                     "calling the current function. This can't work."));
-
           // we will need to interpolate from the children of this cell
           // to the current one. in the hp context, this also means
           // we need to figure out which finite element space to interpolate
@@ -376,18 +365,23 @@ SolutionTransfer<dim, VectorType, DoFHandlerType>::
           std::set<unsigned int> fe_indices_children;
           for (unsigned int child_index = 0; child_index < cell->n_children();
                ++child_index)
-            fe_indices_children.insert(
-              cell->child(child_index)->active_fe_index());
+            {
+              const auto child = cell->child(child_index);
+              Assert(child->active() && child->coarsen_flag_set(),
+                     typename dealii::Triangulation<
+                       dim>::ExcInconsistentCoarseningFlags());
+
+              fe_indices_children.insert(child->active_fe_index());
+            }
+          Assert(!fe_indices_children.empty(), ExcInternalError());
 
           const unsigned int target_fe_index =
             dof_handler->get_fe_collection().find_dominated_fe_extended(
               fe_indices_children, /*codim=*/0);
 
           Assert(target_fe_index != numbers::invalid_unsigned_int,
-                 ExcMessage(
-                   "No FiniteElement has been found in your FECollection "
-                   "that is dominated by all children of a cell you are "
-                   "trying to coarsen!"));
+                 typename dealii::hp::FECollection<
+                   dim>::ExcNoDominatedFiniteElementAmongstChildren());
 
           const unsigned int dofs_per_cell =
             dof_handler->get_fe(target_fe_index).dofs_per_cell;
index 8cbb796c3b14aa3c3882011ef8071c15729ada21..ef7dc0cc36c337d9d949850ae225fb316ea9248c 100644 (file)
 #include "../tests.h"
 
 
-template <int dim, int spacedim>
 std::vector<int>
-get_data_of_first_child(const typename parallel::distributed::
-                          Triangulation<dim, spacedim>::cell_iterator &parent,
-                        const std::vector<std::vector<int>> &input_vector)
+get_data_of_first_child(const std::vector<std::vector<int>> &children_values)
 {
-  return input_vector[parent->child(0)->active_cell_index()];
+  return children_values[0];
 }
 
 
@@ -89,10 +86,9 @@ test()
   // ----- transfer -----
   parallel::distributed::
     CellDataTransfer<dim, spacedim, std::vector<std::vector<int>>>
-      cell_data_transfer(
-        tria,
-        /*transfer_variable_size_data=*/true,
-        /*coarsening_strategy=*/&get_data_of_first_child<dim, spacedim>);
+      cell_data_transfer(tria,
+                         /*transfer_variable_size_data=*/true,
+                         /*coarsening_strategy=*/&get_data_of_first_child);
 
   cell_data_transfer.prepare_for_coarsening_and_refinement(cell_data);
   tria.execute_coarsening_and_refinement();

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.