Upgraded neighbour finding matching the COLIBRE fork

Import the hydro and stars neighbour finding implemented in the colibre fast branch minus the task multiplexing capability. A very similar version is currently running on 20k cores of a galaxy formation test so that should be doing fine.

This introduces significantly faster neighbour finding in particle distributions with strong density gradients. For instance planetary applications or problematic feedback-disturbed galaxies.

Changes:

• Rewrite the recursion logic in the hydro and stars sub-task:
  • The interaction functions have extra parameters to optionally consider particles only between 0.5 * width and width
  • The subtask recursion now continues to lower levels if we reach a level where h is too large. From that level on, we will use the newly introduced feature of only considering particles in the appropriate range of h.
  • When recursing, only the h_max of active particles is considered, not all particles.

Note that this makes use of the functionality whereby every particle carries its own depth in the tree and the cells have information about the min/max allowed h in them. These two things have been merged in the main code already (but do nothing there besides correctly tracking the values).

Supersedes !1353 (closed).

src/cell.h

@@ -809,6 +809,36 @@ cell_can_recurse_in_pair_hydro_task(const struct cell *c) {
                        c->hydro.dx_max_part_old) < 0.5f * c->dmin);
 }
 
+/**
+ * @brief Can a sub-pair hydro task recurse to a lower level based
+ * on the status of the particles in the cell.
+ *
+ * @param c The #cell.
+ */
+__attribute__((always_inline)) INLINE static int
+cell_can_recurse_in_subpair_hydro_task(const struct cell *c) {
+
+  /* If so, is the cut-off radius plus the max distance the parts have moved */
+  /* smaller than the sub-cell sizes ? */
+  return ((kernel_gamma * c->hydro.h_max_active + c->hydro.dx_max_part_old) <
+          0.5f * c->dmin);
+}
+
+/**
+ * @brief Can a sub-pair hydro task recurse to a lower level based
+ * on the status of the particles in the cell.
+ *
+ * @param c The #cell.
+ */
+__attribute__((always_inline)) INLINE static int
+cell_can_recurse_in_subpair2_hydro_task(const struct cell *c) {
+
+  /* If so, is the cut-off radius plus the max distance the parts have moved */
+  /* smaller than the sub-cell sizes ? */
+  return ((kernel_gamma * c->hydro.h_max + c->hydro.dx_max_part) <
+          0.5f * c->dmin);
+}
+
 /**
  * @brief Can a sub-self hydro task recurse to a lower level based
  * on the status of the particles in the cell.
@@ -822,6 +852,32 @@ cell_can_recurse_in_self_hydro_task(const struct cell *c) {
   return c->split && (kernel_gamma * c->hydro.h_max_old < 0.5f * c->dmin);
 }
 
+/**
+ * @brief Can a sub-self hydro task recurse to a lower level based
+ * on the status of the particles in the cell.
+ *
+ * @param c The #cell.
+ */
+__attribute__((always_inline)) INLINE static int
+cell_can_recurse_in_subself_hydro_task(const struct cell *c) {
+
+  /* Is the cell not smaller than the smoothing length? */
+  return (kernel_gamma * c->hydro.h_max_active < 0.5f * c->dmin);
+}
+
+/**
+ * @brief Can a hydro task recurse to a lower level based
+ * on the status of the particles in the cell.
+ *
+ * @param c The #cell.
+ */
+__attribute__((always_inline)) INLINE static int
+cell_can_recurse_in_subself2_hydro_task(const struct cell *c) {
+
+  /* Is the cell split and not smaller than the smoothing length? */
+  return c->split && (kernel_gamma * c->hydro.h_max < 0.5f * c->dmin);
+}
+
 /**
  * @brief Can a sub-pair star task recurse to a lower level based
  * on the status of the particles in the cell.
@@ -830,18 +886,29 @@ cell_can_recurse_in_self_hydro_task(const struct cell *c) {
  * @param cj The #cell with hydro parts.
  */
 __attribute__((always_inline)) INLINE static int
-cell_can_recurse_in_pair_stars_task(const struct cell *ci,
-                                    const struct cell *cj) {
+cell_can_recurse_in_pair_stars_task(const struct cell *c) {
 
   /* Is the cell split ? */
   /* If so, is the cut-off radius plus the max distance the parts have moved */
   /* smaller than the sub-cell sizes ? */
   /* Note: We use the _old values as these might have been updated by a drift */
-  return ci->split && cj->split &&
-         ((kernel_gamma * ci->stars.h_max_old + ci->stars.dx_max_part_old) <
-          0.5f * ci->dmin) &&
-         ((kernel_gamma * cj->hydro.h_max_old + cj->hydro.dx_max_part_old) <
-          0.5f * cj->dmin);
+  return c->split && ((kernel_gamma * c->stars.h_max_old +
+                       c->stars.dx_max_part_old) < 0.5f * c->dmin);
+}
+
+/**
+ * @brief Can a sub-pair stars task recurse to a lower level based
+ * on the status of the particles in the cell.
+ *
+ * @param c The #cell.
+ */
+__attribute__((always_inline)) INLINE static int
+cell_can_recurse_in_subpair_stars_task(const struct cell *c) {
+
+  /* If so, is the cut-off radius plus the max distance the parts have moved */
+  /* smaller than the sub-cell sizes ? */
+  return ((kernel_gamma * c->stars.h_max_active + c->stars.dx_max_part_old) <
+          0.5f * c->dmin);
 }
 
 /**
@@ -858,6 +925,19 @@ cell_can_recurse_in_self_stars_task(const struct cell *c) {
          (kernel_gamma * c->hydro.h_max_old < 0.5f * c->dmin);
 }
 
+/**
+ * @brief Can a sub-self stars task recurse to a lower level based
+ * on the status of the particles in the cell.
+ *
+ * @param c The #cell.
+ */
+__attribute__((always_inline)) INLINE static int
+cell_can_recurse_in_subself_stars_task(const struct cell *c) {
+
+  /* Is the cell not smaller than the smoothing length? */
+  return (kernel_gamma * c->stars.h_max_active < 0.5f * c->dmin);
+}
+
 /**
  * @brief Can a sub-pair sink task recurse to a lower level based
  * on the status of the particles in the cell.