Read external fields from file in restarted simulation (#4547)

* Read external fields from file also in restarted simulation

* [pre-commit.ci] auto fixes from pre-commit.com hooks

for more information, see https://pre-commit.ci

* Add automated test with restart

* Add function to avoid duplication

* Fix restart test

* Move function outside of #ifdef

* Fix restart test

* Update Source/WarpX.H

* Update Source/Initialization/WarpXInitData.cpp

* Clarify documentation on external fields

* Add docstring

---------

Co-authored-by: pre-commit-ci[bot] <66853113+pre-commit-ci[bot]@users.noreply.github.com>

Author: RemiLehe

Docs/source/usage/parameters.rst

@@ -1508,9 +1508,15 @@ Laser initialization
 External fields
 ---------------
 
-Grid initialization
+Applied to the grid
 ^^^^^^^^^^^^^^^^^^^
 
+The external fields defined with input parameters that start with ``warpx.B_ext_grid_init_`` or ``warpx.E_ext_grid_init_``
+are applied to the grid directly. In particular, these fields can be seen in the diagnostics that output the fields on the grid.
+
+    - When using an **electromagnetic** field solver, these fields are applied to the grid at the beginning of the simulation, and serve as initial condition for the Maxwell solver.
+    - When using an **electrostatic** or **magnetostatic** field solver, these fields are added to the fields computed by the Poisson solver, at each timestep.
+
 * ``warpx.B_ext_grid_init_style`` (string) optional
     This parameter determines the type of initialization for the external
     magnetic field. By default, the
@@ -1597,6 +1603,9 @@ Grid initialization
 Applied to Particles
 ^^^^^^^^^^^^^^^^^^^^
 
+The external fields defined with input parameters that start with ``warpx.B_ext_particle_init_`` or ``warpx.E_ext_particle_init_``
+are applied to the particles directly, at each timestep. As a results, these fields **cannot** be seen in the diagnostics that output the fields on the grid.
+
 * ``particles.E_ext_particle_init_style`` & ``particles.B_ext_particle_init_style`` (string) optional (default "none")
     These parameters determine the type of the external electric and
     magnetic fields respectively that are applied directly to the particles at every timestep.
@@ -1657,6 +1666,9 @@ Applied to Particles
 Applied to Cold Relativistic Fluids
 ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
 
+The external fields defined with input parameters that start with ``warpx.B_ext_init_`` or ``warpx.E_ext_init_``
+are applied to the fluids directly, at each timestep. As a results, these fields **cannot** be seen in the diagnostics that output the fields on the grid.
+
 * ``<fluid_species_name>.E_ext_init_style`` & ``<fluid_species_name>.B_ext_init_style`` (string) optional (default "none")
     These parameters determine the type of the external electric and
     magnetic fields respectively that are applied directly to the cold relativistic fluids at every timestep.

Examples/Tests/LoadExternalField/inputs_rz

@@ -46,7 +46,7 @@ algo.particle_shape = 1
 #################################
 ############ PLASMA #############
 #################################
-particles.species_names = proton #electron
+particles.species_names = proton
 proton.injection_style = "SingleParticle"
 proton.single_particle_pos = 0.0  0.2  2.5
 proton.single_particle_u = 9.506735958279367e-05  0.0  0.00013435537232359165
@@ -55,14 +55,11 @@ proton.do_not_deposit = 1
 proton.mass = m_p
 proton.charge = q_e
 
-#electron.injection_style = "SingleParticle"
-#electron.single_particle_pos = 0.0  0.2  2.5
-#electron.single_particle_u = 0.0  0.0  0.0
-#electron.single_particle_weight = 1.0
-#electron.mass = 1.0
-#electron.charge = -q_e*1.0e-20
-
 # Diagnostics
-diagnostics.diags_names = diag1
+diagnostics.diags_names = diag1 chk
 diag1.intervals = 300
 diag1.diag_type = Full
+
+chk.intervals = 150
+chk.diag_type = Full
+chk.format = checkpoint

Regression/WarpX-tests.ini

@@ -2214,11 +2214,12 @@ analysisRoutine = Examples/Tests/resampling/analysis_leveling_thinning.py
 [LoadExternalFieldRZ]
 buildDir = .
 inputFile = Examples/Tests/LoadExternalField/inputs_rz
-runtime_params = warpx.abort_on_warning_threshold=medium
+runtime_params = warpx.abort_on_warning_threshold=medium chk.file_prefix=LoadExternalFieldRZ_chk chk.file_min_digits=5
 dim = 2
 addToCompileString = USE_RZ=TRUE
 cmakeSetupOpts = -DWarpX_DIMS=RZ -DWarpX_OPENPMD=ON
-restartTest = 0
+restartTest = 1
+restartFileNum = 150
 useMPI = 1
 numprocs = 1
 useOMP = 1

Source/Initialization/WarpXInitData.cpp

@@ -728,6 +728,9 @@ void
 WarpX::PostRestart ()
 {
     mypc->PostRestart();
+    for (int lev = 0; lev <= maxLevel(); ++lev) {
+        LoadExternalFieldsFromFile(lev);
+    }
 }
 
 
@@ -901,34 +904,7 @@ WarpX::InitLevelData (int lev, Real /*time*/)
        }
     }
 
-    // Reading external fields from data file
-    if (m_p_ext_field_params->B_ext_grid_type == ExternalFieldType::read_from_file) {
-
-#if defined(WARPX_DIM_RZ)
-        WARPX_ALWAYS_ASSERT_WITH_MESSAGE(n_rz_azimuthal_modes == 1,
-                                         "External field reading is not implemented for more than one RZ mode (see #3829)");
-        ReadExternalFieldFromFile(m_p_ext_field_params->external_fields_path, Bfield_fp_external[lev][0].get(), "B", "r");
-        ReadExternalFieldFromFile(m_p_ext_field_params->external_fields_path, Bfield_fp_external[lev][1].get(), "B", "t");
-        ReadExternalFieldFromFile(m_p_ext_field_params->external_fields_path, Bfield_fp_external[lev][2].get(), "B", "z");
-#else
-        ReadExternalFieldFromFile(m_p_ext_field_params->external_fields_path, Bfield_fp_external[lev][0].get(), "B", "x");
-        ReadExternalFieldFromFile(m_p_ext_field_params->external_fields_path, Bfield_fp_external[lev][1].get(), "B", "y");
-        ReadExternalFieldFromFile(m_p_ext_field_params->external_fields_path, Bfield_fp_external[lev][2].get(), "B", "z");
-#endif
-    }
-    if (m_p_ext_field_params->E_ext_grid_type == ExternalFieldType::read_from_file) {
-#if defined(WARPX_DIM_RZ)
-        WARPX_ALWAYS_ASSERT_WITH_MESSAGE(n_rz_azimuthal_modes == 1,
-                                         "External field reading is not implemented for more than one RZ mode (see #3829)");
-        ReadExternalFieldFromFile(m_p_ext_field_params->external_fields_path, Efield_fp_external[lev][0].get(), "E", "r");
-        ReadExternalFieldFromFile(m_p_ext_field_params->external_fields_path, Efield_fp_external[lev][1].get(), "E", "t");
-        ReadExternalFieldFromFile(m_p_ext_field_params->external_fields_path, Efield_fp_external[lev][2].get(), "E", "z");
-#else
-        ReadExternalFieldFromFile(m_p_ext_field_params->external_fields_path, Efield_fp_external[lev][0].get(), "E", "x");
-        ReadExternalFieldFromFile(m_p_ext_field_params->external_fields_path, Efield_fp_external[lev][1].get(), "E", "y");
-        ReadExternalFieldFromFile(m_p_ext_field_params->external_fields_path, Efield_fp_external[lev][2].get(), "E", "z");
-#endif
-    }
+    LoadExternalFieldsFromFile(lev);
 
     if (costs[lev]) {
         const auto iarr = costs[lev]->IndexArray();
@@ -1324,6 +1300,37 @@ void WarpX::CheckKnownIssues()
 #endif
 }
 
+void
+WarpX::LoadExternalFieldsFromFile (int const lev)
+{
+    if (m_p_ext_field_params->B_ext_grid_type == ExternalFieldType::read_from_file) {
+#if defined(WARPX_DIM_RZ)
+        WARPX_ALWAYS_ASSERT_WITH_MESSAGE(n_rz_azimuthal_modes == 1,
+                                         "External field reading is not implemented for more than one RZ mode (see #3829)");
+        ReadExternalFieldFromFile(m_p_ext_field_params->external_fields_path, Bfield_fp_external[lev][0].get(), "B", "r");
+        ReadExternalFieldFromFile(m_p_ext_field_params->external_fields_path, Bfield_fp_external[lev][1].get(), "B", "t");
+        ReadExternalFieldFromFile(m_p_ext_field_params->external_fields_path, Bfield_fp_external[lev][2].get(), "B", "z");
+#else
+        ReadExternalFieldFromFile(m_p_ext_field_params->external_fields_path, Bfield_fp_external[lev][0].get(), "B", "x");
+        ReadExternalFieldFromFile(m_p_ext_field_params->external_fields_path, Bfield_fp_external[lev][1].get(), "B", "y");
+        ReadExternalFieldFromFile(m_p_ext_field_params->external_fields_path, Bfield_fp_external[lev][2].get(), "B", "z");
+#endif
+    }
+    if (m_p_ext_field_params->E_ext_grid_type == ExternalFieldType::read_from_file) {
+#if defined(WARPX_DIM_RZ)
+        WARPX_ALWAYS_ASSERT_WITH_MESSAGE(n_rz_azimuthal_modes == 1,
+                                         "External field reading is not implemented for more than one RZ mode (see #3829)");
+        ReadExternalFieldFromFile(m_p_ext_field_params->external_fields_path, Efield_fp_external[lev][0].get(), "E", "r");
+        ReadExternalFieldFromFile(m_p_ext_field_params->external_fields_path, Efield_fp_external[lev][1].get(), "E", "t");
+        ReadExternalFieldFromFile(m_p_ext_field_params->external_fields_path, Efield_fp_external[lev][2].get(), "E", "z");
+#else
+        ReadExternalFieldFromFile(m_p_ext_field_params->external_fields_path, Efield_fp_external[lev][0].get(), "E", "x");
+        ReadExternalFieldFromFile(m_p_ext_field_params->external_fields_path, Efield_fp_external[lev][1].get(), "E", "y");
+        ReadExternalFieldFromFile(m_p_ext_field_params->external_fields_path, Efield_fp_external[lev][2].get(), "E", "z");
+#endif
+    }
+}
+
 #if defined(WARPX_USE_OPENPMD) && !defined(WARPX_DIM_1D_Z) && !defined(WARPX_DIM_XZ)
 void
 WarpX::ReadExternalFieldFromFile (

Source/WarpX.H

@@ -1020,6 +1020,16 @@ public:
          char field,
          int lev, PatchType patch_type);
 
+    /**
+     * \brief Load field values from a user-specified openPMD file,
+     * for the fields Ex, Ey, Ez, Bx, By, Bz
+     */
+    void LoadExternalFieldsFromFile (int lev);
+
+    /**
+     * \brief Load field values from a user-specified openPMD file
+     * for a specific field (specified by `F_name`)
+     */
     void ReadExternalFieldFromFile (
          std::string read_fields_from_path, amrex::MultiFab* mf,
          std::string F_name, std::string F_component);