Add plt part restart (#122)

* Add function to restart particles from plotfiles if the data exists

* Fix bugs from last commit

* Add reset for spray dmap and ba in post regrid

* Fix newer code from development, change input file name to spray restart

* Update PeleMP submodule

* Updated PeleMP submodule

* Remove unnecessary create spray data and resetting dmap for sprays

* Update PeleMP submodule

* Removed repeated x_velocity read

Exec/RegTests/SprayTest/input.2d

@@ -48,7 +48,6 @@ peleLM.sdc_iterMax = 1
 peleLM.floor_species = 0
 
 peleLM.do_temporals = 0
-peleLM.do_mass_balance = 0
 peleLM.num_init_iter = 0
 
 amr.plot_int = 100
@@ -58,6 +57,7 @@ amr.stop_time = 100.
 amr.cfl = 0.5
 amr.plot_speciesState = 1
 amr.plot_grad_p = 0
+amr.derive_plot_vars = mass_fractions
 
 # --------------- INPUTS TO CHEMISTRY REACTOR ---------------
 peleLM.chem_integrator = "ReactorNull"

Source/PeleLM.H

@@ -107,19 +107,19 @@ class PeleLM : public amrex::AmrCore {
 
 #ifdef PELELM_USE_SPRAY
    // ReadSprayParameters
-   void readSprayParameters();
+   void SprayReadParameters();
 
-   // SprayParticleSetup
-   void sprayParticleSetup();
+   // SpraySetup
+   void SpraySetup();
 
-   // CreateSprayData
-   void createSprayData();
+   // SprayCreateData
+   void SprayCreateData();
 
-   // InitSprays
-   void initSprays();
+   // SprayInit
+   void SprayInit();
 
    // SprayRestart
-   void sprayRestart();
+   void SprayRestart(const std::string& restart_file);
 #endif
 
    // VariablesSetup
@@ -446,15 +446,15 @@ class PeleLM : public amrex::AmrCore {
    void removeVirtualParticles(const int level);
    void setupGhostParticles(const int ngrow, const int level);
    void removeGhostParticles(const int level);
-   amrex::Real estSprayDt();
-   void sprayMKD(const amrex::Real time, const amrex::Real dt);
-   void sprayMKDLevel(const int level,
+   amrex::Real SprayEstDt();
+   void SprayMKD(const amrex::Real time, const amrex::Real dt);
+   void SprayMKDLevel(const int level,
                       const amrex::Real time,
                       const amrex::Real dt);
-   void sprayInjectRedist();
-   void sprayPostRegrid();
-   void setSprayState(const amrex::Real& a_flow_dt);
-   void addSpraySource(const int level);
+   void SprayInjectRedist();
+   void SprayPostRegrid();
+   void SpraySetState(const amrex::Real& a_flow_dt);
+   void SprayAddSource(const int level);
    static int write_spray_ascii_files;
    static int num_spray_src;
    static bool do_spray_particles;

Source/PeleLMAdvance.cpp

@@ -41,9 +41,9 @@ void PeleLM::Advance(int is_initIter) {
    // Compute time-step size
    m_dt = computeDt(is_initIter,AmrOldTime);
 #ifdef PELELM_USE_SPRAY
-   if (!is_initIter && do_spray_particles) {
+   if (!is_initIter) {
      // Create the state MF used for spray interpolation
-     setSprayState(m_dt);
+     SpraySetState(m_dt);
    }
 #endif
 
@@ -110,7 +110,7 @@ void PeleLM::Advance(int is_initIter) {
 
 #ifdef PELELM_USE_SPRAY
    if (!is_initIter) {
-     sprayMKD(m_cur_time, m_dt);
+     SprayMKD(m_cur_time, m_dt);
    }
 #endif
 #ifdef PELELM_USE_SOOT

Source/PeleLMEvolve.cpp

@@ -32,7 +32,7 @@ void PeleLM::Evolve() {
 #ifdef PELELM_USE_SPRAY
       // Inject and redistribute spray particles
       if (do_spray_particles && regridded) {
-        sprayPostRegrid();
+        SprayPostRegrid();
       }
 #endif
       int is_init = 0;
@@ -42,7 +42,7 @@ void PeleLM::Evolve() {
 
 #ifdef PELELM_USE_SPRAY
       if (do_spray_particles) {
-        sprayInjectRedist();
+        SprayInjectRedist();
       }
 #endif
 

Source/PeleLMInit.cpp

@@ -160,9 +160,7 @@ void PeleLM::initData() {
       updateDiagnostics();
 
 #ifdef PELELM_USE_SPRAY
-      if (do_spray_particles) {
-        initSprays();
-      }
+      SprayInit();
 #endif
 
       //----------------------------------------------------------------
@@ -311,9 +309,7 @@ void PeleLM::initData() {
       ReadCheckPointFile();
 
 #ifdef PELELM_USE_SPRAY
-      if (do_spray_particles) {
-        sprayRestart();
-      }
+      SprayInit();
 #endif
 #ifdef PELE_USE_EFIELD
       // If restarting from a non efield simulation

Source/PeleLMPlot.cpp

@@ -317,7 +317,7 @@ void PeleLM::WritePlotFile() {
    }
 
 #ifdef PELELM_USE_SPRAY
-   if (theSprayPC() != nullptr && do_spray_particles) {
+   if (do_spray_particles) {
      bool is_spraycheck = false;
      for (int lev = 0; lev <= finest_level; ++lev) {
        theSprayPC()->SprayParticleIO(
@@ -430,7 +430,7 @@ void PeleLM::WriteCheckPointFile()
       }
    }
 #ifdef PELELM_USE_SPRAY
-   if (theSprayPC() != nullptr && do_spray_particles) {
+   if (do_spray_particles) {
      int write_ascii = 0; // Not for checkpoints
      bool is_spraycheck = true;
      for (int lev = 0; lev <= finest_level; ++lev) {
@@ -665,7 +665,6 @@ void PeleLM::initLevelDataFromPlt(int a_lev,
         if (plt_vars[i] == "Temp")            idT = i; 
       }
       if (plt_vars[i] == "x_velocity")      idV = i; 
-      if (plt_vars[i] == "x_velocity")      idV = i;
       if (firstChars == "Y(" && idY < 0 ) {  // species might not be ordered in the order of the current mech.
          idY = i;
       }

Source/PeleLMSetup.cpp

@@ -64,9 +64,7 @@ void PeleLM::Setup() {
    taggingSetup();
 
 #ifdef PELELM_USE_SPRAY
-   if (do_spray_particles) {
-     sprayParticleSetup();
-   }
+   SpraySetup();
 #endif
 #ifdef PELELM_USE_SOOT
    if (do_soot_solve) {
@@ -477,7 +475,7 @@ void PeleLM::readParameters() {
 #endif
 
 #ifdef PELELM_USE_SPRAY
-   readSprayParameters();
+   SprayReadParameters();
 #endif
 #ifdef PELELM_USE_SOOT
    do_soot_solve = true;

Source/PeleLMTimestep.cpp

@@ -35,7 +35,7 @@ PeleLM::computeDt(int is_init,
          estdt = std::min(estdt, dtions);
 #endif
 #ifdef PELELM_USE_SPRAY
-         Real dtspray = estSprayDt();
+         Real dtspray = SprayEstDt();
          estdt = std::min(estdt, dtspray);
 #endif
          if ( m_verbose ) {

Source/Spray/PeleLMSprayParticles.cpp

@@ -81,13 +81,13 @@ PeleLM::theGhostPC()
 }
 
 Real
-PeleLM::estSprayDt()
+PeleLM::SprayEstDt()
 {
   Real estdt = 1.0e200;
   if (!do_spray_particles || theSprayPC() == nullptr) {
     return estdt;
   }
-  BL_PROFILE("PeleLM::sprayEstTimeStep()");
+  BL_PROFILE("PeleLM::SprayEstDt()");
   for (int lev = 0; lev <= finest_level; ++lev) {
     Real estdt_lev = theSprayPC()->estTimestep(lev, max_spray_cfl);
     if (estdt_lev > 0. && estdt_lev < estdt) {
@@ -101,7 +101,7 @@ PeleLM::estSprayDt()
 }
 
 void
-PeleLM::readSprayParameters()
+PeleLM::SprayReadParameters()
 {
   ParmParse pp("peleLM");
 
@@ -119,8 +119,11 @@ PeleLM::readSprayParameters()
 }
 
 void
-PeleLM::sprayParticleSetup()
+PeleLM::SpraySetup()
 {
+  if (!do_spray_particles) {
+    return;
+  }
   // There must be at least as many fuel species in the spray as
   // there are species in the fluid
   if (SPRAY_FUEL_NUM > NUM_SPECIES) {
@@ -188,7 +191,6 @@ PeleLM::setupVirtualParticles(const int level)
 void
 PeleLM::removeVirtualParticles(const int level)
 {
-  BL_PROFILE("PeleLM::removeVirtualParticles()");
   if (theVirtPC() != nullptr) {
     theVirtPC()->RemoveParticlesAtLevel(level);
   }
@@ -209,14 +211,13 @@ PeleLM::setupGhostParticles(const int ngrow, const int level)
 void
 PeleLM::removeGhostParticles(const int level)
 {
-  BL_PROFILE("PeleLM::removeGhostParticles()");
   if (GhostPC != nullptr) {
     GhostPC->RemoveParticlesAtLevel(level);
   }
 }
 
 void
-PeleLM::createSprayData()
+PeleLM::SprayCreateData()
 {
   SprayPC = new SprayParticleContainer(
     this, &m_phys_bc, sprayData, scomps, wall_temp, max_spray_cfl);
@@ -228,18 +229,30 @@ PeleLM::createSprayData()
 }
 
 void
-PeleLM::initSprays()
+PeleLM::SprayInit()
 {
-  BL_PROFILE("PeleLM::initSprays()");
+  BL_PROFILE("PeleLM::SprayInit()");
 
   //
   // Make sure to call RemoveParticlesOnExit() on exit.
   //
   amrex::ExecOnFinalize(RemoveParticlesOnExit);
 
-  if (do_spray_particles) {
-    if (theSprayPC() == nullptr) {
-      createSprayData();
+  if (!do_spray_particles) {
+    return;
+  }
+  SprayCreateData();
+  // Check if we are restarting and if we restart file has particles
+  std::string restart_file =
+    (m_restart_chkfile.empty()) ? m_restart_pltfile : m_restart_chkfile;
+  std::string restart_partfile = restart_file + "/particles";
+  if (FileSystem::Exists(restart_partfile)) {
+    SprayRestart(restart_file);
+  } else {
+    if (!m_restart_chkfile.empty() || !m_restart_pltfile.empty()) {
+      std::string warn_msg = "Restart file does not contain particles. "
+                             "Particles are being initialized from scratch.";
+      amrex::Warning(warn_msg);
     }
 
     bool init_part = true;
@@ -251,8 +264,8 @@ PeleLM::initSprays()
     }
     ProbParm const* lprobparm = prob_parm;
     theSprayPC()->InitSprayParticles(init_part, *lprobparm);
-    sprayPostRegrid();
-    sprayInjectRedist();
+    SprayPostRegrid();
+    SprayInjectRedist();
     if (spray_verbose >= 1) {
       amrex::Print() << "Total number of initial particles "
                      << theSprayPC()->TotalNumberOfParticles(false, false)
@@ -262,21 +275,18 @@ PeleLM::initSprays()
 }
 
 void
-PeleLM::sprayRestart()
+PeleLM::SprayRestart(const std::string& restart_file)
 {
   if (do_spray_particles) {
-    AMREX_ASSERT(SprayPC == nullptr);
-    createSprayData();
-
     //
     // Make sure to call RemoveParticlesOnExit() on exit.
     //
     amrex::ExecOnFinalize(RemoveParticlesOnExit);
     {
-      theSprayPC()->Restart(m_restart_chkfile, "particles", true);
+      theSprayPC()->Restart(restart_file, "particles");
       ProbParm const* lprobparm = prob_parm;
       theSprayPC()->InitSprayParticles(false, *lprobparm);
-      sprayPostRegrid();
+      SprayPostRegrid();
       amrex::Gpu::Device::streamSynchronize();
     }
   }
@@ -285,8 +295,11 @@ PeleLM::sprayRestart()
 // Sets the number of ghost cells for the state, source, and ghost particles
 // Also creates and fills the state data
 void
-PeleLM::setSprayState(const Real& a_flow_dt)
+PeleLM::SpraySetState(const Real& a_flow_dt)
 {
+  if (!do_spray_particles) {
+    return;
+  }
   spray_cfl.resize(finest_level + 1);
   // Number of ghost cells needed to interpolate state to spray position
   spray_state_ghosts.resize(finest_level + 1);
@@ -344,7 +357,7 @@ PeleLM::setSprayState(const Real& a_flow_dt)
 }
 
 void
-PeleLM::addSpraySource(const int level)
+PeleLM::SprayAddSource(const int level)
 {
   MultiFab& source = *(m_spraysource[level].get());
   MultiFab& extsource = *(m_extSource[level].get());
@@ -364,7 +377,7 @@ PeleLM::addSpraySource(const int level)
 }
 
 void
-PeleLM::sprayMKD(const Real time, const Real dt)
+PeleLM::SprayMKD(const Real time, const Real dt)
 {
   if (!do_spray_particles) {
     return;
@@ -373,22 +386,23 @@ PeleLM::sprayMKD(const Real time, const Real dt)
     amrex::Print()
       << "moveKickDrift ... updating particle positions and velocity\n";
   }
+  BL_PROFILE("PeleLM::SprayMKD()");
   // Setup the virtual particles that represent particles on finer levels
   setupVirtualParticles(0);
   for (int lev = 0; lev <= finest_level; ++lev) {
     if (spray_verbose > 1) {
-      amrex::Print() << "sprayMKDLevel " << lev << std::endl;
+      amrex::Print() << "SprayMKDLevel " << lev << std::endl;
     }
-    sprayMKDLevel(lev, time, dt);
-    addSpraySource(lev);
+    SprayMKDLevel(lev, time, dt);
+    SprayAddSource(lev);
     removeGhostParticles(lev);
     removeVirtualParticles(lev);
     m_spraysource[lev]->setVal(0.);
   }
 }
 
 void
-PeleLM::sprayMKDLevel(const int level, const Real time, const Real dt)
+PeleLM::SprayMKDLevel(const int level, const Real time, const Real dt)
 {
   if (level < finest_level) {
     // Make a copy of the particles on this level into ghost particles
@@ -430,7 +444,7 @@ PeleLM::sprayMKDLevel(const int level, const Real time, const Real dt)
 }
 
 void
-PeleLM::sprayPostRegrid()
+PeleLM::SprayPostRegrid()
 {
   static Vector<BoxArray> ba_spray;
   static Vector<DistributionMapping> dm_spray;
@@ -465,9 +479,9 @@ PeleLM::sprayPostRegrid()
 }
 
 void
-PeleLM::sprayInjectRedist()
+PeleLM::SprayInjectRedist()
 {
-  BL_PROFILE("PeleLM::sprayInjectRedist");
+  BL_PROFILE("PeleLM::SprayInjectRedist");
   Long prev_count = 0;
   if (spray_verbose >= 3) {
     prev_count = theSprayPC()->TotalNumberOfParticles(true, false);