]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Group all variables in p::d::Tria for data attachment into one structure.
authorWolfgang Bangerth <bangerth@colostate.edu>
Fri, 20 Apr 2018 20:17:12 +0000 (14:17 -0600)
committerWolfgang Bangerth <bangerth@colostate.edu>
Fri, 20 Apr 2018 20:17:12 +0000 (14:17 -0600)
include/deal.II/distributed/tria.h
source/distributed/tria.cc

index 7d3d01c9de94c3f4c61baa673589a1477cec2ef3..ccf57a85af242f89a48916fa6fa2149536ade144 100644 (file)
@@ -32,6 +32,7 @@
 #include <utility>
 #include <functional>
 #include <tuple>
+#include <type_traits>
 
 #ifdef DEAL_II_WITH_MPI
 #  include <mpi.h>
@@ -839,37 +840,46 @@ namespace parallel
       typename dealii::internal::p4est::types<dim>::ghost  *parallel_ghost;
 
       /**
-       * number of bytes that get attached to the Triangulation through
-       * register_data_attach() for example SolutionTransfer.
+       * A structure that stores information about the data that has been, or
+       * will be, attached to cells via the register_data_attach() function
+       * and later retrieved via notify_ready_to_unpack().
        */
-      unsigned int attached_data_size;
+      struct CellAttachedData
+      {
+        /**
+         * number of bytes that get attached to the Triangulation through
+         * register_data_attach() for example SolutionTransfer.
+         */
+        unsigned int attached_data_size;
 
-      /**
-       * number of functions that get attached to the Triangulation through
-       * register_data_attach() for example SolutionTransfer.
-       */
-      unsigned int n_attached_datas;
+        /**
+         * number of functions that get attached to the Triangulation through
+         * register_data_attach() for example SolutionTransfer.
+         */
+        unsigned int n_attached_datas;
 
-      /**
-       * number of functions that need to unpack their data after a call from
-       * load()
-       */
-      unsigned int n_attached_deserialize;
+        /**
+         * number of functions that need to unpack their data after a call from
+         * load()
+         */
+        unsigned int n_attached_deserialize;
 
-      typedef  std::function<
-      void(typename Triangulation<dim,spacedim>::cell_iterator, CellStatus, void *)
-      > pack_callback_t;
+        using pack_callback_t = std::function<void (typename Triangulation<dim,spacedim>::cell_iterator,
+                                                    CellStatus,
+                                                    void *)>;
 
-      typedef std::pair<unsigned int, pack_callback_t> callback_pair_t;
+        using callback_pair_t = std::pair<unsigned int, pack_callback_t>;
 
-      typedef std::list<callback_pair_t> callback_list_t;
+        using callback_list_t = std::list<callback_pair_t>;
 
-      /**
-       * List of callback functions registered by register_data_attach() that
-       * are going to be called for packing data.
-       */
-      callback_list_t attached_data_pack_callbacks;
+        /**
+         * List of callback functions registered by register_data_attach() that
+         * are going to be called for packing data.
+         */
+        callback_list_t attached_data_pack_callbacks;
+      };
 
+      CellAttachedData cell_attached_data;
 
       /**
        * Two arrays that store which p4est tree corresponds to which coarse
index ec773e0f52e0215ea5bc2d08e1f2512871144fb5..d364fa2bd8f65a50406ef93b81cf81d9159fab7d 100644 (file)
@@ -1315,9 +1315,7 @@ namespace parallel
       triangulation_has_content (false),
       connectivity (nullptr),
       parallel_forest (nullptr),
-      attached_data_size(0),
-      n_attached_datas(0),
-      n_attached_deserialize(0)
+      cell_attached_data ({0, 0, 0, {}})
     {
       parallel_ghost = nullptr;
     }
@@ -1795,15 +1793,15 @@ namespace parallel
     Triangulation<dim,spacedim>::
     save(const char *filename) const
     {
-      Assert(n_attached_deserialize==0,
+      Assert(cell_attached_data.n_attached_deserialize==0,
              ExcMessage ("not all SolutionTransfer's got deserialized after the last load()"));
 
       // signal that serialization is going to happen
       this->signals.pre_distributed_save();
 
       int real_data_size = 0;
-      if (attached_data_size>0)
-        real_data_size = attached_data_size+sizeof(CellStatus);
+      if (cell_attached_data.attached_data_size>0)
+        real_data_size = cell_attached_data.attached_data_size+sizeof(CellStatus);
 
       Assert(this->n_cells()>0, ExcMessage("Can not save() an empty Triangulation."));
 
@@ -1815,18 +1813,19 @@ namespace parallel
             << 2 << " "
             << Utilities::MPI::n_mpi_processes (this->mpi_communicator) << " "
             << real_data_size << " "
-            << attached_data_pack_callbacks.size() << " "
+            << cell_attached_data.attached_data_pack_callbacks.size() << " "
             << this->n_cells(0)
             << std::endl;
         }
 
-      if (attached_data_size>0)
+      if (cell_attached_data.attached_data_size>0)
         {
           const_cast<dealii::parallel::distributed::Triangulation<dim, spacedim>*>(this)
           ->attach_mesh_data();
         }
 
-      dealii::internal::p4est::functions<dim>::save(filename, parallel_forest, attached_data_size>0);
+      dealii::internal::p4est::functions<dim>::save(filename, parallel_forest,
+                                                    cell_attached_data.attached_data_size>0);
 
       // clear all of the callback data, as explained in the documentation of
       // register_data_attach()
@@ -1834,9 +1833,9 @@ namespace parallel
         dealii::parallel::distributed::Triangulation<dim, spacedim> *tria
           = const_cast<dealii::parallel::distributed::Triangulation<dim, spacedim>*>(this);
 
-        tria->n_attached_datas = 0;
-        tria->attached_data_size = 0;
-        tria->attached_data_pack_callbacks.clear();
+        tria->cell_attached_data.n_attached_datas = 0;
+        tria->cell_attached_data.attached_data_size = 0;
+        tria->cell_attached_data.attached_data_pack_callbacks.clear();
       }
 
       // and release the data
@@ -1889,9 +1888,9 @@ namespace parallel
 
       // clear all of the callback data, as explained in the documentation of
       // register_data_attach()
-      attached_data_size = 0;
-      n_attached_datas = 0;
-      n_attached_deserialize = attached_count;
+      cell_attached_data.attached_data_size = 0;
+      cell_attached_data.n_attached_datas = 0;
+      cell_attached_data.n_attached_deserialize = attached_count;
 
 #if DEAL_II_P4EST_VERSION_GTE(0,3,4,3)
       parallel_forest = dealii::internal::p4est::functions<dim>::load_ext (
@@ -3245,13 +3244,13 @@ namespace parallel
                                                    void *)> &pack_callback)
     {
       Assert(size>0, ExcMessage("register_data_attach(), size==0"));
-      Assert(attached_data_pack_callbacks.size()==n_attached_datas,
+      Assert(cell_attached_data.attached_data_pack_callbacks.size()==cell_attached_data.n_attached_datas,
              ExcMessage("register_data_attach(), not all data has been unpacked last time?"));
 
-      const unsigned int offset = attached_data_size+sizeof(CellStatus);
-      ++n_attached_datas;
-      attached_data_size += size;
-      attached_data_pack_callbacks.emplace_back(offset, pack_callback);
+      const unsigned int offset = cell_attached_data.attached_data_size+sizeof(CellStatus);
+      ++cell_attached_data.n_attached_datas;
+      cell_attached_data.attached_data_size += size;
+      cell_attached_data.attached_data_pack_callbacks.emplace_back(offset, pack_callback);
 
       return offset;
     }
@@ -3268,9 +3267,9 @@ namespace parallel
     {
       Assert (handle >= sizeof(CellStatus),
               ExcMessage ("Invalid offset."));
-      Assert (handle < sizeof(CellStatus)+attached_data_size,
+      Assert (handle < sizeof(CellStatus)+cell_attached_data.attached_data_size,
               ExcMessage ("Invalid offset."));
-      Assert (n_attached_datas > 0,
+      Assert (cell_attached_data.n_attached_datas > 0,
               ExcMessage ("The notify_ready_to_unpack() has already been called "
                           "once for each registered callback."));
 
@@ -3302,11 +3301,11 @@ namespace parallel
                                                      unpack_callback);
         }
 
-      --n_attached_datas;
-      if (n_attached_deserialize > 0)
+      --cell_attached_data.n_attached_datas;
+      if (cell_attached_data.n_attached_deserialize > 0)
         {
-          --n_attached_deserialize;
-          attached_data_pack_callbacks.pop_front();
+          --cell_attached_data.n_attached_deserialize;
+          cell_attached_data.attached_data_pack_callbacks.pop_front();
         }
 
       // important: only remove data if we are not in the deserialization
@@ -3315,11 +3314,13 @@ namespace parallel
       // would destroy the saved data before the second SolutionTransfer can
       // get it. This created a bug that is documented in
       // tests/mpi/p4est_save_03 with more than one SolutionTransfer.
-      if (n_attached_datas == 0 && n_attached_deserialize == 0)
+      if (cell_attached_data.n_attached_datas == 0
+          &&
+          cell_attached_data.n_attached_deserialize == 0)
         {
           // everybody got their data, time for cleanup!
-          attached_data_size = 0;
-          attached_data_pack_callbacks.clear();
+          cell_attached_data.attached_data_size = 0;
+          cell_attached_data.attached_data_pack_callbacks.clear();
 
           // and release the data
           void *userptr = parallel_forest->user_pointer;
@@ -3585,8 +3586,8 @@ namespace parallel
         + MemoryConsumption::memory_consumption(triangulation_has_content)
         + MemoryConsumption::memory_consumption(connectivity)
         + MemoryConsumption::memory_consumption(parallel_forest)
-        + MemoryConsumption::memory_consumption(attached_data_size)
-        + MemoryConsumption::memory_consumption(n_attached_datas)
+        + MemoryConsumption::memory_consumption(cell_attached_data.attached_data_size)
+        + MemoryConsumption::memory_consumption(cell_attached_data.n_attached_datas)
 //      + MemoryConsumption::memory_consumption(attached_data_pack_callbacks) //TODO[TH]: how?
         + MemoryConsumption::memory_consumption(coarse_cell_to_p4est_tree_permutation)
         + MemoryConsumption::memory_consumption(p4est_tree_to_coarse_cell_permutation)
@@ -3637,8 +3638,8 @@ namespace parallel
         {
           coarse_cell_to_p4est_tree_permutation = other_tria_x->coarse_cell_to_p4est_tree_permutation;
           p4est_tree_to_coarse_cell_permutation = other_tria_x->p4est_tree_to_coarse_cell_permutation;
-          attached_data_size = other_tria_x->attached_data_size;
-          n_attached_datas   = other_tria_x->n_attached_datas;
+          cell_attached_data.attached_data_size = other_tria_x->cell_attached_data.attached_data_size;
+          cell_attached_data.n_attached_datas   = other_tria_x->cell_attached_data.n_attached_datas;
 
           settings           = other_tria_x->settings;
         }
@@ -3673,9 +3674,9 @@ namespace parallel
     {
       // determine size of memory in bytes to attach to each cell. This needs
       // to be constant because of p4est.
-      if (attached_data_size==0)
+      if (cell_attached_data.attached_data_size==0)
         {
-          Assert(n_attached_datas==0, ExcInternalError());
+          Assert(cell_attached_data.n_attached_datas==0, ExcInternalError());
 
           //nothing to do
           return;
@@ -3684,7 +3685,7 @@ namespace parallel
       // realloc user_data in p4est
       void *userptr = parallel_forest->user_pointer;
       dealii::internal::p4est::functions<dim>::reset_data (parallel_forest,
-                                                           attached_data_size+sizeof(CellStatus),
+                                                           cell_attached_data.attached_data_size+sizeof(CellStatus),
                                                            nullptr, nullptr);
       parallel_forest->user_pointer = userptr;
 
@@ -3712,7 +3713,7 @@ namespace parallel
           attach_mesh_data_recursively<dim,spacedim>(*tree,
                                                      cell,
                                                      p4est_coarse_cell,
-                                                     attached_data_pack_callbacks);
+                                                     cell_attached_data.attached_data_pack_callbacks);
         }
     }
 

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.