allocate a separate sort array for each dimensions.

src/cell.c

@@ -1279,7 +1279,11 @@ void cell_check_multipole(struct cell *c, void *data) {
  */
 void cell_clean(struct cell *c) {
 
-  free(c->sort);
+  for (int i = 0; i < 13; i++)
+    if (c->sort[i] != NULL) {
+      free(c->sort[i]);
+      c->sort[i] = NULL;
+      }
 
   /* Recurse */
   for (int k = 0; k < 8; k++)

src/cell.h

@@ -123,7 +123,7 @@ struct cell {
   struct spart *sparts;
 
   /*! Pointer for the sorted indices. */
-  struct entry *sort;
+  struct entry *sort[13];
 
   /*! Pointers to the next level of cells. */
   struct cell *progeny[8];

src/runner.c

@@ -361,13 +361,14 @@ void runner_do_sort(struct runner *r, struct cell *c, int flags, int cleanup,
      were re-set. */
   if (c->sorted == 0) c->ti_sort = r->e->ti_current;
 
-  /* start by allocating the entry arrays. */
-  if (c->sort == NULL) {
-    if ((c->sort = (struct entry *)malloc(sizeof(struct entry) * (count + 1) *
-                                          13)) == NULL)
-      error("Failed to allocate sort memory.");
+  /* start by allocating the entry arrays in the requested dimensions. */
+  for (int j = 0; j < 13; j++) {
+    if ((flags & (1 << j)) && c->sort[j] == NULL) {
+      if ((c->sort[j] = (struct entry *)malloc(sizeof(struct entry) *
+                                               (count + 1))) == NULL)
+        error("Failed to allocate sort memory.");
+    }
   }
-  struct entry *sort = c->sort;
 
   /* Does this cell have any progeny? */
   if (c->split) {
@@ -409,7 +410,7 @@ void runner_do_sort(struct runner *r, struct cell *c, int flags, int cleanup,
       for (int k = 0; k < 8; k++) {
         inds[k] = k;
         if (c->progeny[k] != NULL && c->progeny[k]->count > 0) {
-          fingers[k] = &c->progeny[k]->sort[j * (c->progeny[k]->count + 1)];
+          fingers[k] = c->progeny[k]->sort[j];
           buff[k] = fingers[k]->d;
           off[k] = off[k];
         } else
@@ -426,7 +427,7 @@ void runner_do_sort(struct runner *r, struct cell *c, int flags, int cleanup,
           }
 
       /* For each entry in the new sort list. */
-      finger = &sort[j * (count + 1)];
+      finger = c->sort[j];
       for (int ind = 0; ind < count; ind++) {
 
         /* Copy the minimum into the new sort array. */
@@ -447,8 +448,8 @@ void runner_do_sort(struct runner *r, struct cell *c, int flags, int cleanup,
       } /* Merge. */
 
       /* Add a sentinel. */
-      sort[j * (count + 1) + count].d = FLT_MAX;
-      sort[j * (count + 1) + count].i = 0;
+      c->sort[j][count].d = FLT_MAX;
+      c->sort[j][count].i = 0;
 
       /* Mark as sorted. */
       atomic_or(&c->sorted, 1 << j);
@@ -484,19 +485,19 @@ void runner_do_sort(struct runner *r, struct cell *c, int flags, int cleanup,
       const double px[3] = {parts[k].x[0], parts[k].x[1], parts[k].x[2]};
       for (int j = 0; j < 13; j++)
         if (flags & (1 << j)) {
-          sort[j * (count + 1) + k].i = k;
-          sort[j * (count + 1) + k].d = px[0] * runner_shift[j][0] +
-                                        px[1] * runner_shift[j][1] +
-                                        px[2] * runner_shift[j][2];
+          c->sort[j][k].i = k;
+          c->sort[j][k].d = px[0] * runner_shift[j][0] +
+                            px[1] * runner_shift[j][1] +
+                            px[2] * runner_shift[j][2];
         }
     }
 
     /* Add the sentinel and sort. */
     for (int j = 0; j < 13; j++)
       if (flags & (1 << j)) {
-        sort[j * (count + 1) + count].d = FLT_MAX;
-        sort[j * (count + 1) + count].i = 0;
-        runner_do_sort_ascending(&sort[j * (count + 1)], count);
+        c->sort[j][count].d = FLT_MAX;
+        c->sort[j][count].i = 0;
+        runner_do_sort_ascending(c->sort[j], count);
         atomic_or(&c->sorted, 1 << j);
       }
   }
@@ -505,7 +506,7 @@ void runner_do_sort(struct runner *r, struct cell *c, int flags, int cleanup,
   /* Verify the sorting. */
   for (int j = 0; j < 13; j++) {
     if (!(flags & (1 << j))) continue;
-    finger = &sort[j * (count + 1)];
+    finger = c->sort[j];
     for (int k = 1; k < count; k++) {
       if (finger[k].d < finger[k - 1].d)
         error("Sorting failed, ascending array.");

src/runner_doiact.h

@@ -640,7 +640,7 @@ void DOPAIR_SUBSET(struct runner *r, struct cell *restrict ci,
     error("Interacting unsorted cells.");
 
   /* Pick-out the sorted lists. */
-  const struct entry *restrict sort_j = &cj->sort[sid * (cj->count + 1)];
+  const struct entry *restrict sort_j = cj->sort[sid];
   const float dxj = cj->dx_max_sort;
 
   /* Parts are on the left? */
@@ -906,8 +906,8 @@ void DOPAIR1(struct runner *r, struct cell *ci, struct cell *cj, const int sid,
   for (int k = 0; k < 3; k++) rshift += shift[k] * runner_shift[sid][k];
 
   /* Pick-out the sorted lists. */
-  const struct entry *restrict sort_i = &ci->sort[sid * (ci->count + 1)];
-  const struct entry *restrict sort_j = &cj->sort[sid * (cj->count + 1)];
+  const struct entry *restrict sort_i = ci->sort[sid];
+  const struct entry *restrict sort_j = cj->sort[sid];
 
 #ifdef SWIFT_DEBUG_CHECKS
   /* Check that the dx_max_sort values in the cell are indeed an upper
@@ -1205,8 +1205,8 @@ void DOPAIR2(struct runner *r, struct cell *ci, struct cell *cj) {
   for (int k = 0; k < 3; k++) rshift += shift[k] * runner_shift[sid][k];
 
   /* Pick-out the sorted lists. */
-  struct entry *restrict sort_i = &ci->sort[sid * (ci->count + 1)];
-  struct entry *restrict sort_j = &cj->sort[sid * (cj->count + 1)];
+  struct entry *restrict sort_i = ci->sort[sid];
+  struct entry *restrict sort_j = cj->sort[sid];
 
 #ifdef SWIFT_DEBUG_CHECKS
   /* Check that the dx_max_sort values in the cell are indeed an upper

src/space.c

@@ -233,10 +233,11 @@ void space_rebuild_recycle_mapper(void *map_data, int num_elements,
     c->xparts = NULL;
     c->gparts = NULL;
     c->sparts = NULL;
-    if (c->sort != NULL) {
-      free(c->sort);
-      c->sort = NULL;
-    }
+    for (int i = 0; i < 13; i++)
+      if (c->sort[i] != NULL) {
+        free(c->sort[i]);
+        c->sort[i] = NULL;
+      }
 #if WITH_MPI
     c->recv_xv = NULL;
     c->recv_rho = NULL;
@@ -1804,10 +1805,11 @@ void space_gparts_sort_mapper(void *map_data, int num_elements,
  */
 void space_map_clearsort(struct cell *c, void *data) {
 
-  if (c->sort != NULL) {
-    free(c->sort);
-    c->sort = NULL;
-  }
+  for (int i = 0; i < 13; i++)
+    if (c->sort[i] != NULL) {
+      free(c->sort[i]);
+      c->sort[i] = NULL;
+    }
 }
 
 /**
@@ -2305,7 +2307,8 @@ void space_recycle(struct space *s, struct cell *c) {
     error("Failed to destroy spinlocks.");
 
   /* Clear this cell's sort arrays. */
-  if (c->sort != NULL) free(c->sort);
+  for (int i = 0; i < 13; i++)
+    if (c->sort[i] != NULL) free(c->sort[i]);
 
   /* Lock the space. */
   lock_lock(&s->lock);
@@ -2357,7 +2360,8 @@ void space_recycle_list(struct space *s, struct cell *cell_list_begin,
       error("Failed to destroy spinlocks.");
 
     /* Clear this cell's sort arrays. */
-    if (c->sort != NULL) free(c->sort);
+    for (int i = 0; i < 13; i++)
+      if (c->sort[i] != NULL) free(c->sort[i]);
 
     /* Count this cell. */
     count += 1;