]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Triangulation now manages offsets. Data attaching classes now take handles from the... 6472/head
authorMarc Fehling <marc.fehling@gmx.net>
Wed, 2 May 2018 19:36:01 +0000 (13:36 -0600)
committerMarc Fehling <marc.fehling@gmx.net>
Tue, 8 May 2018 21:47:35 +0000 (15:47 -0600)
include/deal.II/base/quadrature_point_data.h
include/deal.II/distributed/solution_transfer.h
include/deal.II/distributed/tria.h
include/deal.II/particles/particle_handler.h
source/distributed/solution_transfer.cc
source/distributed/tria.cc
source/particles/particle_handler.cc
tests/mpi/attach_data_01.cc
tests/mpi/attach_data_01.mpirun=1.output
tests/mpi/attach_data_01.mpirun=3.output

index a9915939aa9e69adcff7b3ad197db0e562dfef76..f5894f4291d084670b501f463bbb40968b9bc650 100644 (file)
@@ -438,10 +438,10 @@ namespace parallel
       FullMatrix<double> matrix_quadrature;
 
       /**
-       * The offset that the parallel::distributed::Triangulation has assigned
-       * to this object starting at which we are allowed to read.
+       * The handle that the parallel::distributed::Triangulation has assigned
+       * to this object while registering the pack_callback function.
        */
-      unsigned int offset;
+      unsigned int handle;
 
       /**
        * A pointer to the CellDataStorage class whose data will be transferred.
@@ -674,7 +674,7 @@ namespace parallel
       n_q_points(rhs_quadrature.size()),
       project_to_fe_matrix(projection_fe->dofs_per_cell,n_q_points),
       project_to_qp_matrix(n_q_points,projection_fe->dofs_per_cell),
-      offset(0),
+      handle(numbers::invalid_unsigned_int),
       data_storage(nullptr),
       triangulation(nullptr)
     {
@@ -732,7 +732,7 @@ namespace parallel
       matrix_quadrature.reinit(n_q_points,number_of_values);
       data_size_in_bytes = sizeof(double) * dofs_per_cell * number_of_values;
 
-      offset = triangulation->register_data_attach(data_size_in_bytes,
+      handle = triangulation->register_data_attach(data_size_in_bytes,
                                                    std::bind(&ContinuousQuadratureDataTransfer<dim,DataType>::pack_function,
                                                              this,
                                                              std::placeholders::_1,
@@ -745,7 +745,7 @@ namespace parallel
     template <int dim, typename DataType>
     void ContinuousQuadratureDataTransfer<dim,DataType>::interpolate ()
     {
-      triangulation->notify_ready_to_unpack(offset,
+      triangulation->notify_ready_to_unpack(handle,
                                             std::bind(&ContinuousQuadratureDataTransfer<dim,DataType>::unpack_function,
                                                       this,
                                                       std::placeholders::_1,
index 8de9df5d3b4d5ce8b7f6a670e3f3ecdd28109e15..9fa8aab21a3fe7ba2a34ed24730dd87a699a51ab 100644 (file)
@@ -233,10 +233,10 @@ namespace parallel
       std::vector<const VectorType *> input_vectors;
 
       /**
-       * The offset that the Triangulation has assigned to this object
-       * starting at which we are allowed to write.
+       * The handle that the Triangulation has assigned to this object
+       * with which we can access our memory offset and our pack function.
        */
-      unsigned int offset;
+      unsigned int handle;
 
       /**
        * A callback function used to pack the data on the current mesh into
index d742a651f2ad8b47ffdc417a17b1e8de2f99071d..8e6ee330bdd97112e06a0d5e89fc8f50b3c4e95d 100644 (file)
@@ -872,12 +872,12 @@ namespace parallel
 
         /**
          * List of callback functions registered by register_data_attach() that
-         * are going to be called for packing data. The key of this map
-         * variable is the offset at which each callback is allowed to
-         * write into the per-cell buffer (counted in bytes) whereas the
-         * value of the map is the callback function object.
+         * are going to be called for packing data. The callback functions objects
+         * will be stored together with their offset at which each callback is
+         * allowed to write into the per-cell buffer (counted in bytes). These
+         * pairs will be stored in the order on how they have been registered.
          */
-        std::map<unsigned int, pack_callback_t> pack_callbacks;
+        std::vector< std::pair<unsigned int, pack_callback_t> > pack_callbacks;
       };
 
       CellAttachedData cell_attached_data;
index 4a30326bcdc341bb2083dddaac84051d08dd8ec8..b8d01b6f142d23ebdcaf0d6294287286d08c51f4 100644 (file)
@@ -426,9 +426,10 @@ namespace Particles
     /**
      * This variable is set by the register_store_callback_function()
      * function and used by the register_load_callback_function() function
-     * to check where the particle data was stored.
+     * to check where the particle data was registered in the corresponding
+     * triangulation object.
      */
-    unsigned int data_offset;
+    unsigned int handle;
 
 #ifdef DEAL_II_WITH_MPI
     /**
index da4cbb55e2520bf85675447d6a76421822e329dd..b594bc599aeb71841b73f532fe697fdbdafe6650 100644 (file)
@@ -47,7 +47,7 @@ namespace parallel
     SolutionTransfer<dim, VectorType, DoFHandlerType>::SolutionTransfer (const DoFHandlerType &dof)
       :
       dof_handler(&dof, typeid(*this).name()),
-      offset (numbers::invalid_unsigned_int)
+      handle (numbers::invalid_unsigned_int)
     {
       Assert ((dynamic_cast<const parallel::distributed::Triangulation<dim,DoFHandlerType::space_dimension>*>
                (&dof_handler->get_triangulation()) != nullptr),
@@ -80,7 +80,7 @@ namespace parallel
             (&dof_handler->get_triangulation())));
       Assert (tria != nullptr, ExcInternalError());
 
-      offset
+      handle
         = tria->register_data_attach(size,
                                      std::bind(&SolutionTransfer<dim, VectorType,
                                                DoFHandlerType>::pack_callback,
@@ -161,7 +161,7 @@ namespace parallel
             (&dof_handler->get_triangulation())));
       Assert (tria != nullptr, ExcInternalError());
 
-      tria->notify_ready_to_unpack(offset,
+      tria->notify_ready_to_unpack(handle,
                                    std::bind(&SolutionTransfer<dim, VectorType,
                                              DoFHandlerType>::unpack_callback,
                                              this,
index db5d21d2b12a354e17b6303df140c046e0da60f8..bc7bc40137c966c1694b6b381f6463027babe4a6 100644 (file)
@@ -552,11 +552,11 @@ namespace
   attach_mesh_data_recursively (const typename internal::p4est::types<dim>::tree &tree,
                                 const typename Triangulation<dim,spacedim>::cell_iterator &dealii_cell,
                                 const typename internal::p4est::types<dim>::quadrant &p4est_cell,
-                                const std::map<unsigned int, typename std::function<
+                                const std::vector< std::pair<unsigned int, typename std::function<
                                 void(typename parallel::distributed::Triangulation<dim,spacedim>::cell_iterator,
                                      typename parallel::distributed::Triangulation<dim,spacedim>::CellStatus,
                                      void *)
-                                > > &pack_callbacks)
+                                > > &pack_callbacks)
   {
     int idx = sc_array_bsearch(const_cast<sc_array_t *>(&tree.quadrants),
                                &p4est_cell,
@@ -803,7 +803,7 @@ namespace
                               const typename Triangulation<dim,spacedim>::cell_iterator &dealii_cell,
                               const typename Triangulation<dim,spacedim>::cell_iterator &parent_cell,
                               const typename internal::p4est::types<dim>::quadrant &p4est_cell,
-                              const unsigned int handle,
+                              const unsigned int offset,
                               const typename std::function<
                               void(typename parallel::distributed::Triangulation<dim,spacedim>::cell_iterator, typename parallel::distributed::Triangulation<dim,spacedim>::CellStatus, void *)
                               > &unpack_callback)
@@ -850,7 +850,7 @@ namespace
                                                       dealii_cell->child(c),
                                                       dealii_cell,
                                                       p4est_child[c],
-                                                      handle,
+                                                      offset,
                                                       unpack_callback);
           }
       }
@@ -863,7 +863,7 @@ namespace
               sc_array_index (const_cast<sc_array_t *>(&tree.quadrants), idx)
             );
 
-        void *ptr = static_cast<char *>(q->p.user_data) + handle;
+        void *ptr = static_cast<char *>(q->p.user_data) + offset;
         typename parallel::distributed::Triangulation<dim,spacedim>::CellStatus
         status = * static_cast<
                  typename parallel::distributed::Triangulation<dim,spacedim>::CellStatus *
@@ -3226,18 +3226,30 @@ namespace parallel
                                                    const CellStatus,
                                                    void *)> &pack_callback)
     {
+#ifdef DEBUG
       Assert(size>0, ExcMessage("register_data_attach(), size==0"));
-      Assert(cell_attached_data.pack_callbacks.size()==cell_attached_data.n_attached_datas,
+
+      unsigned int remaining_callback_counter = 0;
+      for (auto it : cell_attached_data.pack_callbacks)
+        if (std::get<0>(it) != numbers::invalid_unsigned_int)
+          ++remaining_callback_counter;
+
+      Assert(remaining_callback_counter==cell_attached_data.n_attached_datas,
              ExcMessage("register_data_attach(), not all data has been unpacked last time?"));
+#endif
+
+      // add new pair with offset and pack_callback
+      cell_attached_data.pack_callbacks.push_back(
+        std::make_pair(cell_attached_data.cumulative_fixed_data_size + sizeof(CellStatus),
+                       pack_callback)
+      );
 
-      const unsigned int current_buffer_positition = cell_attached_data.cumulative_fixed_data_size+sizeof(CellStatus);
+      // increase counters
       ++cell_attached_data.n_attached_datas;
       cell_attached_data.cumulative_fixed_data_size += size;
-      cell_attached_data.pack_callbacks[current_buffer_positition] = pack_callback;
 
-      // use the offset as the handle that callers receive and later hand back
-      // to notify_ready_to_unpack()
-      return current_buffer_positition;
+      // return the position of the callback in the register vector as a handle
+      return cell_attached_data.pack_callbacks.size() - 1;
     }
 
 
@@ -3250,14 +3262,16 @@ namespace parallel
                                                       const CellStatus,
                                                       const void *)> &unpack_callback)
     {
-      Assert (handle >= sizeof(CellStatus),
-              ExcMessage ("Invalid offset."));
-      Assert (handle < cell_attached_data.cumulative_fixed_data_size + sizeof(CellStatus),
-              ExcMessage ("Invalid offset."));
+      Assert (handle < cell_attached_data.pack_callbacks.size(),
+              ExcMessage ("Invalid handle."));
       Assert (cell_attached_data.n_attached_datas > 0,
               ExcMessage ("The notify_ready_to_unpack() has already been called "
                           "once for each registered callback."));
 
+      const unsigned int offset = std::get<0> (cell_attached_data.pack_callbacks[handle]);
+      // check if unpack function has been previously called
+      Assert (offset!=dealii::numbers::invalid_unsigned_int, ExcInternalError());
+
       // Recurse over p4est and hand the caller the data back
       for (typename Triangulation<dim, spacedim>::cell_iterator
            cell = this->begin (0);
@@ -3282,24 +3296,21 @@ namespace parallel
                                                      cell,
                                                      cell,
                                                      p4est_coarse_cell,
-                                                     handle,
+                                                     offset,
                                                      unpack_callback);
         }
 
       --cell_attached_data.n_attached_datas;
       if (cell_attached_data.n_attached_deserialize > 0)
-        {
-          --cell_attached_data.n_attached_deserialize;
-
-          // Remove the entry that corresponds to the handle that was
-          // passed to the current function. Double check that such an
-          // entry really existed.
-          const unsigned int n_elements_removed
-            = cell_attached_data.pack_callbacks.erase (handle);
-          (void)n_elements_removed;
-          Assert (n_elements_removed == 1,
-                  ExcInternalError());
-        }
+        --cell_attached_data.n_attached_deserialize;
+
+      // deregister corresponding pack_callback function
+      // first reset with invalid entries to preserve ambiguity of handles,
+      // then free memory when all were unpacked (see below)
+      std::get<0> (cell_attached_data.pack_callbacks[handle]) // offset value
+        = dealii::numbers::invalid_unsigned_int;
+      std::get<1> (cell_attached_data.pack_callbacks[handle]) // pack_callback function
+        = typename CellAttachedData::pack_callback_t ();
 
       // important: only remove data if we are not in the deserialization
       // process. There, each SolutionTransfer registers and unpacks before
index 6d365878132e27899fd72d6b05fdf2d36380faa7..70f5d10a19d826ebe6e43d8cb1fca95776e855df 100644 (file)
@@ -40,7 +40,7 @@ namespace Particles
     size_callback(),
     store_callback(),
     load_callback(),
-    data_offset(numbers::invalid_unsigned_int)
+    handle(numbers::invalid_unsigned_int)
   {}
 
 
@@ -61,7 +61,7 @@ namespace Particles
     size_callback(),
     store_callback(),
     load_callback(),
-    data_offset(numbers::invalid_unsigned_int)
+    handle(numbers::invalid_unsigned_int)
   {}
 
 
@@ -922,7 +922,7 @@ namespace Particles
                                                     :
                                                     std::pow(2,dim));
 
-        data_offset = non_const_triangulation->register_data_attach(transfer_size_per_cell,callback_function);
+        handle = non_const_triangulation->register_data_attach(transfer_size_per_cell,callback_function);
       }
   }
 
@@ -968,11 +968,11 @@ namespace Particles
         // space for the data in case a cell is coarsened
         const std::size_t transfer_size_per_cell = sizeof (unsigned int) +
                                                    (size_per_particle * global_max_particles_per_cell);
-        data_offset = non_const_triangulation->register_data_attach(transfer_size_per_cell,callback_function);
+        handle = non_const_triangulation->register_data_attach(transfer_size_per_cell,callback_function);
       }
 
     // Check if something was stored and load it
-    if (data_offset != numbers::invalid_unsigned_int)
+    if (handle != numbers::invalid_unsigned_int)
       {
         const std::function<void(const typename Triangulation<dim,spacedim>::cell_iterator &,
                                  const typename Triangulation<dim,spacedim>::CellStatus,
@@ -983,11 +983,11 @@ namespace Particles
                       std::placeholders::_2,
                       std::placeholders::_3);
 
-        non_const_triangulation->notify_ready_to_unpack(data_offset,callback_function);
+        non_const_triangulation->notify_ready_to_unpack(handle,callback_function);
 
-        // Reset offset and update global number of particles. The number
+        // Reset handle and update global number of particles. The number
         // can change because of discarded or newly generated particles
-        data_offset = numbers::invalid_unsigned_int;
+        handle = numbers::invalid_unsigned_int;
         update_cached_numbers();
       }
   }
index 308bce201b9bddd648a0507b03b8d0d8fec7b0da..b32ed451fae17aa6c36938f8404fe6f670ed2fae 100644 (file)
@@ -125,9 +125,9 @@ void test()
             }
         }
 
-      unsigned int offset = tr.register_data_attach(sizeof(int), pack_function<dim>);
+      unsigned int handle = tr.register_data_attach(sizeof(int), pack_function<dim>);
 
-      deallog << "offset=" << offset << std::endl;
+      deallog << "handle=" << handle << std::endl;
       tr.execute_coarsening_and_refinement();
 
       deallog << "cells after: " << tr.n_global_active_cells() << std::endl;
@@ -141,7 +141,7 @@ void test()
       deallog << "myid=" << myid << " cellid=" << cell->id() << std::endl;
       }*/
 
-      tr.notify_ready_to_unpack(offset, unpack_function<dim>);
+      tr.notify_ready_to_unpack(handle, unpack_function<dim>);
 
       const unsigned int checksum = tr.get_checksum ();
       deallog << "Checksum: "
index c13bf0b12b6e09175020a3880e92647439432378..85bce57b1932cf5a372f03f3af3d116a3a040049 100644 (file)
@@ -16,7 +16,7 @@ DEAL:0::myid=0 cellid=3_1:0 coarsening
 DEAL:0::myid=0 cellid=3_1:1 coarsening
 DEAL:0::myid=0 cellid=3_1:2 coarsening
 DEAL:0::myid=0 cellid=3_1:3 coarsening
-DEAL:0::offset=4
+DEAL:0::handle=0
 DEAL:0::packing cell 0_1:0 with data=0 status=REFINE
 DEAL:0::packing cell 0_1:1 with data=1 status=PERSIST
 DEAL:0::packing cell 0_1:2 with data=2 status=PERSIST
index ca1c6de2549397e818a7153e20906ee17ac65be7..791e44f0f4e46331b60dd3f19ac8f5db9c369374 100644 (file)
@@ -4,7 +4,7 @@ DEAL:0::myid=0 cellid=0_1:0 refining
 DEAL:0::myid=0 cellid=0_1:1
 DEAL:0::myid=0 cellid=0_1:2
 DEAL:0::myid=0 cellid=0_1:3
-DEAL:0::offset=4
+DEAL:0::handle=0
 DEAL:0::packing cell 0_1:0 with data=0 status=REFINE
 DEAL:0::packing cell 0_1:1 with data=1 status=PERSIST
 DEAL:0::packing cell 0_1:2 with data=2 status=PERSIST
@@ -24,7 +24,7 @@ DEAL:1::myid=1 cellid=2_1:0
 DEAL:1::myid=1 cellid=2_1:1
 DEAL:1::myid=1 cellid=2_1:2
 DEAL:1::myid=1 cellid=2_1:3
-DEAL:1::offset=4
+DEAL:1::handle=0
 DEAL:1::packing cell 1_1:0 with data=0 status=PERSIST
 DEAL:1::packing cell 1_1:1 with data=1 status=PERSIST
 DEAL:1::packing cell 1_1:2 with data=2 status=PERSIST
@@ -49,7 +49,7 @@ DEAL:2::myid=2 cellid=3_1:0 coarsening
 DEAL:2::myid=2 cellid=3_1:1 coarsening
 DEAL:2::myid=2 cellid=3_1:2 coarsening
 DEAL:2::myid=2 cellid=3_1:3 coarsening
-DEAL:2::offset=4
+DEAL:2::handle=0
 DEAL:2::packing cell 3_0: with data=0 status=COARSEN
 DEAL:2::cells after: 16
 DEAL:2::unpacking cell 2_1:0 with data=4 status=PERSIST

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.