From 1cefa9b3bcfc63020d751b55c4769faa255c7034 Mon Sep 17 00:00:00 2001
From: Matthieu Schaller <matthieu.schaller@durham.ac.uk>
Date: Tue, 19 Dec 2017 17:41:59 +0100
Subject: [PATCH] Make the hydro schemes compute the dirfted velocities
themselves when collecting statistics.
---
src/hydro/Default/hydro.h | 17 +++++++++++++++++
src/hydro/Gadget2/hydro.h | 17 +++++++++++++++++
src/hydro/Gizmo/hydro.h | 12 ++++++++++++
src/hydro/Minimal/hydro.h | 17 +++++++++++++++++
src/hydro/PressureEntropy/hydro.h | 17 +++++++++++++++++
src/hydro/Shadowswift/hydro.h | 12 ++++++++++++
src/statistics.c | 21 +++++----------------
7 files changed, 97 insertions(+), 16 deletions(-)
diff --git a/src/hydro/Default/hydro.h b/src/hydro/Default/hydro.h
index 31f0c41720..890380b23a 100644
--- a/src/hydro/Default/hydro.h
+++ b/src/hydro/Default/hydro.h
@@ -97,6 +97,23 @@ __attribute__((always_inline)) INLINE static float hydro_get_mass(
return p->mass;
}
+/**
+ * @brief Returns the velocities drifted to the current time of a particle.
+ *
+ * @param p The particle of interest
+ * @param xp The extended data of the particle.
+ * @param dt The time since the last kick.
+ * @param v (return) The velocities at the current time.
+ */
+__attribute__((always_inline)) INLINE static void hydro_get_drifted_velocities(
+ const struct part *restrict p, const struct xpart *xp, float dt,
+ float v[3]) {
+
+ v[0] = xp->v_full[0] + p->a_hydro[0] * dt;
+ v[1] = xp->v_full[1] + p->a_hydro[1] * dt;
+ v[2] = xp->v_full[2] + p->a_hydro[2] * dt;
+}
+
/**
* @brief Returns the time derivative of internal energy of a particle
*
diff --git a/src/hydro/Gadget2/hydro.h b/src/hydro/Gadget2/hydro.h
index fdd0795cfc..f3c33bc7a0 100644
--- a/src/hydro/Gadget2/hydro.h
+++ b/src/hydro/Gadget2/hydro.h
@@ -106,6 +106,23 @@ __attribute__((always_inline)) INLINE static float hydro_get_mass(
return p->mass;
}
+/**
+ * @brief Returns the velocities drifted to the current time of a particle.
+ *
+ * @param p The particle of interest
+ * @param xp The extended data of the particle.
+ * @param dt The time since the last kick.
+ * @param v (return) The velocities at the current time.
+ */
+__attribute__((always_inline)) INLINE static void hydro_get_drifted_velocities(
+ const struct part *restrict p, const struct xpart *xp, float dt,
+ float v[3]) {
+
+ v[0] = xp->v_full[0] + p->a_hydro[0] * dt;
+ v[1] = xp->v_full[1] + p->a_hydro[1] * dt;
+ v[2] = xp->v_full[2] + p->a_hydro[2] * dt;
+}
+
/**
* @brief Returns the time derivative of internal energy of a particle
*
diff --git a/src/hydro/Gizmo/hydro.h b/src/hydro/Gizmo/hydro.h
index 3b22d46da9..920556114a 100644
--- a/src/hydro/Gizmo/hydro.h
+++ b/src/hydro/Gizmo/hydro.h
@@ -754,6 +754,18 @@ __attribute__((always_inline)) INLINE static float hydro_get_mass(
return p->conserved.mass;
}
+/**
+ * @brief Returns the velocities drifted to the current time of a particle.
+ *
+ * @param p The particle of interest
+ * @param xp The extended data of the particle.
+ * @param dt The time since the last kick.
+ * @param v (return) The velocities at the current time.
+ */
+__attribute__((always_inline)) INLINE static void hydro_get_drifted_velocities(
+ const struct part* restrict p, const struct xpart* xp, float dt,
+ float v[3]) {}
+
/**
* @brief Returns the density of a particle
*
diff --git a/src/hydro/Minimal/hydro.h b/src/hydro/Minimal/hydro.h
index 905412040c..935b5d9684 100644
--- a/src/hydro/Minimal/hydro.h
+++ b/src/hydro/Minimal/hydro.h
@@ -116,6 +116,23 @@ __attribute__((always_inline)) INLINE static float hydro_get_mass(
return p->mass;
}
+/**
+ * @brief Returns the velocities drifted to the current time of a particle.
+ *
+ * @param p The particle of interest
+ * @param xp The extended data of the particle.
+ * @param dt The time since the last kick.
+ * @param v (return) The velocities at the current time.
+ */
+__attribute__((always_inline)) INLINE static void hydro_get_drifted_velocities(
+ const struct part *restrict p, const struct xpart *xp, float dt,
+ float v[3]) {
+
+ v[0] = xp->v_full[0] + p->a_hydro[0] * dt;
+ v[1] = xp->v_full[1] + p->a_hydro[1] * dt;
+ v[2] = xp->v_full[2] + p->a_hydro[2] * dt;
+}
+
/**
* @brief Returns the time derivative of internal energy of a particle
*
diff --git a/src/hydro/PressureEntropy/hydro.h b/src/hydro/PressureEntropy/hydro.h
index 080b796b21..2c10e4fde8 100644
--- a/src/hydro/PressureEntropy/hydro.h
+++ b/src/hydro/PressureEntropy/hydro.h
@@ -106,6 +106,23 @@ __attribute__((always_inline)) INLINE static float hydro_get_mass(
return p->mass;
}
+/**
+ * @brief Returns the velocities drifted to the current time of a particle.
+ *
+ * @param p The particle of interest
+ * @param xp The extended data of the particle.
+ * @param dt The time since the last kick.
+ * @param v (return) The velocities at the current time.
+ */
+__attribute__((always_inline)) INLINE static void hydro_get_drifted_velocities(
+ const struct part *restrict p, const struct xpart *xp, float dt,
+ float v[3]) {
+
+ v[0] = xp->v_full[0] + p->a_hydro[0] * dt;
+ v[1] = xp->v_full[1] + p->a_hydro[1] * dt;
+ v[2] = xp->v_full[2] + p->a_hydro[2] * dt;
+}
+
/**
* @brief Returns the time derivative of internal energy of a particle
*
diff --git a/src/hydro/Shadowswift/hydro.h b/src/hydro/Shadowswift/hydro.h
index eb9938a341..e4ea0d971c 100644
--- a/src/hydro/Shadowswift/hydro.h
+++ b/src/hydro/Shadowswift/hydro.h
@@ -544,6 +544,18 @@ __attribute__((always_inline)) INLINE static float hydro_get_mass(
return p->conserved.mass;
}
+/**
+ * @brief Returns the velocities drifted to the current time of a particle.
+ *
+ * @param p The particle of interest
+ * @param xp The extended data of the particle.
+ * @param dt The time since the last kick.
+ * @param v (return) The velocities at the current time.
+ */
+__attribute__((always_inline)) INLINE static void hydro_get_drifted_velocities(
+ const struct part* restrict p, const struct xpart* xp, float dt,
+ float v[3]) {}
+
/**
* @brief Returns the density of a particle
*
diff --git a/src/statistics.c b/src/statistics.c
index d50143b9eb..f3604a2417 100644
--- a/src/statistics.c
+++ b/src/statistics.c
@@ -135,22 +135,11 @@ void stats_collect_part_mapper(void *map_data, int nr_parts, void *extra_data) {
const integertime_t ti_end = get_integer_time_end(ti_current, p->time_bin);
const float dt = (ti_current - ((ti_begin + ti_end) / 2)) * timeBase;
- /* Get the total acceleration */
- float a_tot[3] = {p->a_hydro[0], p->a_hydro[1], p->a_hydro[2]};
- if (gp != NULL) {
- a_tot[0] += gp->a_grav[0];
- a_tot[1] += gp->a_grav[1];
- a_tot[2] += gp->a_grav[2];
- }
-
- /* Extrapolate velocities to current time */
- const float v[3] = {xp->v_full[0] + a_tot[0] * dt,
- xp->v_full[1] + a_tot[1] * dt,
- xp->v_full[2] + a_tot[2] * dt};
-
- const float m = hydro_get_mass(p);
+ float v[3];
+ hydro_get_drifted_velocities(p, xp, dt, v);
const double x[3] = {p->x[0], p->x[1], p->x[2]};
- const float rho = hydro_get_density(p);
+ const float m = hydro_get_mass(p);
+ const float entropy = hydro_get_entropy(p);
const float u_int = hydro_get_internal_energy(p);
/* Collect mass */
@@ -182,7 +171,7 @@ void stats_collect_part_mapper(void *map_data, int nr_parts, void *extra_data) {
}
/* Collect entropy */
- stats.entropy += m * gas_entropy_from_internal_energy(rho, u_int);
+ stats.entropy += m * entropy;
}
/* Now write back to memory */
--
GitLab