PID mesh adaptation quick fix

applications_tests/lethe-fluid/CMakeLists.txt

@@ -10,13 +10,13 @@ file(COPY box.msh DESTINATION "${CMAKE_CURRENT_BINARY_DIR}")
 file(COPY cylinder_manifold.msh DESTINATION "${CMAKE_CURRENT_BINARY_DIR}")
 file(COPY lpbf_benchmark_mesh.msh DESTINATION "${CMAKE_CURRENT_BINARY_DIR}")
 
-file(COPY cavity_ht_average_velocity_restart_files/restart.averagevelocities DESTINATION "${CMAKE_CURRENT_BINARY_DIR}/cavity_ht_average_velocity_restart.${_build_type}/mpirun=1/")
-file(COPY cavity_ht_average_velocity_restart_files/restart.pvdhandler DESTINATION "${CMAKE_CURRENT_BINARY_DIR}/cavity_ht_average_velocity_restart.${_build_type}/mpirun=1/")
-file(COPY cavity_ht_average_velocity_restart_files/restart.simulationcontrol DESTINATION "${CMAKE_CURRENT_BINARY_DIR}/cavity_ht_average_velocity_restart.${_build_type}/mpirun=1/")
-file(COPY cavity_ht_average_velocity_restart_files/restart.triangulation DESTINATION "${CMAKE_CURRENT_BINARY_DIR}/cavity_ht_average_velocity_restart.${_build_type}/mpirun=1/")
-file(COPY cavity_ht_average_velocity_restart_files/restart.triangulation.info DESTINATION "${CMAKE_CURRENT_BINARY_DIR}/cavity_ht_average_velocity_restart.${_build_type}/mpirun=1/")
-file(COPY cavity_ht_average_velocity_restart_files/restart.triangulation_fixed.data DESTINATION "${CMAKE_CURRENT_BINARY_DIR}/cavity_ht_average_velocity_restart.${_build_type}/mpirun=1/")
-file(COPY cavity_ht_average_velocity_restart_files/temperature_statistics.checkpoint DESTINATION "${CMAKE_CURRENT_BINARY_DIR}/cavity_ht_average_velocity_restart.${_build_type}/mpirun=1/")
+file(COPY cavity_ht_average_velocity_initial_condition_files/restart.averagevelocities DESTINATION "${CMAKE_CURRENT_BINARY_DIR}/cavity_ht_average_velocity_as_initial_condition.${_build_type}/mpirun=1/")
+file(COPY cavity_ht_average_velocity_initial_condition_files/restart.pvdhandler DESTINATION "${CMAKE_CURRENT_BINARY_DIR}/cavity_ht_average_velocity_as_initial_condition.${_build_type}/mpirun=1/")
+file(COPY cavity_ht_average_velocity_initial_condition_files/restart.simulationcontrol DESTINATION "${CMAKE_CURRENT_BINARY_DIR}/cavity_ht_average_velocity_as_initial_condition.${_build_type}/mpirun=1/")
+file(COPY cavity_ht_average_velocity_initial_condition_files/restart.triangulation DESTINATION "${CMAKE_CURRENT_BINARY_DIR}/cavity_ht_average_velocity_as_initial_condition.${_build_type}/mpirun=1/")
+file(COPY cavity_ht_average_velocity_initial_condition_files/restart.triangulation.info DESTINATION "${CMAKE_CURRENT_BINARY_DIR}/cavity_ht_average_velocity_as_initial_condition.${_build_type}/mpirun=1/")
+file(COPY cavity_ht_average_velocity_initial_condition_files/restart.triangulation_fixed.data DESTINATION "${CMAKE_CURRENT_BINARY_DIR}/cavity_ht_average_velocity_as_initial_condition.${_build_type}/mpirun=1/")
+file(COPY cavity_ht_average_velocity_initial_condition_files/temperature_statistics.checkpoint DESTINATION "${CMAKE_CURRENT_BINARY_DIR}/cavity_ht_average_velocity_as_initial_condition.${_build_type}/mpirun=1/")
 
 file(COPY mms-transient-conduction-restart_files/restart.pvdhandler DESTINATION "${CMAKE_CURRENT_BINARY_DIR}/mms-transient-conduction-restart_gls.${_build_type}/mpirun=1/")
 file(COPY mms-transient-conduction-restart_files/restart.simulationcontrol DESTINATION "${CMAKE_CURRENT_BINARY_DIR}/mms-transient-conduction-restart_gls.${_build_type}/mpirun=1/")

applications_tests/lethe-fluid/buoyant_flow_between_parallel_plates.output

@@ -8,7 +8,7 @@ Running on 1 MPI rank(s)...
 Steady iteration:        1/3
 *****************************
 L2 error velocity : 2.97964e-05
-L2 error temperature : 5.62951e-13
+L2 error temperature : 6.63477e-15
 
 *****************************
 Steady iteration:        2/3
@@ -18,7 +18,7 @@ Steady iteration:        2/3
    Volume of triangulation:      0.05
    Number of thermal degrees of freedom: 289
 L2 error velocity : 6.87211e-06
-L2 error temperature : 8.6248e-13
+L2 error temperature : 7.97464e-13
 
 *****************************
 Steady iteration:        3/3
@@ -28,12 +28,12 @@ Steady iteration:        3/3
    Volume of triangulation:      0.05
    Number of thermal degrees of freedom: 1089
 L2 error velocity : 6.84085e-06
-L2 error temperature : 1.46454e-12
+L2 error temperature : 1.45997e-12
 cells  error_velocity     error_pressure   
    64 2.979638e-05    - 6.154347e-01     - 
-  256 6.872113e-06 2.12 3.166237e-01  0.96 
+  256 6.872113e-06 2.12 3.166238e-01  0.96 
  1024 6.840850e-06 0.01 3.166244e-01 0.00 
 cells error_temperature  
-   64 5.629510e-13     - 
-  256 8.624800e-13 -0.62 
- 1024 1.464538e-12 -0.76 
+   64 6.634769e-15     - 
+  256 7.974639e-13 -6.91 
+ 1024 1.459965e-12 -0.87 

applications_tests/lethe-fluid/cavity_ht_average_velocity_as_initial_condition.mpirun=1.output

@@ -0,0 +1,55 @@
+Running on 1 MPI rank(s)...
+   Number of active cells:       4096
+   Number of degrees of freedom: 12675
+   Volume of triangulation:      1
+   Number of thermal degrees of freedom: 4225
+   Number of VOF degrees of freedom: 4225
+-------------------------------------------------
+Initial condition using average velocity profile
+-------------------------------------------------
+   Number of active cells:       4096
+   Number of degrees of freedom: 12675
+   Volume of triangulation:      1
+   Number of thermal degrees of freedom: 4225
+   Number of VOF degrees of freedom: 4225
+Temperature statistics on fluid : 
+	     Min : 0
+	     Max : 0.788675
+	 Average : 0.00769043
+	 Std-Dev : 0.0709978
+
+*******************************************************************************
+Transient iteration: 1        Time: 0.1      Time step: 0.1      CFL: 0       
+*******************************************************************************
+Temperature statistics on fluid : 
+	     Min : 2.7801e-17
+	     Max : 0.946974
+	 Average : 0.0360919
+	 Std-Dev : 0.136556
+
+*******************************************************************************
+Transient iteration: 2        Time: 0.2      Time step: 0.1      CFL: 5.27749 
+*******************************************************************************
+Temperature statistics on fluid : 
+	     Min : 2.98746e-16
+	     Max : 0.976436
+	 Average : 0.0532096
+	 Std-Dev : 0.165839
+
+*******************************************************************************
+Transient iteration: 3        Time: 0.3      Time step: 0.1      CFL: 15.8325 
+*******************************************************************************
+Temperature statistics on fluid : 
+	     Min : 1.99946e-15
+	     Max : 0.982787
+	 Average : 0.0654717
+	 Std-Dev : 0.182144
+
+*******************************************************************************
+Transient iteration: 4        Time: 0.4      Time step: 0.1      CFL: 15.8325 
+*******************************************************************************
+Temperature statistics on fluid : 
+	     Min : 1.0449e-14
+	     Max : 0.985964
+	 Average : 0.0758285
+	 Std-Dev : 0.194116

applications_tests/lethe-fluid/cavity_ht_average_velocity_restart.prm → applications_tests/lethe-fluid/cavity_ht_average_velocity_as_initial_condition.prm

@@ -20,7 +20,8 @@ end
 #---------------------------------------------------
 
 subsection initial conditions
-  set type = nodal
+  set type = average_velocity_profile
+
   subsection VOF
     set Function expression = if (x<0.5, 1, 0)
 
@@ -60,8 +61,6 @@ subsection multiphysics
   set fluid dynamics = false
   set heat transfer  = true
   set VOF            = true
-
-  set use time average velocity field = true
 end
 
 #---------------------------------------------------
@@ -84,7 +83,7 @@ subsection restart
   set checkpoint = false
   set frequency  = 1
   set filename   = restart
-  set restart    = true
+  set restart    = false
 end
 
 #---------------------------------------------------

doc/source/parameters/cfd/initial_conditions.rst

@@ -7,7 +7,7 @@ It is often necessary to set-up complex initial conditions when simulating trans
 .. code-block:: text
 
   subsection initial conditions
-    # Type of initial conditions. Choices are L2projection, viscous or nodal
+    # Type of initial conditions. Choices are L2projection, viscous, nodal or average_velocity_profile
     set type                = nodal
 
     # Kinematic viscosity for viscous initial conditions
@@ -46,10 +46,15 @@ It is often necessary to set-up complex initial conditions when simulating trans
         set alpha      = 0.5
       end
     end
+
+    subsection average velocity profile
+      set checkpoint folder    = ./
+      set checkpoint file name = restart
+    end
   end
 
 
-* The ``type`` parameter indicates which strategy is used to impose the initial conditions. The choices are : ``L2projection``, ``viscous``, ``nodal`` or ``ramp``.
+* The ``type`` parameter indicates which strategy is used to impose the initial conditions. The choices are : ``L2projection``, ``viscous``, ``nodal``, ``ramp`` or ``average_velocity_profile``.
 
 * The ``kinematic viscosity`` parameter controls the kinematic viscosity that is  used to solve the Navier-Stokes equations when the viscous initial condition is chosen.
 
@@ -95,3 +100,11 @@ Likewise, in the ``subection n``, the parameters for ramping on the ``n`` value
 
 .. math::
   n_{i+1} = n_i + \alpha (n_{\text{end}} - n_i)
+
+
+* The subsection ``average velocity profile`` uses the time averaged fluid velocity calculated in a previous simulation as an initial condition. This is useful when the flow dynamics and the subphysics reach a pseudo-steady state at different time scales. Physics can then be run independently, one to solve for the fluid dynamics and one for the subphysics. To use this feature, the user should launch a simulation with the fluid mechanics solver while using the time averaging and checkpointing feature. Once the time average of the velocity field is sufficiently established, the simulation should be stopped and a new simulation can be restarted without the fluid mechanics solver. The subphysics can then be solved using a larger time step.
+
+.. important::
+   * If only an auxiliary physic must be solved without the fluid dynamics, ``set fluid dynamics = false`` needs to be specified in the ``multiphysics`` section. The average velocity field will then be used for the whole duration of the simulation.
+   * This feature uses the checkpoint mechanism to load the time averaged velocity field. Make sure to activate checkpointing in the restart section of the first simulation. 
+   * The same mesh needs to be used for the fluid dynamics and the auxiliary physics simulations. The mesh should not be modified between the two simulations.

doc/source/parameters/cfd/multiphysics.rst

@@ -23,7 +23,6 @@ This subsection defines the multiphysics interface of Lethe and enables the solu
     # Cahn-Hilliard equations
     set cahn hilliard                   = false
 
-    set use time average velocity field = false
   end
 
 
@@ -59,9 +58,5 @@ This subsection defines the multiphysics interface of Lethe and enables the solu
 
   The VOF solver is used in the example :doc:`../../examples/multiphysics/dam-break/dam-break`.
 
-* ``use time average velocity field``:  controls if the advection of the subphysics is done using the average velocity field instead of the current velocity field. This is useful when the flow dynamics and the subphysics reach a pseudo-steady state at different time scales. To use this feature, the user should launch a simulation with the fluid mechanics solver and use the time averaging feature. Once the time average of the velocity field is sufficiently established, the simulation should be stopped and restarted without the fluid mechanics solver. The subphysics can then be solved using a larger time step.
-
-.. important::
-   The subphysic must also be enabled in the first part of the simulation (for example, ``set heat transfer = true``).
 
 

include/core/mesh_controller.h

@@ -31,6 +31,8 @@ class MeshController
    */
   MeshController(const unsigned int target_number_of_elements)
     : target_number_of_elements(target_number_of_elements)
+    , previous_number_of_elements(0)
+    , previous_mesh_control_error(0.0)
   {}
   virtual ~MeshController()
   {}

include/core/parameters_multiphysics.h

@@ -232,8 +232,6 @@ namespace Parameters
     bool VOF;
     bool cahn_hilliard;
 
-    bool use_time_average_velocity_field;
-
     // subparameters for heat_transfer
     bool viscous_dissipation;
     bool buoyancy_force;

include/core/simulation_control.h

@@ -475,8 +475,20 @@ class SimulationControl
 
   void
   save(std::string filename);
+
+  /**
+   * @brief Reads the simulation control information from the checkpoint file and updates the time step vector, the CFL value, the time and the iteration number.
+   */
   void
   read(std::string filename);
+
+  /**
+   * @brief Reads and returns the simulation control information from the checkpoint file filename without updating the simulation control information.
+   *
+   * @return A vector containing the last checkpointed file and time step.
+   */
+  std::vector<double>
+  get_checkpointed_simulation_control_info(std::string filename);
 };
 
 

include/solvers/initial_conditions.h

@@ -34,7 +34,8 @@ namespace Parameters
     L2projection,
     viscous,
     nodal,
-    ramp
+    ramp,
+    average_velocity_profile
   };
 
   struct Ramp_n
@@ -108,6 +109,10 @@ namespace Parameters
     Functions::ParsedFunction<dim> cahn_hilliard =
       Functions::ParsedFunction<dim>(2);
 
+    // Path to the checkpointed average velocity profile
+    std::string average_velocity_folder;
+    std::string average_velocity_file_name;
+
     void
     declare_parameters(ParameterHandler &prm);
     void
@@ -120,11 +125,13 @@ namespace Parameters
   {
     prm.enter_subsection("initial conditions");
     {
-      prm.declare_entry("type",
-                        "nodal",
-                        Patterns::Selection("L2projection|viscous|nodal|ramp"),
-                        "Type of initial condition"
-                        "Choices are <L2projection|viscous|nodal|ramp>.");
+      prm.declare_entry(
+        "type",
+        "nodal",
+        Patterns::Selection(
+          "L2projection|viscous|nodal|ramp|average_velocity_profile"),
+        "Type of initial condition"
+        "Choices are <L2projection|viscous|nodal|ramp|average_velocity_profile>.");
       prm.enter_subsection("uvwp");
       uvwp.declare_parameters(prm, dim + 1);
       prm.leave_subsection();
@@ -163,6 +170,18 @@ namespace Parameters
       cahn_hilliard.declare_parameters(prm, 2);
       prm.leave_subsection();
 
+      prm.enter_subsection("average velocity profile");
+      prm.declare_entry(
+        "checkpoint folder",
+        "./",
+        Patterns::FileName(),
+        "the path leading to the checkpointed average velocity profile");
+      prm.declare_entry("checkpoint file name",
+                        "restart",
+                        Patterns::FileName(),
+                        "checkpoint file name");
+      prm.leave_subsection();
+
       ramp.declare_parameters(prm);
     }
     prm.leave_subsection();
@@ -183,6 +202,8 @@ namespace Parameters
         type = InitialConditionType::nodal;
       else if (op == "ramp")
         type = InitialConditionType::ramp;
+      else if (op == "average_velocity_profile")
+        type = InitialConditionType::average_velocity_profile;
 
       kinematic_viscosity = prm.get_double("kinematic viscosity");
       prm.enter_subsection("uvwp");
@@ -214,6 +235,12 @@ namespace Parameters
       prm.enter_subsection("cahn hilliard");
       cahn_hilliard.parse_parameters(prm);
       prm.leave_subsection();
+
+      prm.enter_subsection("average velocity profile");
+      average_velocity_folder    = prm.get("checkpoint folder");
+      average_velocity_file_name = prm.get("checkpoint file name");
+
+      prm.leave_subsection();
     }
     prm.leave_subsection();
   }

include/solvers/navier_stokes_base.h

@@ -382,11 +382,17 @@ class NavierStokesBase : public PhysicsSolver<VectorType>
   refine_mesh_uniform();
 
   /**
-   * @brief read_checkpoint
+   * @brief Restart a previous simulation from a checkpoint file.
    */
   virtual void
   read_checkpoint();
 
+  /**
+   * @brief Read the solution from a checkpoint file and set it as the current solution.
+   */
+  virtual void
+  set_solution_from_checkpoint(std::string checkpoint_file_prefix);
+
   /**
    * @brief set_nodal_values
    * Set the nodal values of velocity and pressure

source/core/parameters_multiphysics.cc

@@ -55,12 +55,6 @@ Parameters::Multiphysics::declare_parameters(ParameterHandler &prm)
                       Patterns::Bool(),
                       "Cahn-Hilliard calculation <true|false>");
 
-    prm.declare_entry(
-      "use time average velocity field",
-      "false",
-      Patterns::Bool(),
-      "Use the average velocity field in subphysics, instead of the present velocity field <true|false>");
-
     // subparameters for heat_transfer
     prm.declare_entry("viscous dissipation",
                       "false",
@@ -88,8 +82,6 @@ Parameters::Multiphysics::parse_parameters(ParameterHandler &prm)
     tracer         = prm.get_bool("tracer");
     VOF            = prm.get_bool("VOF");
     cahn_hilliard  = prm.get_bool("cahn hilliard");
-    use_time_average_velocity_field =
-      prm.get_bool("use time average velocity field");
 
     // subparameter for heat_transfer
     viscous_dissipation = prm.get_bool("viscous dissipation");