Stars: optimized dopair1 works

src/runner_doiact_stars.h

@@ -259,7 +259,6 @@ void DO_SYM_PAIR1_STARS(struct runner *r, struct cell *ci, struct cell *cj, cons
 
   if (do_ci_stars) {
     /* Pick-out the sorted lists. */
-    const struct entry *restrict sort_i = ci->stars.sort[sid];
     const struct entry *restrict sort_j = cj->hydro.sort[sid];
 
 #ifdef SWIFT_DEBUG_CHECKS
@@ -285,21 +284,23 @@ void DO_SYM_PAIR1_STARS(struct runner *r, struct cell *ci, struct cell *cj, cons
       (ci->stars.dx_max_sort + cj->hydro.dx_max_sort) - rshift;
 
     /* Loop over the sparts in ci. */
-    for (int pid = count_i - 1; pid >= 0; pid--) {
+    for (int pid = 0; pid < count_i; pid++) {
 
       /* Get a hold of the ith part in ci. */
-      struct spart *restrict spi = &sparts_i[sort_i[pid].i];
-
-      if (sort_i[pid].d + spi->h * kernel_gamma + dx_max_rshift < dj_min)
-	continue;
+      struct spart *restrict spi = &sparts_i[pid];
+      const float hi = spi->h;
 
       /* Skip inactive particles */
       if (!spart_is_active(spi, e) || spart_is_inhibited(spi, e)) continue;
 
-      const float hi = spi->h;
-
+      /* Compute distance from the other cell. */
+      const double px[3] = {spi->x[0], spi->x[1], spi->x[2]};
+      float dist = px[0] * runner_shift[sid][0] +
+	           px[1] * runner_shift[sid][1] +
+                   px[2] * runner_shift[sid][2];
+      
       /* Is there anything we need to interact with ? */
-      const double di = sort_i[pid].d + hi * kernel_gamma + dx_max_rshift;
+      const double di = dist + hi * kernel_gamma + dx_max_rshift;
       if (di < dj_min) continue;
 
       /* Get some additional information about pi */
@@ -374,7 +375,6 @@ void DO_SYM_PAIR1_STARS(struct runner *r, struct cell *ci, struct cell *cj, cons
   if (do_cj_stars) {
     /* Pick-out the sorted lists. */
     const struct entry *restrict sort_i = ci->hydro.sort[sid];
-    const struct entry *restrict sort_j = cj->stars.sort[sid];
 
 #ifdef SWIFT_DEBUG_CHECKS
     /* Some constants used to checks that the parts are in the right frame */
@@ -395,25 +395,28 @@ void DO_SYM_PAIR1_STARS(struct runner *r, struct cell *ci, struct cell *cj, cons
     struct part *restrict parts_i = ci->hydro.parts;
     struct spart *restrict sparts_j = cj->stars.parts;
     const double di_max = sort_i[count_i - 1].d - rshift;
-    const float dx_max = (ci->hydro.dx_max_sort + cj->stars.dx_max_sort);
+    const float dx_max_rshift =
+      (ci->hydro.dx_max_sort + cj->stars.dx_max_sort) + rshift;
 
     /* Loop over the parts in cj. */
     for (int pjd = 0; pjd < count_j; pjd++) {
       
       /* Get a hold of the jth part in cj. */
-      struct spart *spj = &sparts_j[sort_j[pjd].i];
+      struct spart *spj = &sparts_j[pjd];
+      const float hj = spj->h;
 
-      /* Is it possible to interact with a particle? */
-      if (sort_j[pjd].d - spj->h * kernel_gamma - dx_max > di_max)
-	continue;
       /* Skip inactive particles */
       if (!spart_is_active(spj, e) || spart_is_inhibited(spj, e)) continue;
 
-      const float hj = spj->h;
 
+      /* Compute distance from the other cell. */
+      const double px[3] = {spj->x[0], spj->x[1], spj->x[2]};
+      float dist = px[0] * runner_shift[sid][0] +
+	           px[1] * runner_shift[sid][1] +
+                   px[2] * runner_shift[sid][2];
 
       /* Is there anything we need to interact with ? */
-      const double dj = sort_j[pjd].d - hj * kernel_gamma - dx_max + rshift;
+      const double dj = dist - hj * kernel_gamma - dx_max_rshift;
       if (dj - rshift > di_max) continue;
 
       /* Get some additional information about pj */