Commit dd7f2e51 authored by Matthieu Schaller's avatar Matthieu Schaller
Browse files

Make the particles carry their own softening to speed up the cache construction

parent 9ba6efdb
......@@ -4717,6 +4717,7 @@ void cell_drift_gpart(struct cell *c, const struct engine *e, int force) {
const integertime_t ti_old_gpart = c->grav.ti_old_part;
const integertime_t ti_current = e->ti_current;
struct gpart *const gparts = c->grav.parts;
const struct gravity_props *grav_props = e->gravity_properties;
/* Drift irrespective of cell flags? */
force = (force || cell_get_flag(c, cell_flag_do_grav_drift));
......@@ -4778,7 +4779,7 @@ void cell_drift_gpart(struct cell *c, const struct engine *e, int force) {
if (gpart_is_inhibited(gp, e)) continue;
/* Drift... */
drift_gpart(gp, dt_drift, ti_old_gpart, ti_current);
drift_gpart(gp, dt_drift, ti_old_gpart, ti_current, grav_props);
#ifdef SWIFT_DEBUG_CHECKS
/* Make sure the particle does not drift by more than a box length. */
......
......@@ -43,7 +43,7 @@
*/
__attribute__((always_inline)) INLINE static void drift_gpart(
struct gpart *restrict gp, double dt_drift, integertime_t ti_old,
integertime_t ti_current) {
integertime_t ti_current, const struct gravity_props *grav_props) {
#ifdef SWIFT_DEBUG_CHECKS
if (gp->ti_drift != ti_old)
......@@ -60,6 +60,8 @@ __attribute__((always_inline)) INLINE static void drift_gpart(
gp->x[0] += gp->v_full[0] * dt_drift;
gp->x[1] += gp->v_full[1] * dt_drift;
gp->x[2] += gp->v_full[2] * dt_drift;
gravity_predict_extra(gp, grav_props);
}
/**
......
......@@ -212,6 +212,17 @@ __attribute__((always_inline)) INLINE static void gravity_end_force(
#endif
}
/**
* @brief Update the #gpart after a drift step.
*
* This is typically used to update the softening lengths.
*
* @param gp The particle to act upon
* @param grav_props The global properties of the gravity calculation.
*/
__attribute__((always_inline)) INLINE static void gravity_predict_extra(
struct gpart* gp, const struct gravity_props* grav_props) {}
/**
* @brief Kick the additional variables
*
......
......@@ -50,23 +50,7 @@ __attribute__((always_inline)) INLINE static float gravity_get_mass(
__attribute__((always_inline)) INLINE static float gravity_get_softening(
const struct gpart* gp, const struct gravity_props* restrict grav_props) {
switch (gp->type) {
case swift_type_dark_matter:
return grav_props->epsilon_DM_cur;
case swift_type_stars:
return grav_props->epsilon_baryon_cur;
case swift_type_gas:
return grav_props->epsilon_baryon_cur;
case swift_type_black_hole:
return grav_props->epsilon_baryon_cur;
case swift_type_dark_matter_background:
return grav_props->epsilon_background_fac * cbrtf(gp->mass);
default:
#ifdef SWIFT_DEBUG_CHECKS
error("Invalid gpart type!");
#endif
return 0.f;
}
return gp->epsilon;
}
/**
......@@ -239,6 +223,41 @@ __attribute__((always_inline)) INLINE static void gravity_end_force(
#endif
}
/**
* @brief Update the #gpart after a drift step.
*
* This is typically used to update the softening lengths.
*
* @param gp The particle to act upon
* @param grav_props The global properties of the gravity calculation.
*/
__attribute__((always_inline)) INLINE static void gravity_predict_extra(
struct gpart* gp, const struct gravity_props* grav_props) {
switch (gp->type) {
case swift_type_dark_matter:
gp->epsilon = grav_props->epsilon_DM_cur;
break;
case swift_type_stars:
gp->epsilon = grav_props->epsilon_baryon_cur;
break;
case swift_type_gas:
gp->epsilon = grav_props->epsilon_baryon_cur;
break;
case swift_type_black_hole:
gp->epsilon = grav_props->epsilon_baryon_cur;
break;
case swift_type_dark_matter_background:
gp->epsilon = grav_props->epsilon_background_fac * cbrtf(gp->mass);
break;
default:
#ifdef SWIFT_DEBUG_CHECKS
error("Invalid gpart type!");
#endif
break;
}
}
/**
* @brief Kick the additional variables
*
......@@ -272,6 +291,29 @@ __attribute__((always_inline)) INLINE static void gravity_first_init_gpart(
gp->time_bin = 0;
gp->old_a_grav_norm = 0.f;
switch (gp->type) {
case swift_type_dark_matter:
gp->epsilon = grav_props->epsilon_DM_cur;
break;
case swift_type_stars:
gp->epsilon = grav_props->epsilon_baryon_cur;
break;
case swift_type_gas:
gp->epsilon = grav_props->epsilon_baryon_cur;
break;
case swift_type_black_hole:
gp->epsilon = grav_props->epsilon_baryon_cur;
break;
case swift_type_dark_matter_background:
gp->epsilon = grav_props->epsilon_background_fac * cbrtf(gp->mass);
break;
default:
#ifdef SWIFT_DEBUG_CHECKS
error("Invalid gpart type!");
#endif
break;
}
gravity_init_gpart(gp);
}
......
......@@ -46,6 +46,9 @@ struct gpart {
/*! Norm of the acceleration at the previous step. */
float old_a_grav_norm;
/*! Current co-moving spline softening of the particle */
float epsilon;
/*! Particle FoF properties (group ID, group size, ...) */
struct fof_gpart_data fof_data;
......
......@@ -195,6 +195,17 @@ __attribute__((always_inline)) INLINE static void gravity_end_force(
#endif
}
/**
* @brief Update the #gpart after a drift step.
*
* This is typically used to update the softening lengths.
*
* @param gp The particle to act upon
* @param grav_props The global properties of the gravity calculation.
*/
__attribute__((always_inline)) INLINE static void gravity_predict_extra(
struct gpart* gp, const struct gravity_props* grav_props) {}
/**
* @brief Kick the additional variables
*
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment