diff --git a/src/cell.c b/src/cell.c
index cc77f9c94bcf37db3e99ec9ab19a58c938279949..2c4afe81a060103bb3dc239e23233c75c0458ec0 100644
--- a/src/cell.c
+++ b/src/cell.c
@@ -1419,25 +1419,64 @@ void cell_make_multipoles(struct cell *c, integertime_t ti_current) {
   c->ti_old_multipole = ti_current;
 }
 
+/**
+ * @brief Recursively verify that the multipoles are the sum of their progenies.
+ *
+ * This function does not check whether the multipoles match the particle
+ * content as we may not have received the particles.
+ *
+ * @param c The #cell to recursively search and verify.
+ */
+void cell_check_foreign_multipole(const struct cell *c) {
+
+#ifdef SWIFT_DEBUG_CHECKS
+
+  if (c->split) {
+
+    double M_000 = 0.;
+    long long num_gpart = 0;
+
+    for (int k = 0; k < 8; k++) {
+      const struct cell *cp = c->progeny[k];
+
+      if (cp != NULL) {
+
+        /* Check the mass */
+        M_000 += cp->multipole->m_pole.M_000;
+
+        /* Check the number of particles */
+        num_gpart += cp->multipole->m_pole.num_gpart;
+
+        /* Now recurse */
+        cell_check_foreign_multipole(cp);
+      }
+    }
+
+    if (num_gpart != c->multipole->m_pole.num_gpart)
+      error("Sum of particles in progenies does not match");
+  }
+
+#else
+  error("Calling debugging code without debugging flag activated.");
+#endif
+}
+
 /**
  * @brief Computes the multi-pole brutally and compare to the
  * recursively computed one.
  *
  * @param c Cell to act upon
- * @param data Unused parameter
  */
-void cell_check_multipole(struct cell *c, void *data) {
+void cell_check_multipole(struct cell *c) {
 
 #ifdef SWIFT_DEBUG_CHECKS
   struct gravity_tensors ma;
   const double tolerance = 1e-3; /* Relative */
 
-  return;
-
   /* First recurse */
   if (c->split)
     for (int k = 0; k < 8; k++)
-      if (c->progeny[k] != NULL) cell_check_multipole(c->progeny[k], NULL);
+      if (c->progeny[k] != NULL) cell_check_multipole(c->progeny[k]);
 
   if (c->gcount > 0) {
 
@@ -1456,7 +1495,7 @@ void cell_check_multipole(struct cell *c, void *data) {
     }
 
     /* Check that the upper limit of r_max is good enough */
-    if (!(c->multipole->r_max >= ma.r_max)) {
+    if (!(1.1 * c->multipole->r_max >= ma.r_max)) {
       error("Upper-limit r_max=%e too small. Should be >=%e.",
             c->multipole->r_max, ma.r_max);
     } else if (c->multipole->r_max * c->multipole->r_max >
diff --git a/src/cell.h b/src/cell.h
index 747b79f05a799b9401c67c23b72b722d385c5ced..bd450d3b307bdfcab6b16ab7b6cc7f954e9674b4 100644
--- a/src/cell.h
+++ b/src/cell.h
@@ -515,7 +515,8 @@ int cell_link_gparts(struct cell *c, struct gpart *gparts);
 int cell_link_sparts(struct cell *c, struct spart *sparts);
 void cell_clean_links(struct cell *c, void *data);
 void cell_make_multipoles(struct cell *c, integertime_t ti_current);
-void cell_check_multipole(struct cell *c, void *data);
+void cell_check_multipole(struct cell *c);
+void cell_check_foreign_multipole(const struct cell *c);
 void cell_clean(struct cell *c);
 void cell_check_part_drift_point(struct cell *c, void *data);
 void cell_check_gpart_drift_point(struct cell *c, void *data);
diff --git a/src/engine.c b/src/engine.c
index 6c001e8279ec0eefa88f6ff1294a91357e51c87a..c3a024edc1651f9617420af5d2c2b9208f8600fa 100644
--- a/src/engine.c
+++ b/src/engine.c
@@ -3969,6 +3969,18 @@ void engine_rebuild(struct engine *e, int clean_smoothing_length_values) {
   engine_exchange_cells(e);
 #endif
 
+#ifdef SWIFT_DEBUG_CHECKS
+  long long counter = 0;
+
+  /* Let's check that what we received makes sense */
+  for (int i = 0; i < e->s->nr_cells; ++i) {
+    const struct gravity_tensors *m = &e->s->multipoles_top[i];
+    counter += m->m_pole.num_gpart;
+  }
+  if (counter != e->total_nr_gparts)
+    error("Total particles in multipoles inconsistent with engine");
+#endif
+
   /* Re-build the tasks. */
   engine_maketasks(e);
 
@@ -3981,6 +3993,9 @@ void engine_rebuild(struct engine *e, int clean_smoothing_length_values) {
    * previously been active on this rank. */
   space_check_drift_point(e->s, e->ti_current,
                           e->policy & engine_policy_self_gravity);
+
+  for (int k = 0; k < e->s->nr_local_cells; k++)
+    cell_check_foreign_multipole(&e->s->cells_top[e->s->local_cells_top[k]]);
 #endif
 
   /* Run through the tasks and mark as skip or not. */
diff --git a/src/space.c b/src/space.c
index 9a9dd05c237a69d034bf618490c834ade5ee33f7..c3d079bf621dec9e9945f005d1abee8f716ac7bd 100644
--- a/src/space.c
+++ b/src/space.c
@@ -960,7 +960,7 @@ void space_rebuild(struct space *s, int verbose) {
   /* Check that the multipole construction went OK */
   if (s->gravity)
     for (int k = 0; k < s->nr_cells; k++)
-      cell_check_multipole(&s->cells_top[k], NULL);
+      cell_check_multipole(&s->cells_top[k]);
 #endif
 
   /* Clean up any stray sort indices in the cell buffer. */
@@ -2311,6 +2311,7 @@ void space_getcells(struct space *s, int nr_cells, struct cell **cells) {
       if (cells[j]->sort[k] != NULL) free(cells[j]->sort[k]);
     struct gravity_tensors *temp = cells[j]->multipole;
     bzero(cells[j], sizeof(struct cell));
+    bzero(temp, sizeof(struct gravity_tensors));
     cells[j]->multipole = temp;
     cells[j]->nodeID = -1;
     if (lock_init(&cells[j]->lock) != 0 || lock_init(&cells[j]->glock) != 0 ||