]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Fix p:f:T::load() 16940/head
authorPeter Munch <peterrmuench@gmail.com>
Mon, 29 Apr 2024 21:17:43 +0000 (23:17 +0200)
committerPeter Munch <peterrmuench@gmail.com>
Wed, 1 May 2024 14:57:47 +0000 (16:57 +0200)
include/deal.II/distributed/shared_tria.h
include/deal.II/distributed/tria_base.h
source/distributed/fully_distributed_tria.cc
source/distributed/tria_base.cc
tests/sharedtria/serialization_01.cc
tests/sharedtria/serialization_01.mpirun=2.output

index d10b2ea256f25087363d997b99cb3ff9ef4304fa..14b626734cf6e42ffb157a214f790061141375cf 100644 (file)
@@ -319,6 +319,21 @@ namespace parallel
       copy_triangulation(
         const dealii::Triangulation<dim, spacedim> &other_tria) override;
 
+#  ifdef DOXYGEN
+      /**
+       * Write and read the data of this object from a stream for the purpose
+       * of serialization. using the [BOOST serialization
+       * library](https://www.boost.org/doc/libs/1_74_0/libs/serialization/doc/index.html).
+       */
+      template <class Archive>
+      void
+      serialize(Archive &archive, const unsigned int version);
+#  else
+      // This macro defines the serialize() method that is compatible with
+      // the templated save() and load() method that have been implemented.
+      BOOST_SERIALIZATION_SPLIT_MEMBER()
+#  endif
+
       /**
        * Save the triangulation into the given file. This file needs to be
        * reachable from all nodes in the computation on a shared network file
index 6bd1ae0b9f33d4777a248b93ba181d3685a46216..0156954533fa3b1932372bc394f2a1fa678276ed 100644 (file)
@@ -304,6 +304,15 @@ namespace parallel
     virtual types::coarse_cell_id
     n_global_coarse_cells() const override;
 
+    /**
+     * Reset this triangulation to an empty state by deleting all data.
+     *
+     * Note that this operation is only allowed if no subscriptions to this
+     * object exist any more, such as DoFHandler objects using it.
+     */
+    virtual void
+    clear() override;
+
   protected:
     /**
      * MPI communicator to be used for the triangulation. We create a unique
@@ -460,15 +469,6 @@ namespace parallel
                  smooth_grid = (dealii::Triangulation<dim, spacedim>::none),
       const bool check_for_distorted_cells = false);
 
-    /**
-     * Reset this triangulation into a virgin state by deleting all data.
-     *
-     * Note that this operation is only allowed if no subscriptions to this
-     * object exist any more, such as DoFHandler objects using it.
-     */
-    virtual void
-    clear() override;
-
     using cell_iterator =
       typename dealii::Triangulation<dim, spacedim>::cell_iterator;
 
index 89434c72386f71eb575e2dcec7c4c52cc04274fc..33b78780a616cfa06108075a63d7518b843f4cbc 100644 (file)
@@ -598,21 +598,6 @@ namespace parallel
       Assert(this->n_cells() == 0,
              ExcMessage("load() only works if the Triangulation is empty!"));
 
-      // Compute global offset for each rank.
-      unsigned int n_locally_owned_cells = this->n_locally_owned_active_cells();
-
-      unsigned int global_first_cell = 0;
-
-      int ierr = MPI_Exscan(&n_locally_owned_cells,
-                            &global_first_cell,
-                            1,
-                            MPI_UNSIGNED,
-                            MPI_SUM,
-                            this->mpi_communicator);
-      AssertThrowMPI(ierr);
-
-      global_first_cell *= sizeof(unsigned int);
-
 
       unsigned int version, numcpus, attached_count_fixed,
         attached_count_variable, n_global_active_cells;
@@ -628,8 +613,6 @@ namespace parallel
 
       AssertThrow(version == 4,
                   ExcMessage("Incompatible version found in .info file."));
-      Assert(this->n_global_active_cells() == n_global_active_cells,
-             ExcMessage("Number of global active cells differ!"));
 
       // Load description and construct the triangulation.
       {
@@ -709,6 +692,24 @@ namespace parallel
         this->create_triangulation(construction_data);
       }
 
+      // Compute global offset for each rank.
+      unsigned int n_locally_owned_cells = this->n_locally_owned_active_cells();
+
+      unsigned int global_first_cell = 0;
+
+      int ierr = MPI_Exscan(&n_locally_owned_cells,
+                            &global_first_cell,
+                            1,
+                            MPI_UNSIGNED,
+                            MPI_SUM,
+                            this->mpi_communicator);
+      AssertThrowMPI(ierr);
+
+      global_first_cell *= sizeof(unsigned int);
+
+      Assert(this->n_global_active_cells() == n_global_active_cells,
+             ExcMessage("Number of global active cells differ!"));
+
       // clear all of the callback data, as explained in the documentation of
       // register_data_attach()
       this->cell_attached_data.n_attached_data_sets = 0;
index cdc17609a7f907c432a7527ef5285f8e93a46086..e044ba36d1c6f3c21661a41c31278d891a2cfdf0 100644 (file)
@@ -684,9 +684,11 @@ namespace parallel
 
   template <int dim, int spacedim>
   DEAL_II_CXX20_REQUIRES((concepts::is_valid_dim_spacedim<dim, spacedim>))
-  void DistributedTriangulationBase<dim, spacedim>::clear()
+  void TriangulationBase<dim, spacedim>::clear()
   {
     dealii::Triangulation<dim, spacedim>::clear();
+
+    number_cache = {};
   }
 
 
index 2fc47b4ed4b05a920a1065526b71d7f8f14b0de8..4e2183e6bfe8736362bcb6d98164437d9fd04bb4 100644 (file)
@@ -50,7 +50,7 @@ public:
 
 template <int dim, typename TriangulationType>
 void
-test(TriangulationType &triangulation)
+test(TriangulationType &triangulation, TriangulationType &triangulation_clean)
 {
   DoFHandler<dim> dof_handler(triangulation);
   dof_handler.distribute_dofs(FE_Q<dim>(2));
@@ -65,27 +65,49 @@ test(TriangulationType &triangulation)
 
   print_statistics(triangulation, false);
 
-  std::stringstream stream;
+  std::stringstream stream_0, stream_1;
   {
-    boost::archive::text_oarchive oa(stream);
-    oa << triangulation;
-    oa << vector;
+    boost::archive::text_oarchive oa_0(stream_0);
+    boost::archive::text_oarchive oa_1(stream_1);
+    oa_0 << triangulation;
+    oa_0 << vector;
+    oa_1 << triangulation;
+    oa_1 << vector;
   }
 
   triangulation.clear();
 
   {
-    boost::archive::text_iarchive ia(stream);
-    ia >> triangulation;
-    ia >> vector_loaded;
+    // load into existing tringulation
+    {
+      boost::archive::text_iarchive ia(stream_0);
+      ia >> triangulation;
+      ia >> vector_loaded;
+    }
+    print_statistics(triangulation, false);
+
+    // Verify that error is 0.
+    VectorType error(vector);
+    error.add(-1, vector_loaded);
+
+    deallog << (error.linfty_norm() < 1e-16 ? "PASSED" : "FAILED") << std::endl;
   }
-  print_statistics(triangulation, false);
-
-  // Verify that error is 0.
-  VectorType error(vector);
-  error.add(-1, vector_loaded);
 
-  deallog << (error.linfty_norm() < 1e-16 ? "PASSED" : "FAILED") << std::endl;
+  {
+    // load into new tringulation
+    {
+      boost::archive::text_iarchive ia(stream_1);
+      ia >> triangulation_clean;
+      ia >> vector_loaded;
+    }
+    print_statistics(triangulation_clean, false);
+
+    // Verify that error is 0.
+    VectorType error(vector);
+    error.add(-1, vector_loaded);
+
+    deallog << (error.linfty_norm() < 1e-16 ? "PASSED" : "FAILED") << std::endl;
+  }
 }
 
 int
@@ -106,7 +128,14 @@ main(int argc, char *argv[])
     GridGenerator::hyper_cube(triangulation);
     triangulation.refine_global(3);
 
-    test<dim>(triangulation);
+
+    parallel::shared::Triangulation<dim> triangulation_other(
+      MPI_COMM_WORLD,
+      ::Triangulation<dim>::none,
+      false,
+      parallel::shared::Triangulation<dim>::partition_zorder);
+
+    test<dim>(triangulation, triangulation_other);
   }
   deallog.pop();
 
@@ -122,7 +151,13 @@ main(int argc, char *argv[])
     GridGenerator::hyper_cube(triangulation);
     triangulation.refine_global(3);
 
-    test<dim>(triangulation);
+    parallel::shared::Triangulation<dim> triangulation_other(
+      MPI_COMM_WORLD,
+      ::Triangulation<dim>::none,
+      false,
+      parallel::shared::Triangulation<dim>::partition_zorder);
+
+    test<dim>(triangulation, triangulation_other);
   }
   deallog.pop();
 }
index 985cf7deff47fa3fcb1e3640acd805cdedf7aaba..4edb900d08a9e7d9a20e536f3555286b6e695790 100644 (file)
@@ -10,6 +10,12 @@ DEAL:0:2d::n_active_cells:            64
 DEAL:0:2d::has_hanging_nodes:         false
 DEAL:0:2d::
 DEAL:0:2d::PASSED
+DEAL:0:2d::n_levels:                  4
+DEAL:0:2d::n_cells:                   85
+DEAL:0:2d::n_active_cells:            64
+DEAL:0:2d::has_hanging_nodes:         false
+DEAL:0:2d::
+DEAL:0:2d::PASSED
 DEAL:0:3d::n_levels:                  4
 DEAL:0:3d::n_cells:                   585
 DEAL:0:3d::n_active_cells:            512
@@ -21,6 +27,12 @@ DEAL:0:3d::n_active_cells:            512
 DEAL:0:3d::has_hanging_nodes:         false
 DEAL:0:3d::
 DEAL:0:3d::PASSED
+DEAL:0:3d::n_levels:                  4
+DEAL:0:3d::n_cells:                   585
+DEAL:0:3d::n_active_cells:            512
+DEAL:0:3d::has_hanging_nodes:         false
+DEAL:0:3d::
+DEAL:0:3d::PASSED
 
 DEAL:1:2d::n_levels:                  4
 DEAL:1:2d::n_cells:                   85
@@ -33,6 +45,12 @@ DEAL:1:2d::n_active_cells:            64
 DEAL:1:2d::has_hanging_nodes:         false
 DEAL:1:2d::
 DEAL:1:2d::PASSED
+DEAL:1:2d::n_levels:                  4
+DEAL:1:2d::n_cells:                   85
+DEAL:1:2d::n_active_cells:            64
+DEAL:1:2d::has_hanging_nodes:         false
+DEAL:1:2d::
+DEAL:1:2d::PASSED
 DEAL:1:3d::n_levels:                  4
 DEAL:1:3d::n_cells:                   585
 DEAL:1:3d::n_active_cells:            512
@@ -44,4 +62,10 @@ DEAL:1:3d::n_active_cells:            512
 DEAL:1:3d::has_hanging_nodes:         false
 DEAL:1:3d::
 DEAL:1:3d::PASSED
+DEAL:1:3d::n_levels:                  4
+DEAL:1:3d::n_cells:                   585
+DEAL:1:3d::n_active_cells:            512
+DEAL:1:3d::has_hanging_nodes:         false
+DEAL:1:3d::
+DEAL:1:3d::PASSED
 

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.