Separated the calculation of new time-steps from the second kick

src/cell.c

@@ -974,6 +974,7 @@ int cell_unskip_tasks(struct cell *c, struct scheduler *s) {
   if (c->drift != NULL) scheduler_activate(s, c->drift);
   if (c->kick1 != NULL) scheduler_activate(s, c->kick1);
   if (c->kick2 != NULL) scheduler_activate(s, c->kick2);
+  if (c->timestep != NULL) scheduler_activate(s, c->timestep);
   if (c->cooling != NULL) scheduler_activate(s, c->cooling);
   if (c->sourceterms != NULL) scheduler_activate(s, c->sourceterms);
 

src/cell.h

@@ -156,6 +156,9 @@ struct cell {
   /*! The second kick task */
   struct task *kick2;
 
+  /*! The task to compute time-steps */
+  struct task *timestep;
+
   /*! Task constructing the multipole from the particles */
   struct task *grav_up;
 

src/engine.c

@@ -147,6 +147,12 @@ void engine_make_hierarchical_tasks(struct engine *e, struct cell *c) {
       c->kick2 = scheduler_addtask(s, task_type_kick2, task_subtype_none, 0, 0,
                                    c, NULL, 0);
 
+      /* Add the time-step calculation task and its dependency */
+      c->timestep = scheduler_addtask(s, task_type_timestep, task_subtype_none,
+                                      0, 0, c, NULL, 0);
+
+      scheduler_addunlock(s, c->kick2, c->timestep);
+
       /* Add the drift task and its dependencies. */
       c->drift = scheduler_addtask(s, task_type_drift, task_subtype_none, 0, 0,
                                    c, NULL, 0);
@@ -2131,9 +2137,12 @@ void engine_marktasks_mapper(void *map_data, int num_elements,
     /* Kick? */
     else if (t->type == task_type_kick1) {
       if (t->ci->ti_end_min <= ti_end) scheduler_activate(s, t);
+    } else if (t->type == task_type_kick2) {
+      if (t->ci->ti_end_min <= ti_end) scheduler_activate(s, t);
     }
 
-    else if (t->type == task_type_kick2) {
+    /* Time-step? */
+    else if (t->type == task_type_timestep) {
       t->ci->updated = 0;
       t->ci->g_updated = 0;
       if (t->ci->ti_end_min <= ti_end) scheduler_activate(s, t);
@@ -2365,7 +2374,7 @@ void engine_barrier(struct engine *e, int tid) {
 void engine_collect_kick(struct cell *c) {
 
   /* Skip super-cells (Their values are already set) */
-  if (c->kick2 != NULL) return;
+  if (c->timestep != NULL) return;
 
   /* Counters for the different quantities. */
   int updated = 0, g_updated = 0;

src/runner.c

@@ -896,8 +896,7 @@ void runner_do_kick1(struct runner *r, struct cell *c, int timer) {
 /**
  * @brief Perform the second half-kick on all the active particles in a cell.
  *
- * Also computes the next time-step of all active particles, prepare them to be
- * drifted and update the cell's statistics.
+ * Also prepares particles to be drifted.
  *
  * @param r The runner thread.
  * @param c The cell.
@@ -918,17 +917,13 @@ void runner_do_kick2(struct runner *r, struct cell *c, int timer) {
   TIMER_TIC;
 
   /* Anything to do here? */
-  if (!cell_is_active(c, e)) {
-    c->updated = 0;
-    c->g_updated = 0;
-    return;
-  }
-
-  int updated = 0, g_updated = 0;
-  integertime_t ti_end_min = max_nr_timesteps, ti_end_max = 0;
+  if (!cell_is_active(c, e)) return;
 
-  /* No children? */
-  if (!c->split) {
+  /* Recurse? */
+  if (c->split) {
+    for (int k = 0; k < 8; k++)
+      if (c->progeny[k] != NULL) runner_do_kick2(r, c->progeny[k], 0);
+  } else {
 
     /* Loop over the particles in this cell. */
     for (int k = 0; k < count; k++) {
@@ -947,8 +942,6 @@ void runner_do_kick2(struct runner *r, struct cell *c, int timer) {
         const integertime_t ti_step = get_integer_timestep(p->time_bin);
         const integertime_t ti_begin =
             get_integer_time_begin(ti_current, p->time_bin);
-        const integertime_t ti_end =
-            get_integer_time_end(ti_current, p->time_bin);
 
         /* Finish the time-step with a second half-kick */
         kick_part(p, xp, ti_begin + ti_step / 2, ti_begin + ti_step, ti_current,
@@ -956,6 +949,82 @@ void runner_do_kick2(struct runner *r, struct cell *c, int timer) {
 
         /* Prepare the values to be drifted */
         hydro_reset_predicted_values(p, xp);
+      }
+    }
+
+    /* Loop over the g-particles in this cell. */
+    for (int k = 0; k < gcount; k++) {
+
+      /* Get a handle on the part. */
+      struct gpart *restrict gp = &gparts[k];
+
+      /* If the g-particle has no counterpart */
+      if (gp->id_or_neg_offset > 0) {
+
+        /* need to be kicked ? */
+        if (gpart_is_active(gp, e)) {
+
+          /* First, finish the force loop */
+          gravity_end_force(gp, const_G);
+
+          const integertime_t ti_step = get_integer_timestep(gp->time_bin);
+          const integertime_t ti_begin =
+              get_integer_time_begin(ti_current, gp->time_bin);
+
+          /* Finish the time-step with a second half-kick */
+          kick_gpart(gp, ti_begin + ti_step / 2, ti_begin + ti_step, ti_current,
+                     timeBase);
+        }
+      }
+    }
+  }
+  if (timer) TIMER_TOC(timer_kick2);
+}
+
+/**
+ * @brief Computes the next time-step of all active particles in this cell
+ * and update the cell's statistics.
+ *
+ * @param r The runner thread.
+ * @param c The cell.
+ * @param timer Are we timing this ?
+ */
+void runner_do_timestep(struct runner *r, struct cell *c, int timer) {
+
+  const struct engine *e = r->e;
+  const integertime_t ti_current = e->ti_current;
+  const int count = c->count;
+  const int gcount = c->gcount;
+  struct part *restrict parts = c->parts;
+  struct xpart *restrict xparts = c->xparts;
+  struct gpart *restrict gparts = c->gparts;
+
+  TIMER_TIC;
+
+  int updated = 0, g_updated = 0;
+  integertime_t ti_end_min = max_nr_timesteps, ti_end_max = 0;
+
+  /* No children? */
+  if (!c->split) {
+
+    /* Loop over the particles in this cell. */
+    for (int k = 0; k < count; k++) {
+
+      /* Get a handle on the part. */
+      struct part *restrict p = &parts[k];
+      struct xpart *restrict xp = &xparts[k];
+
+      /* If particle needs updating */
+      if (part_is_active(p, e)) {
+
+#ifdef SWIFT_DEBUG_CHECKS
+        /* Current end of time-step */
+        const integertime_t ti_end =
+            get_integer_time_end(ti_current, p->time_bin);
+
+        if (ti_end != ti_current)
+          error("Computing time-step of rogue particle.");
+#endif
 
         /* Get new time-step */
         const integertime_t ti_new_step = get_part_timestep(p, xp, e);
@@ -969,8 +1038,8 @@ void runner_do_kick2(struct runner *r, struct cell *c, int timer) {
         if (p->gpart != NULL) g_updated++;
 
         /* What is the next sync-point ? */
-        ti_end_min = min(ti_end + ti_new_step, ti_end_min);
-        ti_end_max = max(ti_end + ti_new_step, ti_end_max);
+        ti_end_min = min(ti_current + ti_new_step, ti_end_min);
+        ti_end_max = max(ti_current + ti_new_step, ti_end_max);
       }
 
       else { /* part is inactive */
@@ -993,21 +1062,17 @@ void runner_do_kick2(struct runner *r, struct cell *c, int timer) {
       /* If the g-particle has no counterpart */
       if (gp->id_or_neg_offset > 0) {
 
-        /* need to be kicked ? */
+        /* need to be updated ? */
         if (gpart_is_active(gp, e)) {
 
-          /* First, finish the force loop */
-          gravity_end_force(gp, const_G);
-
-          const integertime_t ti_step = get_integer_timestep(gp->time_bin);
-          const integertime_t ti_begin =
-              get_integer_time_begin(ti_current, gp->time_bin);
+#ifdef SWIFT_DEBUG_CHECKS
+          /* Current end of time-step */
           const integertime_t ti_end =
               get_integer_time_end(ti_current, gp->time_bin);
 
-          /* Finish the time-step with a second half-kick */
-          kick_gpart(gp, ti_begin + ti_step / 2, ti_begin + ti_step, ti_current,
-                     timeBase);
+          if (ti_end != ti_current)
+            error("Computing time-step of rogue particle.");
+#endif
 
           /* Get new time-step */
           const integertime_t ti_new_step = get_gpart_timestep(gp, e);
@@ -1019,8 +1084,8 @@ void runner_do_kick2(struct runner *r, struct cell *c, int timer) {
           g_updated++;
 
           /* What is the next sync-point ? */
-          ti_end_min = min(ti_end + ti_new_step, ti_end_min);
-          ti_end_max = max(ti_end + ti_new_step, ti_end_max);
+          ti_end_min = min(ti_current + ti_new_step, ti_end_min);
+          ti_end_max = max(ti_current + ti_new_step, ti_end_max);
         } else { /* gpart is inactive */
 
           const integertime_t ti_end =
@@ -1032,10 +1097,7 @@ void runner_do_kick2(struct runner *r, struct cell *c, int timer) {
         }
       }
     }
-  }
-
-  /* Otherwise, aggregate data from children. */
-  else {
+  } else {
 
     /* Loop over the progeny. */
     for (int k = 0; k < 8; k++)
@@ -1043,7 +1105,7 @@ void runner_do_kick2(struct runner *r, struct cell *c, int timer) {
         struct cell *restrict cp = c->progeny[k];
 
         /* Recurse */
-        runner_do_kick2(r, cp, 0);
+        runner_do_timestep(r, cp, 0);
 
         /* And aggregate */
         updated += cp->updated;
@@ -1059,7 +1121,7 @@ void runner_do_kick2(struct runner *r, struct cell *c, int timer) {
   c->ti_end_min = ti_end_min;
   c->ti_end_max = ti_end_max;
 
-  if (timer) TIMER_TOC(timer_kick2);
+  if (timer) TIMER_TOC(timer_timestep);
 }
 
 /**
@@ -1297,6 +1359,9 @@ void *runner_main(void *data) {
         case task_type_kick2:
           runner_do_kick2(r, ci, 1);
           break;
+        case task_type_timestep:
+          runner_do_timestep(r, ci, 1);
+          break;
 #ifdef WITH_MPI
         case task_type_send:
           if (t->subtype == task_subtype_tend) {

src/space.c

@@ -215,6 +215,7 @@ void space_rebuild_recycle_mapper(void *map_data, int num_elements,
     c->ghost = NULL;
     c->kick1 = NULL;
     c->kick2 = NULL;
+    c->timestep = NULL;
     c->drift = NULL;
     c->cooling = NULL;
     c->sourceterms = NULL;

src/task.c

@@ -48,10 +48,11 @@
 
 /* Task type names. */
 const char *taskID_names[task_type_count] = {
-    "none",     "sort",    "self",    "pair",        "sub_self",
-    "sub_pair", "init",    "ghost",   "extra_ghost", "drift",
-    "kick1",    "kick2",   "send",    "recv",        "grav_gather_m",
-    "grav_fft", "grav_mm", "grav_up", "cooling",     "sourceterms"};
+    "none",          "sort",     "self",     "pair",        "sub_self",
+    "sub_pair",      "init",     "ghost",    "extra_ghost", "drift",
+    "kick1",         "kick2",    "timestep", "send",        "recv",
+    "grav_gather_m", "grav_fft", "grav_mm",  "grav_up",     "cooling",
+    "sourceterms"};
 
 const char *subtaskID_names[task_subtype_count] = {
     "none", "density", "gradient", "force", "grav", "external_grav", "tend"};
@@ -149,6 +150,7 @@ __attribute__((always_inline)) INLINE static enum task_actions task_acts_on(
     case task_type_init:
     case task_type_kick1:
     case task_type_kick2:
+    case task_type_timestep:
     case task_type_send:
     case task_type_recv:
     case task_type_drift:

src/task.h

@@ -48,6 +48,7 @@ enum task_types {
   task_type_drift,
   task_type_kick1,
   task_type_kick2,
+  task_type_timestep,
   task_type_send,
   task_type_recv,
   task_type_grav_gather_m,

src/timers.h

@@ -35,6 +35,7 @@ enum {
   timer_drift,
   timer_kick1,
   timer_kick2,
+  timer_timestep,
   timer_dosort,
   timer_doself_density,
   timer_doself_gradient,