Reinstate the debugging checks

examples/EAGLE_low_z/EAGLE_6/eagle_6.yml

@@ -53,12 +53,14 @@ Statistics:
 
 # Parameters for the self-gravity scheme
 Gravity:
-  eta:                    0.025    # Constant dimensionless multiplier for time integration.
-  MAC:                    adaptive
-  epsilon_fmm:            0.01
-  theta_cr:               0.7      # Opening angle (Multipole acceptance criterion)
-  mesh_side_length:       16
-  dithering:              0
+  eta:                      0.025    # Constant dimensionless multiplier for time integration.
+  MAC:                      adaptive
+  epsilon_fmm:              0.001
+  theta_cr:                 0.7      # Opening angle (Multipole acceptance criterion)
+  use_tree_below_softening: 0
+  allow_truncation_in_MAC:  1
+  mesh_side_length:         16
+  dithering:                0
   comoving_DM_softening:         0.0026994 # Comoving DM softening length (in internal units).
   max_physical_DM_softening:     0.0007    # Max physical DM softening length (in internal units).
   comoving_baryon_softening:     0.0026994 # Comoving DM softening length (in internal units).

src/multipole_accept.h

@@ -59,8 +59,8 @@ __attribute__((nonnull, pure)) INLINE static int gravity_M2L_accept(
   const float rho_B = use_rebuild_sizes ? B->r_max_rebuild : B->r_max;
 
 #ifdef SWIFT_DEBUG_CHECKS
-  // if (rho_A == 0.) error("Size of multipole A is 0!");
-  // if (rho_B == 0.) error("Size of multipole B is 0!");
+  if (rho_A == 0.) error("Size of multipole A is 0!");
+  if (rho_B == 0.) error("Size of multipole B is 0!");
 #endif
 
   /* Compute the error estimator (without the 1/M_B term that cancels out) */

src/runner_doiact_grav.c

@@ -538,7 +538,11 @@ static INLINE void runner_dopair_grav_pm_full(
 
     const float r2 = dx * dx + dy * dy + dz * dz;
 
-    // MATTHIEU: Do we want a check that the particle can indeed use M2P?
+#ifdef SWIFT_DEBUG_CHECKS
+    if (!gravity_M2P_accept(e->gravity_properties, &gparts_i[pid],
+                            cj->grav.multipole, r2 * 1.01, periodic))
+      error("use_mpole[i] set when M2P accept fails");
+#endif
 
     /* Interact! */
     float f_x, f_y, f_z, pot_ij;
@@ -670,7 +674,11 @@ static INLINE void runner_dopair_grav_pm_truncated(
 
     const float r2 = dx * dx + dy * dy + dz * dz;
 
-    // MATTHIEU: Do we want a check that the particle can indeed use M2P?
+#ifdef SWIFT_DEBUG_CHECKS
+    if (!gravity_M2P_accept(e->gravity_properties, &gparts_i[pid],
+                            cj->grav.multipole, r2 * 1.01, /*periodic=*/1))
+      error("use_mpole[i] set when M2P accept fails");
+#endif
 
     /* Interact! */
     float f_x, f_y, f_z, pot_ij;