diff --git a/src/hydro/PressureEnergy/hydro_iact.h b/src/hydro/PressureEnergy/hydro_iact.h index 28ed5c32fc22a12f2c6d732191a8618de83ceeb0..fa973ad9cd91a8c080c9428e7a729fbe1e58cfe6 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 } /**