Add more global time step control

src/driver/driver.cpp

@@ -205,25 +205,60 @@ void EvolutionDriver::InitializeBlockTimeSteps() {
 // \brief function that loops over all MeshBlocks and find new timestep
 
 void EvolutionDriver::SetGlobalTimeStep() {
-  // don't allow dt to grow by more than 2x
-  // consider making this configurable in the input
-  if (tm.dt < 0.1 * std::numeric_limits<Real>::max()) {
-    tm.dt *= 2.0;
-  }
-  Real big = std::numeric_limits<Real>::max();
-  for (auto const &pmb : pmesh->block_list) {
-    tm.dt = std::min(tm.dt, pmb->NewDt());
-    pmb->SetAllowedDt(big);
+  // Check if user wants to force the value
+  if (dt_force > 0.0) {
+    tm.dt = dt_force;
+  } else {
+    // Check for special first cycle value
+    if (tm.ncycle == 0) {
+      tm.dt = std::min(tm.dt, dt_init);
+    }
+    // don't allow dt to grow by more than 2x
+    // consider making this configurable in the input
+    if (tm.dt < 0.1 * std::numeric_limits<Real>::max()) {
+      tm.dt *= 2.0;
+    }
+    // Allow the meshblocks to vote
+    Real big = std::numeric_limits<Real>::max();
+    for (auto const &pmb : pmesh->block_list) {
+      tm.dt = std::min(tm.dt, pmb->NewDt());
+      pmb->SetAllowedDt(big);
+    }
+    // Allow the user to vote
+    tm.dt = std::min(tm.dt, dt_user);
   }
 
 #ifdef MPI_PARALLEL
   PARTHENON_MPI_CHECK(MPI_Allreduce(MPI_IN_PLACE, &tm.dt, 1, MPI_PARTHENON_REAL, MPI_MIN,
                                     MPI_COMM_WORLD));
 #endif
 
-  if (tm.time < tm.tlim &&
-      (tm.tlim - tm.time) < tm.dt) // timestep would take us past desired endpoint
+  // Check that we have not gone off the rails
+  if (tm.dt <= dt_min) {
+    if (++dt_min_count >= dt_min_count_max) {
+      std::stringstream msg;
+      msg << "Timesetep has fallen bellow minimum (parthenon/time/dt_min=" << dt_min
+          << ") for more than " << dt_min_count_max << " steps";
+      PARTHENON_FAIL(msg);
+    }
+  } else {
+    dt_min_count = 0;
+  }
+  if (tm.dt >= dt_max) {
+    if (++dt_max_count >= dt_max_count_max) {
+      std::stringstream msg;
+      msg << "Timesetep has risen above maximum (parthenon/time/dt_max=" << dt_max
+          << ") for more than " << dt_max_count_max << " steps";
+      PARTHENON_FAIL(msg);
+    }
+  } else {
+    dt_max_count = 0;
+  }
+
+  // Limit timestep if it would take us past desired endpoint
+  if (tm.time < tm.tlim && (tm.tlim - tm.time) < tm.dt) {
     tm.dt = tm.tlim - tm.time;
+  }
 }
 
 void EvolutionDriver::DumpInputParameters() {

src/driver/driver.hpp

@@ -69,6 +69,21 @@ class EvolutionDriver : public Driver {
                                       std::numeric_limits<Real>::infinity());
     Real dt =
         pinput->GetOrAddReal("parthenon/time", "dt", std::numeric_limits<Real>::max());
+    dt_min = pinput->GetOrAddReal("parthenon/time", "dt_min",
+                                  std::numeric_limits<Real>::min());
+    dt_max = pinput->GetOrAddReal("parthenon/time", "dt_max",
+                                  std::numeric_limits<Real>::max());
+    dt_init = pinput->GetOrAddReal("parthenon/time", "dt_init",
+                                   std::numeric_limits<Real>::max());
+    dt_user = pinput->GetOrAddReal("parthenon/time", "dt_user",
+                                   std::numeric_limits<Real>::max());
+    dt_force = pinput->GetOrAddReal("parthenon/time", "dt_force",
+                                    std::numeric_limits<Real>::lowest());
+    dt_min_count_max =
+        pinput->GetOrAddInteger("parthenon/time", "dt_min_cycle_limit", 10);
+    dt_max_count_max = pinput->GetOrAddInteger("parthenon/time", "dt_max_cycle_limit", 1);
+    dt_min_count = 0;
+    dt_max_count = 0;
     const auto ncycle = pinput->GetOrAddInteger("parthenon/time", "ncycle", 0);
     const auto nmax = pinput->GetOrAddInteger("parthenon/time", "nlim", -1);
     const auto nout = pinput->GetOrAddInteger("parthenon/time", "ncycle_out", 1);
@@ -88,6 +103,9 @@ class EvolutionDriver : public Driver {
 
  protected:
   void PostExecute(DriverStatus status) override;
+  Real dt_user, dt_force, dt_init, dt_min, dt_max;
+  int dt_min_count, dt_max_count;
+  int dt_min_count_max, dt_max_count_max;
 
  private:
   void InitializeBlockTimeSteps();