Nitsche checkpoint particles

CHANGELOG.md

@@ -3,6 +3,12 @@
 All notable changes to the Lethe project will be documented in this file.
 The format is based on [Keep a Changelog](http://keepachangelog.com/).
 
+## [Master] - 2024-07-09
+
+### Fixed
+
+- MINOR The lethe-fluid-nitsche solver could not checkpoint once it had been restarted. This has been fixed by ensuring that the particles are properly checkpointed and read from the restart files. [#1192](https://github.com/chaos-polymtl/lethe/pull/1192)
+
 ## [Master] - 2024-07-08
 
 ### Changed

applications/lethe-fluid-nitsche/nitsche_navier_stokes.cc

@@ -43,9 +43,6 @@ main(int argc, char *argv[])
           prm.parse_input(argv[1]);
           NSparam.parse(prm);
 
-          AssertThrow(NSparam.nitsche->number_solids > 0,
-                      NoSolidWarning("gls_nitsche_navier_stokes_22"));
-
           GLSNitscheNavierStokesSolver<2> problem_22(NSparam);
           problem_22.solve();
         }
@@ -59,9 +56,6 @@ main(int argc, char *argv[])
           prm.parse_input(argv[1]);
           NSparam.parse(prm);
 
-          AssertThrow(NSparam.nitsche->number_solids > 0,
-                      NoSolidWarning("gls_nitsche_navier_stokes_22"));
-
           GLSNitscheNavierStokesSolver<3> problem_33(NSparam);
           problem_33.solve();
         }

applications_tests/lethe-fluid-nitsche/CMakeLists.txt

@@ -13,6 +13,8 @@ file(COPY two-bar-mixer-restart/restart.triangulation_variable.data DESTINATION
 file(COPY two-bar-mixer-restart/restart.triangulation.info DESTINATION "${CMAKE_CURRENT_BINARY_DIR}/two-bar-mixer-restart.${_build_type}/mpirun=1/")
 file(COPY two-bar-mixer-restart/restart_solid_particles_00.pvdhandler DESTINATION "${CMAKE_CURRENT_BINARY_DIR}/two-bar-mixer-restart.${_build_type}/mpirun=1/")
 file(COPY two-bar-mixer-restart/restart_solid_particles_01.pvdhandler DESTINATION "${CMAKE_CURRENT_BINARY_DIR}/two-bar-mixer-restart.${_build_type}/mpirun=1/")
+file(COPY two-bar-mixer-restart/restart_solid_00.particles DESTINATION "${CMAKE_CURRENT_BINARY_DIR}/two-bar-mixer-restart.${_build_type}/mpirun=1/")
+file(COPY two-bar-mixer-restart/restart_solid_01.particles DESTINATION "${CMAKE_CURRENT_BINARY_DIR}/two-bar-mixer-restart.${_build_type}/mpirun=1/")
 file(COPY two-bar-mixer-restart/restart_solid_00.triangulation DESTINATION "${CMAKE_CURRENT_BINARY_DIR}/two-bar-mixer-restart.${_build_type}/mpirun=1/")
 file(COPY two-bar-mixer-restart/restart_solid_00.triangulation.info DESTINATION "${CMAKE_CURRENT_BINARY_DIR}/two-bar-mixer-restart.${_build_type}/mpirun=1/")
 file(COPY two-bar-mixer-restart/restart_solid_00.triangulation_fixed.data DESTINATION "${CMAKE_CURRENT_BINARY_DIR}/two-bar-mixer-restart.${_build_type}/mpirun=1/")

applications_tests/lethe-fluid-nitsche/two-bar-mixer-start.prm → applications_tests/lethe-fluid-nitsche/generators/two-bar-mixer-start.prm

@@ -72,6 +72,14 @@ subsection initial conditions
   end
 end
 
+#---------------------------------------------------
+# Boundary Conditions
+#---------------------------------------------------
+
+subsection boundary conditions
+  set number = 0
+end
+
 #---------------------------------------------------
 # Physical Properties
 #---------------------------------------------------

applications_tests/lethe-fluid-nitsche/two-bar-mixer-restart.prm

@@ -33,7 +33,7 @@ end
 
 subsection restart
   # Enable checkpointing
-  set checkpoint = false
+  set checkpoint = true
 
   # Enable restart
   set restart = true

applications_tests/lethe-fluid-nitsche/two-bar-mixer-restart/restart_solid_00.particles

@@ -0,0 +1 @@
+0 0 16 4 16

applications_tests/lethe-fluid-nitsche/two-bar-mixer-restart/restart_solid_01.particles

@@ -0,0 +1 @@
+0 0 16 4 16

source/core/solid_base.cc

@@ -722,14 +722,28 @@ template <int dim, int spacedim>
 void
 SolidBase<dim, spacedim>::write_checkpoint(std::string prefix)
 {
-#if (DEAL_II_VERSION_MAJOR < 10 && DEAL_II_VERSION_MINOR < 4)
-  parallel::distributed::
-    SolutionTransfer<dim, GlobalVectorType, DoFHandler<dim, spacedim>>
-      system_trans_vectors(this->displacement_dh);
-#else
+  // First, manage the serialization of the particle information.
+  // The particles are attached to the background triangulation and not
+  // the solid triangulation. Consequently, they only need to be prepared
+  // for serialization and will be serialized when the background triangulation
+  // itself is serialized by the main solver.
+  std::string                   prefix_particles = prefix + "_particles";
+  std::ostringstream            oss;
+  boost::archive::text_oarchive oa(oss, boost::archive::no_header);
+  oa << *(this->solid_particle_handler);
+
+  // Write additional particle information for deserialization
+  std::string   particle_filename = prefix + ".particles";
+  std::ofstream output(particle_filename.c_str());
+  output << oss.str() << std::endl;
+
+  //   Prepare particle handler for serialization
+  this->solid_particle_handler->prepare_for_serialization();
+
+  // Prepare serialization of the solid triangulation which is
+  // only used for visualization purposes
   parallel::distributed::SolutionTransfer<dim, GlobalVectorType, spacedim>
     system_trans_vectors(this->displacement_dh);
-#endif
 
   std::vector<const GlobalVectorType *> sol_set_transfer;
   displacement_relevant = displacement;
@@ -773,14 +787,8 @@ SolidBase<dim, spacedim>::read_checkpoint(std::string prefix)
   std::vector<GlobalVectorType *> x_system(1);
   x_system[0] = &(displacement);
 
-#if (DEAL_II_VERSION_MAJOR < 10 && DEAL_II_VERSION_MINOR < 4)
-  parallel::distributed::
-    SolutionTransfer<dim, GlobalVectorType, DoFHandler<dim, spacedim>>
-      system_trans_vectors(this->displacement_dh);
-#else
   parallel::distributed::SolutionTransfer<dim, GlobalVectorType, spacedim>
     system_trans_vectors(this->displacement_dh);
-#endif
 
   system_trans_vectors.deserialize(x_system);
   displacement_relevant = displacement;
@@ -791,8 +799,28 @@ SolidBase<dim, spacedim>::read_checkpoint(std::string prefix)
   if (!param->solid_mesh.simplex)
     move_solid_triangulation_with_displacement();
 
-  // We did not checkpoint particles, we re-create them from scratch
-  setup_particles();
+  // Gather particle serialization information and deserialize them
+  std::string   particle_filename = prefix + ".particles";
+  std::ifstream input(particle_filename.c_str());
+  AssertThrow(input, ExcFileNotOpen(particle_filename));
+
+  std::string buffer;
+  std::getline(input, buffer);
+  std::istringstream            iss(buffer);
+  boost::archive::text_iarchive ia(iss, boost::archive::no_header);
+
+
+  // At this stage, the pointer for the particle handler has not been allocated
+  // it thus needs to be done by calling the setup
+  // TODO, this function should be called in the constructor
+  this->setup_particles_handler();
+  ia >> *(this->solid_particle_handler);
+
+  this->solid_particle_handler->deserialize();
+
+  // Number of particles is reset to the number of particles in the
+  // solid_particle_handler after deserialization
+  initial_number_of_particles = solid_particle_handler->n_global_particles();
 }
 
 

source/solvers/gls_nitsche_navier_stokes.cc

@@ -114,13 +114,8 @@ GLSNitscheNavierStokesSolver<dim, spacedim>::assemble_nitsche_restriction()
           local_rhs    = 0;
 
 
-#if (DEAL_II_VERSION_MAJOR < 10 && DEAL_II_VERSION_MINOR < 4)
-          const auto &cell =
-            particle->get_surrounding_cell(*this->triangulation);
-#else
-          const auto &cell = particle->get_surrounding_cell();
-#endif
-          double h_cell = 0;
+          const auto &cell   = particle->get_surrounding_cell();
+          double      h_cell = 0;
           if (dim == 2)
             h_cell = std::sqrt(4. * cell->measure() / M_PI) /
                      this->velocity_fem_degree;
@@ -977,9 +972,16 @@ template <int dim, int spacedim>
 void
 GLSNitscheNavierStokesSolver<dim, spacedim>::write_checkpoint()
 {
-  // Prepare solid particles for serialization
-  for (auto solid : solids)
-    solid->get_solid_particle_handler()->prepare_for_serialization();
+  // Write solid base checkpoint
+  for (unsigned int i_solid = 0; i_solid < solids.size(); ++i_solid)
+    {
+      std::string filename =
+        this->simulation_parameters.simulation_control.output_folder +
+        this->simulation_parameters.restart_parameters.filename + "_solid_" +
+        Utilities::int_to_string(i_solid, 2);
+
+      this->solids[i_solid]->write_checkpoint(filename);
+    }
 
   // Call regular checkpointing routine
   this->GLSNavierStokesSolver<spacedim>::write_checkpoint();
@@ -1005,16 +1007,6 @@ GLSNitscheNavierStokesSolver<dim, spacedim>::write_checkpoint()
             Utilities::int_to_string(i_solid, 2));
         }
     }
-
-  // Write solid base checkpoint
-  for (unsigned int i_solid = 0; i_solid < solids.size(); ++i_solid)
-    {
-      std::string filename =
-        this->simulation_parameters.simulation_control.output_folder +
-        this->simulation_parameters.restart_parameters.filename + "_solid_" +
-        Utilities::int_to_string(i_solid, 2);
-      this->solids[i_solid]->write_checkpoint(filename);
-    }
 }
 
 template <int dim, int spacedim>
@@ -1024,10 +1016,12 @@ GLSNitscheNavierStokesSolver<dim, spacedim>::read_checkpoint()
   this->GLSNavierStokesSolver<spacedim>::read_checkpoint();
 
   TimerOutput::Scope t(this->computing_timer,
-                       "Reset Nitsche solid mesh and particles");
-
+                       "Reload Nitsche solid mesh and particles");
 
   // Reload initial configurations
+  // This must be done after the background triangulation is read
+  // because the solid particle_handler is deserialized physically
+  // after the triangulation is read.
   for (unsigned int i_solid = 0; i_solid < solids.size(); ++i_solid)
     {
       std::string prefix =
@@ -1038,6 +1032,7 @@ GLSNitscheNavierStokesSolver<dim, spacedim>::read_checkpoint()
       solids[i_solid]->read_checkpoint(prefix);
     }
 
+
   // Reload particle and solid pvd handlers
   std::string prefix =
     this->simulation_parameters.simulation_control.output_folder +