]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Triangulation: pack/unpack fixed/variable data 17273/head
authorPeter Munch <peterrmuench@gmail.com>
Wed, 17 Jul 2024 07:43:39 +0000 (09:43 +0200)
committerPeter Munch <peterrmuench@gmail.com>
Wed, 17 Jul 2024 15:48:32 +0000 (17:48 +0200)
include/deal.II/distributed/fully_distributed_tria.h
include/deal.II/distributed/tria.h
include/deal.II/grid/tria.h
source/grid/tria.cc
source/particles/particle_handler.cc

index a54abc6d04b0225c01c2e069b3c03c5aadbc6f61..5480ee3b2b9d2262615a49aaaa929d50910f7be9 100644 (file)
@@ -294,8 +294,8 @@ namespace parallel
        * The stored vector will have a size equal to the number of locally owned
        * active cells and will be ordered by the occurrence of those cells.
        */
-      virtual void
-      update_cell_relations() override;
+      void
+      update_cell_relations();
 
       virtual void
       update_number_cache() override;
index f9640b3b98325f25d4933a1379f92061757c142c..414c745b63055f10bce05465e397c035aee55514 100644 (file)
@@ -780,8 +780,8 @@ namespace parallel
        * the refinement process. With this information, we can prepare all
        * buffers for data transfer accordingly.
        */
-      virtual void
-      update_cell_relations() override;
+      void
+      update_cell_relations();
 
       /**
        * Two arrays that store which p4est tree corresponds to which coarse
@@ -964,8 +964,8 @@ namespace parallel
        * This function is not implemented, but needs to be present for the
        * compiler.
        */
-      virtual void
-      update_cell_relations() override;
+      void
+      update_cell_relations();
 
       /**
        * Dummy arrays. This class isn't usable but the compiler wants to see
@@ -1104,8 +1104,8 @@ namespace parallel
        * Dummy replacement to allow for better error messages when compiling
        * this class.
        */
-      virtual void
-      update_cell_relations() override
+      void
+      update_cell_relations()
       {}
     };
   } // namespace distributed
index 86f3f2946d27426396646c083d8f1dc55ddc99a8..9f187889d864242ce44287557b4479853148f94b 100644 (file)
@@ -367,7 +367,7 @@ namespace internal
    * will be, attached to cells via the register_data_attach() function
    * and later retrieved via notify_ready_to_unpack().
    *
-   * This internalclass is dedicated to the data serialization and transfer
+   * This internal class is dedicated to the data serialization and transfer
    * across repartitioned meshes and to/from the file system.
    *
    * It is designed to store all data buffers intended for serialization.
@@ -379,6 +379,12 @@ namespace internal
   public:
     using cell_iterator = TriaIterator<CellAccessor<dim, spacedim>>;
 
+    /**
+     * Version number stored in the .info file written by
+     * Triangulation::save().
+     */
+    static inline constexpr unsigned int version_number = 5;
+
     /**
      * Auxiliary data structure for assigning a CellStatus to a deal.II cell
      * iterator. For an extensive description of the former, see the
@@ -3881,18 +3887,16 @@ protected:
                      const unsigned int n_attached_deserialize_variable);
 
   /**
-   * A function to record the CellStatus of currently active cells that
-   * are locally owned. This information is mandatory to transfer data
-   * between meshes during adaptation or serialization, e.g., using
-   * parallel::distributed::SolutionTransfer.
+   * A function to record the CellStatus of currently active cells.
+   * This information is mandatory to transfer data between meshes
+   * during adaptation or serialization, e.g., using SolutionTransfer.
    *
    * Relations will be stored in the private member local_cell_relations. For
    * an extensive description of CellStatus, see the documentation for the
    * member function register_data_attach().
    */
-  virtual void
-  update_cell_relations()
-  {}
+  void
+  update_cell_relations();
 
   /**
    * Vector of pairs, each containing a deal.II cell iterator and its
index 22979763d34ccd40a5745e162275600b30d97f7e..fa228791d3bae242b364773aae21b6d5800d338b 100644 (file)
@@ -12701,6 +12701,8 @@ void Triangulation<dim, spacedim>::create_triangulation(
       clear_user_flags();
     }
 
+  this->update_cell_relations();
+
   // inform all listeners that the triangulation has been created
   signals.create();
 }
@@ -13691,20 +13693,76 @@ template <int dim, int spacedim>
 DEAL_II_CXX20_REQUIRES((concepts::is_valid_dim_spacedim<dim, spacedim>))
 void Triangulation<dim, spacedim>::save(const std::string &filename) const
 {
-  // Create boost archive then call alternative version of the save function
-  std::ofstream                 ofs(filename);
+  // Save triangulation information.
+  std::ofstream                 ofs(filename + "_triangulation.data");
   boost::archive::text_oarchive oa(ofs, boost::archive::no_header);
   save(oa, 0);
+
+  // Save attached data.
+  {
+    std::ofstream ifs(filename + ".info");
+    ifs
+      << "version nproc n_attached_fixed_size_objs n_attached_variable_size_objs n_active_cells"
+      << std::endl
+      << internal::CellAttachedDataSerializer<dim, spacedim>::version_number
+      << " " << 1 << " " << this->cell_attached_data.pack_callbacks_fixed.size()
+      << " " << this->cell_attached_data.pack_callbacks_variable.size() << " "
+      << this->n_global_active_cells() << std::endl;
+  }
+
+  this->save_attached_data(0, this->n_global_active_cells(), filename);
 }
 
 template <int dim, int spacedim>
 DEAL_II_CXX20_REQUIRES((concepts::is_valid_dim_spacedim<dim, spacedim>))
 void Triangulation<dim, spacedim>::load(const std::string &filename)
 {
-  // Create boost archive then call alternative version of the load function
-  std::ifstream                 ifs(filename);
+  // Load triangulation information.
+  std::ifstream                 ifs(filename + "_triangulation.data");
   boost::archive::text_iarchive ia(ifs, boost::archive::no_header);
   load(ia, 0);
+
+  // Load attached data.
+  unsigned int version, numcpus, attached_count_fixed, attached_count_variable,
+    n_global_active_cells;
+  {
+    std::ifstream ifs(std::string(filename) + ".info");
+    AssertThrow(ifs.fail() == false, ExcIO());
+    std::string firstline;
+    getline(ifs, firstline);
+    ifs >> version >> numcpus >> attached_count_fixed >>
+      attached_count_variable >> n_global_active_cells;
+  }
+
+  AssertThrow(numcpus == 1,
+              ExcMessage("Incompatible number of CPUs found in .info file."));
+
+  const auto expected_version =
+    ::dealii::internal::CellAttachedDataSerializer<dim,
+                                                   spacedim>::version_number;
+  AssertThrow(version == expected_version,
+              ExcMessage(
+                "The information saved in the file you are trying "
+                "to read the triangulation from was written with an "
+                "incompatible file format version and cannot be read."));
+  Assert(this->n_global_active_cells() == n_global_active_cells,
+         ExcMessage("The number of cells of the triangulation differs "
+                    "from the number of cells written into the .info file."));
+
+  // Clear all of the callback data, as explained in the documentation of
+  // register_data_attach().
+  this->cell_attached_data.n_attached_data_sets = 0;
+  this->cell_attached_data.n_attached_deserialize =
+    attached_count_fixed + attached_count_variable;
+
+  this->load_attached_data(0,
+                           this->n_global_active_cells(),
+                           this->n_active_cells(),
+                           filename,
+                           attached_count_fixed,
+                           attached_count_variable);
+
+  this->update_cell_relations();
 }
 
 #endif
@@ -15563,6 +15621,26 @@ const typename std::map<
 }
 
 
+template <int dim, int spacedim>
+DEAL_II_CXX20_REQUIRES((concepts::is_valid_dim_spacedim<dim, spacedim>))
+void Triangulation<dim, spacedim>::update_cell_relations()
+{
+  // We only update the cell relations here for serial triangulations.
+  // For other triangulations, this is done at other stages of
+  // mesh creation and mesh refinement.
+  if (dynamic_cast<parallel::DistributedTriangulationBase<dim, spacedim> *>(
+        this))
+    return;
+
+  this->local_cell_relations.clear();
+  this->local_cell_relations.reserve(this->n_active_cells());
+
+  for (const auto &cell : this->active_cell_iterators())
+    this->local_cell_relations.emplace_back(
+      cell, ::dealii::CellStatus::cell_will_persist);
+}
+
+
 
 template <int dim, int spacedim>
 DEAL_II_CXX20_REQUIRES((concepts::is_valid_dim_spacedim<dim, spacedim>))
@@ -15611,6 +15689,9 @@ void Triangulation<dim, spacedim>::execute_coarsening_and_refinement()
 
   update_periodic_face_map();
 
+  if (this->cell_attached_data.n_attached_data_sets == 0)
+    this->update_cell_relations();
+
 #  ifdef DEBUG
 
   // In debug mode, we want to check for some consistency of the
index 1e9d71ef6783adb69595e9796dd51a208c900c94..5e760fbec1dd1714fabc9254b46ccd4839f978e9 100644 (file)
@@ -2194,21 +2194,6 @@ namespace Particles
   void
   ParticleHandler<dim, spacedim>::register_data_attach()
   {
-    parallel::DistributedTriangulationBase<dim, spacedim>
-      *distributed_triangulation =
-        const_cast<parallel::DistributedTriangulationBase<dim, spacedim> *>(
-          dynamic_cast<
-            const parallel::DistributedTriangulationBase<dim, spacedim> *>(
-            &(*triangulation)));
-    (void)distributed_triangulation;
-
-    Assert(
-      distributed_triangulation != nullptr,
-      ExcMessage(
-        "Mesh refinement in a non-distributed triangulation is not supported "
-        "by the ParticleHandler class. Either insert particles after mesh "
-        "creation and do not refine afterwards, or use a distributed triangulation."));
-
     const auto callback_function =
       [this](const typename Triangulation<dim, spacedim>::cell_iterator
                              &cell_iterator,
@@ -2216,8 +2201,9 @@ namespace Particles
         return this->pack_callback(cell_iterator, cell_status);
       };
 
-    handle = distributed_triangulation->register_data_attach(
-      callback_function, /*returns_variable_size_data=*/true);
+    handle = const_cast<Triangulation<dim, spacedim> *>(&*triangulation)
+               ->register_data_attach(callback_function,
+                                      /*returns_variable_size_data=*/true);
   }
 
 
@@ -2256,21 +2242,6 @@ namespace Particles
   ParticleHandler<dim, spacedim>::notify_ready_to_unpack(
     const bool serialization)
   {
-    parallel::DistributedTriangulationBase<dim, spacedim>
-      *distributed_triangulation =
-        const_cast<parallel::DistributedTriangulationBase<dim, spacedim> *>(
-          dynamic_cast<
-            const parallel::DistributedTriangulationBase<dim, spacedim> *>(
-            &(*triangulation)));
-    (void)distributed_triangulation;
-
-    Assert(
-      distributed_triangulation != nullptr,
-      ExcMessage(
-        "Mesh refinement in a non-distributed triangulation is not supported "
-        "by the ParticleHandler class. Either insert particles after mesh "
-        "creation and do not refine afterwards, or use a distributed triangulation."));
-
     // First prepare container for insertion
     clear();
 
@@ -2293,8 +2264,8 @@ namespace Particles
             this->unpack_callback(cell_iterator, cell_status, range_iterator);
           };
 
-        distributed_triangulation->notify_ready_to_unpack(handle,
-                                                          callback_function);
+        const_cast<Triangulation<dim, spacedim> *>(&*triangulation)
+          ->notify_ready_to_unpack(handle, callback_function);
 
         // Reset handle and update global numbers.
         handle = numbers::invalid_unsigned_int;

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.