diff --git a/src/engine.c b/src/engine.c
index 799b5d40247b6f38482aa72bebe0f11f120a4469..800b44f8dff89b043952aa4f4507ba033f112d85 100644
--- a/src/engine.c
+++ b/src/engine.c
@@ -2455,6 +2455,7 @@ void engine_init(struct engine *e, struct space *s,
e->ti_current = 0;
e->timeStep = 0.;
e->timeBase = 0.;
+ e->timeBase_inv = 0.;
e->timeFirstSnapshot =
parser_get_param_double(params, "Snapshots:time_first");
e->deltaTimeSnapshot =
@@ -2596,6 +2597,7 @@ void engine_init(struct engine *e, struct space *s,
/* Deal with timestep */
e->timeBase = (e->timeEnd - e->timeBegin) / max_nr_timesteps;
+ e->timeBase_inv = 1.0 / e->timeBase;
e->ti_current = 0;
/* Fixed time-step case */
diff --git a/src/engine.h b/src/engine.h
index 9ef7d57599d30aad5e8000e64148812493299d23..f2219202be1559ae3e41c7e8e6cb41e31b3e17dd 100644
--- a/src/engine.h
+++ b/src/engine.h
@@ -132,6 +132,7 @@ struct engine {
/* Time base */
double timeBase;
+ double timeBase_inv;
/* Snapshot information */
double timeFirstSnapshot;
diff --git a/src/runner.c b/src/runner.c
index 721a3214e847b2b9400f0916ed804ff2ba3ea759..aeb871b4e7046d87bfad9cd491363293498999d7 100644
--- a/src/runner.c
+++ b/src/runner.c
@@ -1079,20 +1079,14 @@ void runner_do_kick_fixdt(struct runner *r, struct cell *c, int timer) {
*/
void runner_do_kick(struct runner *r, struct cell *c, int timer) {
- const double global_dt_min = r->e->dt_min;
- const double global_dt_max = r->e->dt_max;
+ const struct engine *e = r->e;
+ const double timeBase = e->timeBase;
const int ti_current = r->e->ti_current;
- const double timeBase = r->e->timeBase;
- const double timeBase_inv = 1.0 / r->e->timeBase;
const int count = c->count;
const int gcount = c->gcount;
struct part *const parts = c->parts;
struct xpart *const xparts = c->xparts;
struct gpart *const gparts = c->gparts;
- const struct external_potential *potential = r->e->external_potential;
- const struct hydro_props *hydro_properties = r->e->hydro_properties;
- const struct phys_const *constants = r->e->physical_constants;
- const float ln_max_h_change = hydro_properties->log_max_h_change;
int updated = 0, g_updated = 0;
int ti_end_min = max_nr_timesteps, ti_end_max = 0;
@@ -1124,9 +1118,7 @@ void runner_do_kick(struct runner *r, struct cell *c, int timer) {
gravity_end_force(gp);
/* Compute the next timestep */
- const int new_dti =
- get_gpart_timestep(gp, potential, constants, global_dt_min,
- global_dt_max, timeBase_inv);
+ const int new_dti = get_gpart_timestep(gp, e);
/* Now we have a time step, proceed with the kick */
kick_gpart(gp, new_dti, timeBase);
@@ -1161,9 +1153,7 @@ void runner_do_kick(struct runner *r, struct cell *c, int timer) {
if (p->gpart != NULL) gravity_end_force(p->gpart);
/* Compute the next timestep (hydro condition) */
- const int new_dti = get_part_timestep(
- p, xp, hydro_properties, potential, constants, global_dt_min,
- global_dt_max, timeBase_inv, ln_max_h_change);
+ const int new_dti = get_part_timestep(p, xp, e);
/* Now we have a time step, proceed with the kick */
kick_part(p, xp, new_dti, timeBase);
diff --git a/src/timestep.h b/src/timestep.h
index 794e4de47ea35f57c6108bcd521f78eaf33b4a99..c0697c82cdb084e0ccdadfc8c7b19249d9cfd7c7 100644
--- a/src/timestep.h
+++ b/src/timestep.h
@@ -26,26 +26,14 @@
#include "const.h"
#include "debug.h"
-__attribute__((always_inline)) INLINE static int get_gpart_timestep(
- const struct gpart *gp, const struct external_potential *potential,
- const struct phys_const *constants, double global_dt_min,
- double global_dt_max, double timeBase_inv) {
-
- const float new_dt_external =
- gravity_compute_timestep_external(potential, constants, gp);
- const float new_dt_self = gravity_compute_timestep_self(constants, gp);
-
- float new_dt = fminf(new_dt_external, new_dt_self);
-
- /* Limit timestep within the allowed range */
- new_dt = fminf(new_dt, global_dt_max);
- new_dt = fmaxf(new_dt, global_dt_min);
+__attribute__((always_inline)) INLINE static int get_integer_timestep(
+ float new_dt, int ti_begin, int ti_end, double timeBase_inv) {
/* Convert to integer time */
- int new_dti = new_dt * timeBase_inv;
+ int new_dti = (int)(new_dt * timeBase_inv);
/* Recover the current timestep */
- const int current_dti = gp->ti_end - gp->ti_begin;
+ const int current_dti = ti_end - ti_begin;
/* Limit timestep increase */
if (current_dti > 0) new_dti = min(new_dti, 2 * current_dti);
@@ -58,30 +46,47 @@ __attribute__((always_inline)) INLINE static int get_gpart_timestep(
/* Make sure we are allowed to increase the timestep size */
if (new_dti > current_dti) {
- if ((max_nr_timesteps - gp->ti_end) % new_dti > 0) new_dti = current_dti;
+ if ((max_nr_timesteps - ti_end) % new_dti > 0) new_dti = current_dti;
}
return new_dti;
}
+__attribute__((always_inline)) INLINE static int get_gpart_timestep(
+ const struct gpart *gp, const struct engine *e) {
+
+ const float new_dt_external = gravity_compute_timestep_external(
+ e->external_potential, e->physical_constants, gp);
+ const float new_dt_self =
+ gravity_compute_timestep_self(e->physical_constants, gp);
+
+ float new_dt = fminf(new_dt_external, new_dt_self);
+
+ /* Limit timestep within the allowed range */
+ new_dt = fminf(new_dt, e->dt_max);
+ new_dt = fmaxf(new_dt, e->dt_min);
+
+ /* Convert to integer time */
+ const int new_dti =
+ get_integer_timestep(new_dt, gp->ti_begin, gp->ti_end, e->timeBase_inv);
+
+ return new_dti;
+}
+
__attribute__((always_inline)) INLINE static int get_part_timestep(
- const struct part *p, const struct xpart *xp,
- const struct hydro_props *hydro_properties,
- const struct external_potential *potential,
- const struct phys_const *constants, double global_dt_min,
- double global_dt_max, double timeBase_inv, double ln_max_h_change) {
+ const struct part *p, const struct xpart *xp, const struct engine *e) {
/* Compute the next timestep (hydro condition) */
- const float new_dt_hydro = hydro_compute_timestep(p, xp, hydro_properties);
+ const float new_dt_hydro = hydro_compute_timestep(p, xp, e->hydro_properties);
/* Compute the next timestep (gravity condition) */
float new_dt_grav = FLT_MAX;
if (p->gpart != NULL) {
- const float new_dt_external =
- gravity_compute_timestep_external(potential, constants, p->gpart);
+ const float new_dt_external = gravity_compute_timestep_external(
+ e->external_potential, e->physical_constants, p->gpart);
const float new_dt_self =
- gravity_compute_timestep_self(constants, p->gpart);
+ gravity_compute_timestep_self(e->physical_constants, p->gpart);
new_dt_grav = fminf(new_dt_external, new_dt_self);
}
@@ -91,33 +96,19 @@ __attribute__((always_inline)) INLINE static int get_part_timestep(
/* Limit change in h */
const float dt_h_change =
- (p->h_dt != 0.0f) ? fabsf(ln_max_h_change * p->h / p->h_dt) : FLT_MAX;
+ (p->h_dt != 0.0f)
+ ? fabsf(e->hydro_properties->log_max_h_change * p->h / p->h_dt)
+ : FLT_MAX;
new_dt = fminf(new_dt, dt_h_change);
/* Limit timestep within the allowed range */
- new_dt = fminf(new_dt, global_dt_max);
- new_dt = fmaxf(new_dt, global_dt_min);
+ new_dt = fminf(new_dt, e->dt_max);
+ new_dt = fmaxf(new_dt, e->dt_min);
/* Convert to integer time */
- int new_dti = new_dt * timeBase_inv;
-
- /* Recover the current timestep */
- const int current_dti = p->ti_end - p->ti_begin;
-
- /* Limit timestep increase */
- if (current_dti > 0) new_dti = min(new_dti, 2 * current_dti);
-
- /* Put this timestep on the time line */
- int dti_timeline = max_nr_timesteps;
- while (new_dti < dti_timeline) dti_timeline /= 2;
-
- new_dti = dti_timeline;
-
- /* Make sure we are allowed to increase the timestep size */
- if (new_dti > current_dti) {
- if ((max_nr_timesteps - p->ti_end) % new_dti > 0) new_dti = current_dti;
- }
+ const int new_dti =
+ get_integer_timestep(new_dt, p->ti_begin, p->ti_end, e->timeBase_inv);
return new_dti;
}