From 200c35cc916310307aedae2a8c9f1cf974b43bf7 Mon Sep 17 00:00:00 2001
From: Josh Borrow <joshua.borrow@durham.ac.uk>
Date: Mon, 21 May 2018 07:38:07 +0100
Subject: [PATCH] Added debugging checks for negative values

---
 src/hydro/PressureEnergy/hydro_iact.h | 41 +++++++++++++++++++++++++++
 1 file changed, 41 insertions(+)

diff --git a/src/hydro/PressureEnergy/hydro_iact.h b/src/hydro/PressureEnergy/hydro_iact.h
index 28ed5c32fc..fa973ad9cd 100644
--- a/src/hydro/PressureEnergy/hydro_iact.h
+++ b/src/hydro/PressureEnergy/hydro_iact.h
@@ -63,10 +63,12 @@ __attribute__((always_inline)) INLINE static void runner_iact_density(
   /* Compute density of pi. */
   const float hi_inv = 1.f / hi;
   const float ui = r * hi_inv;
+
   kernel_deval(ui, &wi, &wi_dx);
 
   pi->rho += mj * wi;
   pi->density.rho_dh -= mj * (hydro_dimension * wi + ui * wi_dx);
+
   pi->pressure_bar += mj * wi * pj->u;
   pi->density.pressure_bar_dh -=
       mj * pj->u * (hydro_dimension * wi + ui * wi_dx);
@@ -113,6 +115,26 @@ __attribute__((always_inline)) INLINE static void runner_iact_density(
   pj->density.rot_v[0] += facj * curlvr[0];
   pj->density.rot_v[1] += facj * curlvr[1];
   pj->density.rot_v[2] += facj * curlvr[2];
+
+#ifdef SWIFT_DEBUG_CHECKS
+  if (pi->u < -FLT_MIN)
+    error("Particle %lld has negative u=%.3e\n", pi->id, pi->u);
+
+  if (pj->u < -FLT_MIN)
+    error("Particle %lld has negative u=%.3e\n", pj->id, pj->u);
+
+  if (pi->pressure_bar < -FLT_MIN)
+    error("Particle %lld has negative P_bar=%.3e\n", pi->id, pi->pressure_bar);
+
+  if (pi->rho < -FLT_MIN)
+    error("Particle %lld has negative rho=%.3e\n", pi->id, pi->rho);
+
+  if (mj < -FLT_MIN)
+    error("Particle %lld has negative m=%.3e\n", pj->id, mj);
+
+  if (wi < -FLT_MIN)
+    error("Particle %lld has negative rho=%.3e\n", pi->id, wi);
+#endif
 }
 
 /**
@@ -146,7 +168,9 @@ __attribute__((always_inline)) INLINE static void runner_iact_nonsym_density(
 
   pi->rho += mj * wi;
   pi->density.rho_dh -= mj * (hydro_dimension * wi + ui * wi_dx);
+
   pi->pressure_bar += mj * wi * pj->u;
+
   pi->density.pressure_bar_dh -=
       mj * pj->u * (hydro_dimension * wi + ui * wi_dx);
   pi->density.wcount += wi;
@@ -171,6 +195,23 @@ __attribute__((always_inline)) INLINE static void runner_iact_nonsym_density(
   pi->density.rot_v[0] += faci * curlvr[0];
   pi->density.rot_v[1] += faci * curlvr[1];
   pi->density.rot_v[2] += faci * curlvr[2];
+
+#ifdef SWIFT_DEBUG_CHECKS
+  if (pi->u < -FLT_MIN)
+    error("Particle %lld has negative u=%.3e\n", pi->id, pi->u);
+
+  if (pi->pressure_bar < -FLT_MIN)
+    error("Particle %lld has negative P_bar=%.3e\n", pi->id, pi->pressure_bar);
+
+  if (pi->rho < -FLT_MIN)
+    error("Particle %lld has negative rho=%.3e\n", pi->id, pi->rho);
+
+  if (mj < -FLT_MIN)
+    error("Particle %lld has negative m=%.3e\n", pj->id, mj);
+
+  if (wi < -FLT_MIN)
+    error("Particle %lld has negative rho=%.3e\n", pi->id, wi);
+#endif
 }
 
 /**
-- 
GitLab