Merge branch 'multi_time_step_warnings' into 'master'

Multi time step warnings and stop condition

Solves the issues #111 flagged by @alepper 

We now explicitly check the time-step sizes given in input are sensible.

We now have a function to test whether we have reached the end of the simulation or not. 

Branch can be removed.

See merge request !104

examples/main.c

@@ -171,7 +171,7 @@ int main(int argc, char *argv[]) {
       case 'd':
         if (sscanf(optarg, "%f", &dt_min) != 1)
           error("Error parsing minimal timestep.");
-        if (myrank == 0) message("dt_min set to %e.", dt_max);
+        if (myrank == 0) message("dt_min set to %e.", dt_min);
         fflush(stdout);
         break;
       case 'e':
@@ -192,8 +192,8 @@ int main(int argc, char *argv[]) {
         with_outputs = 0;
         break;
       case 'P':
-        /* Partition type is one of "g", "m", "w", or "v"; "g" can be
-         * followed by three numbers defining the grid. */
+/* Partition type is one of "g", "m", "w", or "v"; "g" can be
+ * followed by three numbers defining the grid. */
 #ifdef WITH_MPI
         switch (optarg[0]) {
           case 'g':
@@ -224,8 +224,8 @@ int main(int argc, char *argv[]) {
           error("Error parsing number of queues.");
         break;
       case 'R':
-        /* Repartition type "n", "b", "v", "e" or "x". 
-         * Note only none is available without METIS. */
+/* Repartition type "n", "b", "v", "e" or "x".
+ * Note only none is available without METIS. */
 #ifdef WITH_MPI
         switch (optarg[0]) {
           case 'n':
@@ -323,10 +323,6 @@ int main(int argc, char *argv[]) {
             aFactor(&us, UNIT_CONV_ENTROPY), hFactor(&us, UNIT_CONV_ENTROPY));
   }
 
-  /* Check we have sensible time step bounds */
-  if (dt_min > dt_max)
-    error("Minimal time step size must be large than maximal time step size ");
-
   /* Check whether an IC file has been provided */
   if (strcmp(ICfileName, "") == 0)
     error("An IC file name must be provided via the option -f");
@@ -468,7 +464,7 @@ int main(int argc, char *argv[]) {
         "[ms]\n");
 
   /* Let loose a runner on the space. */
-  for (j = 0; e.time < time_end; j++) {
+  for (j = 0; !engine_is_done(&e); j++) {
 
 /* Repartition the space amongst the nodes? */
 #ifdef WITH_MPI

src/engine.c

@@ -285,7 +285,6 @@ void engine_redistribute(struct engine *e) {
 #endif
 }
 
-
 /**
  * @brief Repartition the cells amongst the nodes.
  *
@@ -1393,7 +1392,7 @@ void engine_init_particles(struct engine *e) {
 
   struct space *s = e->s;
 
-  if(e->nodeID == 0) message("Initialising particles");
+  if (e->nodeID == 0) message("Initialising particles");
 
   /* Make sure all particles are ready to go */
   /* i.e. clean-up any stupid state in the ICs */
@@ -1625,6 +1624,13 @@ void engine_step(struct engine *e) {
   // printParticle(e->s->parts, e->s->xparts,515050, e->s->nr_parts);
 }
 
+/**
+ * @brief Returns 1 if the simulation has reached its end point, 0 otherwise
+ */
+int engine_is_done(struct engine *e) {
+  return !(e->ti_current < max_nr_timesteps);
+}
+
 /**
  * @brief Create and fill the proxies.
  *
@@ -1831,7 +1837,7 @@ void engine_init(struct engine *e, struct space *s, float dt, int nr_threads,
       int home = numa_node_of_cpu(sched_getcpu()), half = nr_cores / 2;
       bool done = false, swap_hyperthreads = hyperthreads_present();
       if (swap_hyperthreads && nodeID == 0)
-	message("prefer physical cores to hyperthreads");
+        message("prefer physical cores to hyperthreads");
 
       while (!done) {
         done = true;
@@ -1928,15 +1934,25 @@ void engine_init(struct engine *e, struct space *s, float dt, int nr_threads,
   engine_print_policy(e);
 
   /* Print information about the hydro scheme */
-  if (e->nodeID == 0)
-    message("Hydrodynamic scheme: %s", SPH_IMPLEMENTATION);
+  if (e->nodeID == 0) message("Hydrodynamic scheme: %s", SPH_IMPLEMENTATION);
+
+  /* Check we have sensible time bounds */
+  if (timeBegin >= timeEnd)
+    error(
+        "Final simulation time (t_end = %e) must be larger than the start time "
+        "(t_beg = %e)",
+        timeEnd, timeBegin);
+
+  /* Check we have sensible time step bounds */
+  if (e->dt_min > e->dt_max)
+    error(
+        "Minimal time step size must be smaller than maximal time step size ");
 
   /* Deal with timestep */
   e->timeBase = (timeEnd - timeBegin) / max_nr_timesteps;
   e->ti_current = 0;
-  if (e->nodeID == 0)
-    message("Absolute minimal timestep size: %e", e->timeBase);
 
+  /* Fixed time-step case */
   if ((e->policy & engine_policy_fixdt) == engine_policy_fixdt) {
     e->dt_min = e->dt_max;
 
@@ -1947,11 +1963,32 @@ void engine_init(struct engine *e, struct space *s, float dt, int nr_threads,
     e->dt_min = e->dt_max = dti_timeline * e->timeBase;
 
     if (e->nodeID == 0) message("Timestep set to %e", e->dt_max);
+  } else {
+
+    if (e->nodeID == 0) {
+      message("Absolute minimal timestep size: %e", e->timeBase);
+
+      float dt_min = timeEnd - timeBegin;
+      while (dt_min > e->dt_min) dt_min /= 2.f;
+
+      message("Minimal timestep size (on time-line): %e", dt_min);
+
+      float dt_max = timeEnd - timeBegin;
+      while (dt_max > e->dt_max) dt_max /= 2.f;
+
+      message("Maximal timestep size (on time-line): %e", dt_max);
+    }
   }
 
   if (e->dt_min < e->timeBase && e->nodeID == 0)
-    error("Minimal timestep smaller than the absolue possible minimum dt=%e",
-          e->timeBase);
+    error(
+        "Minimal time-step size smaller than the absolute possible minimum "
+        "dt=%e",
+        e->timeBase);
+
+  if (e->dt_max > (e->timeEnd - e->timeBegin) && e->nodeID == 0)
+    error("Maximal time-step size larger than the simulation run time t=%e",
+          e->timeEnd - e->timeBegin);
 
 /* Construct types for MPI communications */
 #ifdef WITH_MPI

src/engine.h

@@ -63,7 +63,6 @@ extern const char *engine_policy_names[];
 #define engine_maxproxies 64
 #define engine_tasksreweight 10
 
-
 /* The rank of the engine as a global variable (for messages). */
 extern int engine_rank;
 
@@ -186,5 +185,6 @@ void engine_makeproxies(struct engine *e);
 void engine_redistribute(struct engine *e);
 struct link *engine_addlink(struct engine *e, struct link *l, struct task *t);
 void engine_print_policy(struct engine *e);
+int engine_is_done(struct engine *e);
 
 #endif /* SWIFT_ENGINE_H */