Construct a list of local top-level cells that are then used in the unskipping

examples/main.c

@@ -544,8 +544,8 @@ int main(int argc, char *argv[]) {
   /* Initialize the space with these data. */
   if (myrank == 0) clocks_gettime(&tic);
   struct space s;
-  space_init(&s, params, dim, parts, gparts, sparts, Ngas, Ngpart, Nspart,
-             periodic, replicate, with_self_gravity, talking, dry_run);
+  space_init(&s, myrank, params, dim, parts, gparts, sparts, Ngas, Ngpart,
+             Nspart, periodic, replicate, with_self_gravity, talking, dry_run);
   if (myrank == 0) {
     clocks_gettime(&toc);
     message("space_init took %.3f %s.", clocks_diff(&tic, &toc),

src/engine.c

@@ -3807,8 +3807,8 @@ void engine_unskip(struct engine *e) {
   const ticks tic = getticks();
 
   /* Activate all the regular tasks */
-  threadpool_map(&e->threadpool, runner_do_unskip_mapper, e->s->cells_top,
-                 e->s->nr_cells, sizeof(struct cell), 1, e);
+  threadpool_map(&e->threadpool, runner_do_unskip_mapper, e->s->local_cells_top,
+                 e->s->nr_local_cells, sizeof(int), 1, e);
 
   /* And the top level gravity FFT one */
   if (e->s->periodic && (e->policy & engine_policy_self_gravity))

src/runner.c

@@ -868,10 +868,11 @@ void runner_do_unskip_mapper(void *map_data, int num_elements,
                              void *extra_data) {
 
   struct engine *e = (struct engine *)extra_data;
-  struct cell *cells = (struct cell *)map_data;
+  struct space *s = e->s;
+  int *local_cells = (int *)map_data;
 
   for (int ind = 0; ind < num_elements; ind++) {
-    struct cell *c = &cells[ind];
+    struct cell *c = &s->cells_top[local_cells[ind]];
     if (c != NULL) runner_do_unskip(c, e);
   }
 }

src/space.c

@@ -389,6 +389,7 @@ void space_regrid(struct space *s, int verbose) {
     /* Free the old cells, if they were allocated. */
     if (s->cells_top != NULL) {
       space_free_cells(s);
+      free(s->local_cells_top);
       free(s->cells_top);
       free(s->multipoles_top);
     }
@@ -420,6 +421,12 @@ void space_regrid(struct space *s, int verbose) {
       bzero(s->multipoles_top, s->nr_cells * sizeof(struct gravity_tensors));
     }
 
+    /* Allocate the indices of local cells */
+    if (posix_memalign((void *)&s->local_cells_top, SWIFT_STRUCT_ALIGNMENT,
+                       s->nr_cells * sizeof(int)) != 0)
+      error("Failed to allocate indices of local top-level cells.");
+    bzero(s->local_cells_top, s->nr_cells * sizeof(int));
+
     /* Set the cells' locks */
     for (int k = 0; k < s->nr_cells; k++) {
       if (lock_init(&s->cells_top[k].lock) != 0)
@@ -496,6 +503,16 @@ void space_regrid(struct space *s, int verbose) {
     }
 #endif /* WITH_MPI */
 
+    /* Let's rebuild the list of local top-level cells */
+    s->nr_local_cells = 0;
+    for (int i = 0; i < s->nr_cells; ++i)
+      if (s->cells_top[i].nodeID == s->nodeID) {
+        s->local_cells_top[s->nr_local_cells] = i;
+        ++s->nr_local_cells;
+      }
+    if (verbose)
+      message("Have %d local cells (total=%d)", s->nr_local_cells, s->nr_cells);
+
     // message( "rebuilding upper-level cells took %.3f %s." ,
     // clocks_from_ticks(double)(getticks() - tic), clocks_getunit());
 
@@ -2652,7 +2669,7 @@ void space_init_sparts(struct space *s) {
  * parts with a cutoff below half the cell width are then split
  * recursively.
  */
-void space_init(struct space *s, const struct swift_params *params,
+void space_init(struct space *s, int nodeID, const struct swift_params *params,
                 double dim[3], struct part *parts, struct gpart *gparts,
                 struct spart *sparts, size_t Npart, size_t Ngpart,
                 size_t Nspart, int periodic, int replicate, int gravity,
@@ -2662,6 +2679,7 @@ void space_init(struct space *s, const struct swift_params *params,
   bzero(s, sizeof(struct space));
 
   /* Store everything in the space. */
+  s->nodeID = nodeID;
   s->dim[0] = dim[0];
   s->dim[1] = dim[1];
   s->dim[2] = dim[2];
@@ -3039,6 +3057,7 @@ void space_clean(struct space *s) {
   for (int i = 0; i < s->nr_cells; ++i) cell_clean(&s->cells_top[i]);
   free(s->cells_top);
   free(s->multipoles_top);
+  free(s->local_cells_top);
   free(s->parts);
   free(s->xparts);
   free(s->gparts);

src/space.h

@@ -74,6 +74,9 @@ struct space {
   /*! Extra space information needed for some hydro schemes. */
   struct hydro_space hs;
 
+  /*! The MPI rank of this space */
+  int nodeID;
+
   /*! Are we doing gravity? */
   int gravity;
 
@@ -101,6 +104,9 @@ struct space {
   /*! Total number of cells (top- and sub-) */
   int tot_cells;
 
+  /*! Number of *local* top-level cells */
+  int nr_local_cells;
+
   /*! The (level 0) cells themselves. */
   struct cell *cells_top;
 
@@ -113,6 +119,9 @@ struct space {
   /*! Buffer of unused multipoles for the sub-cells. */
   struct gravity_tensors *multipoles_sub;
 
+  /*! The indices of the *local* top-level cells */
+  int *local_cells_top;
+
   /*! The total number of parts in the space. */
   size_t nr_parts, size_parts;
 
@@ -173,7 +182,7 @@ void space_sparts_sort(struct space *s, int *ind, size_t N, int min, int max,
 void space_getcells(struct space *s, int nr_cells, struct cell **cells);
 int space_getsid(struct space *s, struct cell **ci, struct cell **cj,
                  double *shift);
-void space_init(struct space *s, const struct swift_params *params,
+void space_init(struct space *s, int nodeID, const struct swift_params *params,
                 double dim[3], struct part *parts, struct gpart *gparts,
                 struct spart *sparts, size_t Npart, size_t Ngpart,
                 size_t Nspart, int periodic, int replicate, int gravity,