From c33df916c368410e413675fc2979a39d099c2fba Mon Sep 17 00:00:00 2001
From: Bert Vandenbroucke <bert.vandenbroucke@ugent.be>
Date: Tue, 28 Feb 2017 10:49:50 +0000
Subject: [PATCH] Cleaned up Gizmo fixes. Added a hydro_timestep_extra method
and added an empty version for all hydro schemes.
---
src/hydro/Default/hydro.h | 10 +++++++++
src/hydro/Gadget2/hydro.h | 10 +++++++++
src/hydro/Gizmo/hydro.h | 35 ++++++++++++++++++++++++-------
src/hydro/Gizmo/hydro_iact.h | 8 +------
src/hydro/Gizmo/hydro_part.h | 4 ++--
src/hydro/Minimal/hydro.h | 10 +++++++++
src/hydro/PressureEntropy/hydro.h | 10 +++++++++
src/kick.h | 2 +-
src/runner.c | 9 +++++---
9 files changed, 77 insertions(+), 21 deletions(-)
diff --git a/src/hydro/Default/hydro.h b/src/hydro/Default/hydro.h
index bfb5cd1ce3..a614d08c30 100644
--- a/src/hydro/Default/hydro.h
+++ b/src/hydro/Default/hydro.h
@@ -148,6 +148,16 @@ __attribute__((always_inline)) INLINE static float hydro_compute_timestep(
return min(dt_cfl, dt_u_change);
}
+/**
+ * @brief Does some extra hydro operations once the actual physical time step
+ * for the particle is known.
+ *
+ * @param p The particle to act upon.
+ * @param dt Physical time step of the particle during the next step.
+ */
+__attribute__((always_inline)) INLINE static void hydro_timestep_extra(
+ struct part *p, float dt) {}
+
/**
* @brief Prepares a particle for the density calculation.
*
diff --git a/src/hydro/Gadget2/hydro.h b/src/hydro/Gadget2/hydro.h
index 160a2d8b5d..cc7b422ccb 100644
--- a/src/hydro/Gadget2/hydro.h
+++ b/src/hydro/Gadget2/hydro.h
@@ -152,6 +152,16 @@ __attribute__((always_inline)) INLINE static float hydro_compute_timestep(
return dt_cfl;
}
+/**
+ * @brief Does some extra hydro operations once the actual physical time step
+ * for the particle is known.
+ *
+ * @param p The particle to act upon.
+ * @param dt Physical time step of the particle during the next step.
+ */
+__attribute__((always_inline)) INLINE static void hydro_timestep_extra(
+ struct part *p, float dt) {}
+
/**
* @brief Prepares a particle for the density calculation.
*
diff --git a/src/hydro/Gizmo/hydro.h b/src/hydro/Gizmo/hydro.h
index 8da949feb2..c29c76bed0 100644
--- a/src/hydro/Gizmo/hydro.h
+++ b/src/hydro/Gizmo/hydro.h
@@ -40,6 +40,27 @@ __attribute__((always_inline)) INLINE static float hydro_compute_timestep(
return CFL_condition * p->h / fabsf(p->timestepvars.vmax);
}
+/**
+ * @brief Does some extra hydro operations once the actual physical time step
+ * for the particle is known.
+ *
+ * We use this to store the physical time step, since it is used for the flux
+ * exchange during the force loop.
+ *
+ * We also set the active flag of the particle to inactive. It will be set to
+ * active in hydro_init_part, which is called the next time the particle becomes
+ * active.
+ *
+ * @param p The particle to act upon.
+ * @param dt Physical time step of the particle during the next step.
+ */
+__attribute__((always_inline)) INLINE static void hydro_timestep_extra(
+ struct part* p, float dt) {
+
+ p->force.dt = dt;
+ p->force.active = 0;
+}
+
/**
* @brief Initialises the particles for the first time
*
@@ -89,6 +110,9 @@ __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;
+
+ /* Set the active flag to active. */
+ p->force.active = 1;
}
/**
@@ -178,10 +202,7 @@ __attribute__((always_inline)) INLINE static void hydro_end_density(
* @param timeBase Conversion factor between integer time and physical time.
*/
__attribute__((always_inline)) INLINE static void hydro_prepare_force(
- struct part* restrict p, struct xpart* restrict xp, double timeBase) {
-
- /* Set the physical time step */
- p->force.dt = get_timestep(p->time_bin, timeBase); // MATTHIEU 0
+ struct part* restrict p, struct xpart* restrict xp) {
/* Initialize time step criterion variables */
p->timestepvars.vmax = 0.0f;
@@ -375,7 +396,7 @@ __attribute__((always_inline)) INLINE static void hydro_end_force(
* @param half_dt Half the physical time step.
*/
__attribute__((always_inline)) INLINE static void hydro_kick_extra(
- struct part* p, struct xpart* xp, float dt, integertime_t ti_current) {
+ struct part* p, struct xpart* xp, float dt) {
float oldm, oldp[3], anew[3];
const float half_dt = 0.5f * dt; // MATTHIEU
@@ -429,9 +450,9 @@ __attribute__((always_inline)) INLINE static void hydro_kick_extra(
/* Store gravitational acceleration and mass flux for next step */
p->gravity.old_a[0] = anew[0];
+ p->gravity.old_a[1] = anew[1];
p->gravity.old_a[2] = anew[2];
p->gravity.old_mflux[0] = p->gravity.mflux[0];
- p->gravity.old_a[1] = anew[1];
p->gravity.old_mflux[1] = p->gravity.mflux[1];
p->gravity.old_mflux[2] = p->gravity.mflux[2];
}
@@ -444,8 +465,6 @@ __attribute__((always_inline)) INLINE static void hydro_kick_extra(
p->conserved.flux.momentum[1] = 0.0f;
p->conserved.flux.momentum[2] = 0.0f;
p->conserved.flux.energy = 0.0f;
-
- p->force.ti_end = get_integer_time_end(ti_current, p->time_bin);
}
/**
diff --git a/src/hydro/Gizmo/hydro_iact.h b/src/hydro/Gizmo/hydro_iact.h
index a0d8ee4c0b..611f0449a5 100644
--- a/src/hydro/Gizmo/hydro_iact.h
+++ b/src/hydro/Gizmo/hydro_iact.h
@@ -412,13 +412,7 @@ __attribute__((always_inline)) INLINE static void runner_iact_fluxes_common(
==> we update particle j if (MODE IS 1) OR (j IS INACTIVE)
*/
- // MATTHIEU
- // const integertime_t pj_ti_end = 0; // get_integer_time_end(pj->time_bin);
- // const integertime_t pi_ti_end = 0; // get_integer_time_end(pi->time_bin);
- integertime_t pi_ti_end = pi->force.ti_end;
- integertime_t pj_ti_end = pj->force.ti_end;
-
- if (mode == 1 || pj_ti_end > pi_ti_end) {
+ if (mode == 1 || pj->force.active == 0) {
/* Store mass flux */
mflux = dtj * Anorm * totflux[0];
pj->gravity.mflux[0] -= mflux * dx[0];
diff --git a/src/hydro/Gizmo/hydro_part.h b/src/hydro/Gizmo/hydro_part.h
index 88b4daef6f..928011d820 100644
--- a/src/hydro/Gizmo/hydro_part.h
+++ b/src/hydro/Gizmo/hydro_part.h
@@ -178,8 +178,8 @@ struct part {
/* Physical time step of the particle. */
float dt;
- /* Integer end time of the particle's time step. */
- integertime_t ti_end;
+ /* Flag keeping track of whether this is an active or inactive particle. */
+ char active;
/* Actual velocity of the particle. */
float v_full[3];
diff --git a/src/hydro/Minimal/hydro.h b/src/hydro/Minimal/hydro.h
index 20856b7e03..56078a8256 100644
--- a/src/hydro/Minimal/hydro.h
+++ b/src/hydro/Minimal/hydro.h
@@ -165,6 +165,16 @@ __attribute__((always_inline)) INLINE static float hydro_compute_timestep(
return dt_cfl;
}
+/**
+ * @brief Does some extra hydro operations once the actual physical time step
+ * for the particle is known.
+ *
+ * @param p The particle to act upon.
+ * @param dt Physical time step of the particle during the next step.
+ */
+__attribute__((always_inline)) INLINE static void hydro_timestep_extra(
+ struct part *p, float dt) {}
+
/**
* @brief Prepares a particle for the density calculation.
*
diff --git a/src/hydro/PressureEntropy/hydro.h b/src/hydro/PressureEntropy/hydro.h
index f22bb8a13a..20238896f1 100644
--- a/src/hydro/PressureEntropy/hydro.h
+++ b/src/hydro/PressureEntropy/hydro.h
@@ -152,6 +152,16 @@ __attribute__((always_inline)) INLINE static float hydro_compute_timestep(
return dt_cfl;
}
+/**
+ * @brief Does some extra hydro operations once the actual physical time step
+ * for the particle is known.
+ *
+ * @param p The particle to act upon.
+ * @param dt Physical time step of the particle during the next step.
+ */
+__attribute__((always_inline)) INLINE static void hydro_timestep_extra(
+ struct part *p, float dt) {}
+
/**
* @brief Prepares a particle for the density calculation.
*
diff --git a/src/kick.h b/src/kick.h
index 1c6f5fb68e..7ccea7d269 100644
--- a/src/kick.h
+++ b/src/kick.h
@@ -107,7 +107,7 @@ __attribute__((always_inline)) INLINE static void kick_part(
}
/* Extra kick work */
- hydro_kick_extra(p, xp, dt, ti_end);
+ hydro_kick_extra(p, xp, dt);
if (p->gpart != NULL) gravity_kick_extra(p->gpart, dt);
}
diff --git a/src/runner.c b/src/runner.c
index ad078f204e..f0970ae701 100644
--- a/src/runner.c
+++ b/src/runner.c
@@ -599,8 +599,6 @@ void runner_do_ghost(struct runner *r, struct cell *c, int timer) {
target_wcount - e->hydro_properties->delta_neighbours;
const int max_smoothing_iter = e->hydro_properties->max_smoothing_iterations;
int redo = 0, count = 0;
- const double timeBase = e->timeBase;
- integertime_t ti_current = e->ti_current;
TIMER_TIC;
@@ -681,7 +679,7 @@ void runner_do_ghost(struct runner *r, struct cell *c, int timer) {
/* As of here, particle force variables will be set. */
/* Compute variables required for the force loop */
- hydro_prepare_force(p, xp, ti_current, timeBase);
+ hydro_prepare_force(p, xp);
/* The particle force values are now set. Do _NOT_
try to read any particle density variables! */
@@ -1098,6 +1096,7 @@ void runner_do_timestep(struct runner *r, struct cell *c, int timer) {
struct xpart *restrict xparts = c->xparts;
struct gpart *restrict gparts = c->gparts;
struct spart *restrict sparts = c->sparts;
+ const double timeBase = e->timeBase;
TIMER_TIC;
@@ -1133,6 +1132,10 @@ void runner_do_timestep(struct runner *r, struct cell *c, int timer) {
p->time_bin = get_time_bin(ti_new_step);
if (p->gpart != NULL) p->gpart->time_bin = get_time_bin(ti_new_step);
+ /* Tell the particle what the new physical time step is */
+ float dt = get_timestep(p->time_bin, timeBase);
+ hydro_timestep_extra(p, dt);
+
/* Number of updated particles */
updated++;
if (p->gpart != NULL) g_updated++;
--
GitLab