From 2f408b256ce8e0a59fe71ff13018bf155f050374 Mon Sep 17 00:00:00 2001
From: Bert Vandenbroucke <bert.vandenbroucke@ugent.be>
Date: Wed, 10 Aug 2016 16:53:54 +0100
Subject: [PATCH] Tried to fix test125cells. Found a few bugs in the process.
---
src/hydro/Gizmo/hydro.h | 5 +++--
src/hydro/Gizmo/hydro_gradients_sph.h | 6 ++++++
src/hydro/Gizmo/hydro_iact.h | 8 ++++++--
tests/test125cells.c | 19 +++++++++++++++++++
4 files changed, 34 insertions(+), 4 deletions(-)
diff --git a/src/hydro/Gizmo/hydro.h b/src/hydro/Gizmo/hydro.h
index a2fd75e9c8..7aa50d96b6 100644
--- a/src/hydro/Gizmo/hydro.h
+++ b/src/hydro/Gizmo/hydro.h
@@ -62,6 +62,9 @@ __attribute__((always_inline)) INLINE static void hydro_first_init_part(
__attribute__((always_inline)) INLINE static void hydro_init_part(
struct part* p) {
+ /* We need to do this before we reset the volume */
+ hydro_gradients_init_density_loop(p);
+
p->density.wcount = 0.0f;
p->density.wcount_dh = 0.0f;
p->geometry.volume = 0.0f;
@@ -74,8 +77,6 @@ __attribute__((always_inline)) INLINE static void hydro_init_part(
p->geometry.matrix_E[2][0] = 0.0f;
p->geometry.matrix_E[2][1] = 0.0f;
p->geometry.matrix_E[2][2] = 0.0f;
-
- hydro_gradients_init_density_loop(p);
}
/**
diff --git a/src/hydro/Gizmo/hydro_gradients_sph.h b/src/hydro/Gizmo/hydro_gradients_sph.h
index 8c9d7be052..29e36ec190 100644
--- a/src/hydro/Gizmo/hydro_gradients_sph.h
+++ b/src/hydro/Gizmo/hydro_gradients_sph.h
@@ -38,6 +38,12 @@ hydro_gradients_init_density_loop(struct part *p) {
p->conserved.momentum[2] * p->conserved.momentum[2]) /
p->conserved.mass) /
p->geometry.volume;
+ } else {
+ p->primitives.rho = 0.0f;
+ p->primitives.v[0] = 0.0f;
+ p->primitives.v[1] = 0.0f;
+ p->primitives.v[2] = 0.0f;
+ p->primitives.P = 0.0f;
}
p->primitives.gradients.rho[0] = 0.0f;
diff --git a/src/hydro/Gizmo/hydro_iact.h b/src/hydro/Gizmo/hydro_iact.h
index 0f847e1b85..82b26f3c78 100644
--- a/src/hydro/Gizmo/hydro_iact.h
+++ b/src/hydro/Gizmo/hydro_iact.h
@@ -155,8 +155,12 @@ __attribute__((always_inline)) INLINE static void runner_iact_fluxes_common(
dtj = pj->force.dt;
/* calculate the maximal signal velocity */
- vmax =
- sqrtf(hydro_gamma * Wi[4] / Wi[0]) + sqrtf(hydro_gamma * Wj[4] / Wj[0]);
+ if (Wi[0] && Wj[0]) {
+ vmax =
+ sqrtf(hydro_gamma * Wi[4] / Wi[0]) + sqrtf(hydro_gamma * Wj[4] / Wj[0]);
+ } else {
+ vmax = 0.0f;
+ }
dvdotdx = (Wi[1] - Wj[1]) * dx[0] + (Wi[2] - Wj[2]) * dx[1] +
(Wi[3] - Wj[3]) * dx[2];
if (dvdotdx > 0.) {
diff --git a/tests/test125cells.c b/tests/test125cells.c
index f2cbe6606f..a839c84f73 100644
--- a/tests/test125cells.c
+++ b/tests/test125cells.c
@@ -253,6 +253,25 @@ struct cell *make_cell(size_t n, const double offset[3], double size, double h,
xpart->v_full[1] = part->v[1];
xpart->v_full[2] = part->v[2];
hydro_first_init_part(part, xpart);
+
+#if defined(GIZMO_SPH)
+ part->geometry.volume = volume;
+ part->primitives.rho = density;
+ part->primitives.v[0] = part->v[0];
+ part->primitives.v[1] = part->v[1];
+ part->primitives.v[2] = part->v[2];
+ part->conserved.mass = part->mass;
+ part->conserved.momentum[0] = part->conserved.mass * part->v[0];
+ part->conserved.momentum[1] = part->conserved.mass * part->v[1];
+ part->conserved.momentum[2] = part->conserved.mass * part->v[2];
+ part->conserved.energy =
+ part->primitives.P / hydro_gamma_minus_one * volume +
+ 0.5f * (part->conserved.momentum[0] * part->conserved.momentum[0] +
+ part->conserved.momentum[1] * part->conserved.momentum[1] +
+ part->conserved.momentum[2] * part->conserved.momentum[2]) /
+ part->conserved.mass;
+#endif
+
++part;
++xpart;
}
--
GitLab