From a0931f6111c6df55db435c6e45348e75b66f0f5e Mon Sep 17 00:00:00 2001
From: Matthieu Schaller <matthieu.schaller@durham.ac.uk>
Date: Wed, 12 Apr 2017 09:46:48 +0100
Subject: [PATCH] Zero the multipoles when rebuilding if there are no gparts.

---
 src/multipole.h |  5 +++++
 src/space.c     | 11 ++++++++++-
 2 files changed, 15 insertions(+), 1 deletion(-)

diff --git a/src/multipole.h b/src/multipole.h
index 0daef6de24..47f4764daf 100644
--- a/src/multipole.h
+++ b/src/multipole.h
@@ -354,6 +354,11 @@ INLINE static void gravity_field_tensors_print(const struct grav_tensor *l) {
 #endif
 }
 
+INLINE static void gravity_multipole_init(struct multipole *m) {
+
+  bzero(m, sizeof(struct multipole));
+}
+
 /**
  * @brief Prints the content of a #multipole to stdout.
  *
diff --git a/src/space.c b/src/space.c
index 72ab8baa44..04bb3a513f 100644
--- a/src/space.c
+++ b/src/space.c
@@ -2150,7 +2150,16 @@ void space_split_recursive(struct space *s, struct cell *c,
     ti_beg_max = get_integer_time_begin(e->ti_current + 1, time_bin_max);
 
     /* Construct the multipole and the centre of mass*/
-    if (s->gravity && gcount > 0) gravity_P2M(c->multipole, c->gparts, c->gcount);
+    if (s->gravity) {
+      if(gcount > 0) 
+	gravity_P2M(c->multipole, c->gparts, c->gcount);
+      else {
+	gravity_multipole_init(&c->multipole->m_pole);
+	c->multipole->CoM[0] = c->loc[0] + c->width[0] / 2.;
+	c->multipole->CoM[1] = c->loc[1] + c->width[1] / 2.;
+	c->multipole->CoM[2] = c->loc[2] + c->width[2] / 2.;
+      }
+    }
   }
 
   /* Set the values for this cell. */
-- 
GitLab