From 3d12f3828823e511ce41157053bd031ef77c4700 Mon Sep 17 00:00:00 2001
From: "Peter W. Draper" <p.w.draper@durham.ac.uk>
Date: Mon, 17 Jul 2017 14:01:19 +0100
Subject: [PATCH] Add checks for particles that are at same position, this
 isn't allowed, and update the cell h_max to that of the particles after the
 initial ghost tasks have ran. This is needed before regridding if the h
 values of the particles in the IC are set to be smaller than reality (ends up
 looking like rebuilding failed as the particles are too far outside their
 cells).

---
 src/engine.c | 32 ++++++++++++++++++++++++++++++++
 1 file changed, 32 insertions(+)

diff --git a/src/engine.c b/src/engine.c
index 2680080ee6..2b6f0f2f42 100644
--- a/src/engine.c
+++ b/src/engine.c
@@ -3359,6 +3359,38 @@ void engine_init_particles(struct engine *e, int flag_entropy_ICs,
   /* Recover the (integer) end of the next time-step */
   engine_collect_timestep_and_rebuild(e, 1);
 
+  /* Check if any particles have the same position and see if we need to
+   * correct the cell h_max to match possible particle updates in the ghost
+   * tasks. */
+  if (s->cells_top != NULL) {
+
+    /* Sorting should put the same positions next to each other... */
+    int failed = 0;
+    double *prev_x = s->parts[0].x;
+    for (size_t k = 1; k < s->nr_parts; k++) {
+      if (prev_x[0] == s->parts[k].x[0] &&
+          prev_x[1] == s->parts[k].x[1] &&
+          prev_x[2] == s->parts[k].x[2]) {
+        message("Two particles occupy location: %f %f %f",
+              prev_x[0], prev_x[1], prev_x[2]);
+        failed = 1;
+      }
+      prev_x = s->parts[k].x;
+    }
+    if (failed)
+      error("Cannot have particles with the same positions");
+
+    for (int i = 0; i < s->nr_cells; i++) {
+      struct cell *c = &s->cells_top[i];
+      if (c->nodeID == engine_rank) {
+        float part_h_max = c->parts[0].h;
+        for (int k = 1; k < c->count; k++) {
+          if (c->parts[k].h > part_h_max) part_h_max = c->parts[k].h;
+        }
+        c->h_max = max(part_h_max, c->h_max);
+      }
+    }
+  }
   clocks_gettime(&time2);
 
 #ifdef SWIFT_DEBUG_CHECKS
-- 
GitLab