Commit 28deacb4 authored by Matthieu Schaller's avatar Matthieu Schaller
Browse files

Make sure the inhibited particles do not take part in the gravity calculation.

parent eec60d76
......@@ -618,8 +618,10 @@ void cell_set_super_mapper(void *map_data, int num_elements, void *extra_data);
int cell_has_tasks(struct cell *c);
void cell_remove_part(const struct engine *e, struct cell *c, struct part *p,
struct xpart *xp);
void cell_remove_gpart(const struct engine *e, struct cell *c, struct gpart *gp);
void cell_remove_spart(const struct engine *e, struct cell *c, struct spart *sp);
void cell_remove_gpart(const struct engine *e, struct cell *c,
struct gpart *gp);
void cell_remove_spart(const struct engine *e, struct cell *c,
struct spart *sp);
int cell_can_use_pair_mm(const struct cell *ci, const struct cell *cj,
const struct engine *e, const struct space *s);
int cell_can_use_pair_mm_rebuild(const struct cell *ci, const struct cell *cj,
......
......@@ -208,12 +208,20 @@ __attribute__((always_inline)) INLINE static void gravity_cache_populate(
/* Fill the input caches */
for (int i = 0; i < gcount; ++i) {
x[i] = (float)(gparts[i].x[0] - shift[0]);
y[i] = (float)(gparts[i].x[1] - shift[1]);
z[i] = (float)(gparts[i].x[2] - shift[2]);
epsilon[i] = gravity_get_softening(&gparts[i], grav_props);
m[i] = gparts[i].mass;
active[i] = (int)(gparts[i].time_bin <= max_active_bin);
/* Make a dummy particle out of the inhibted ones */
if (gparts[i].time_bin == time_bin_inhibited) {
m[i] = 0.f;
active[i] = 0;
} else {
m[i] = gparts[i].mass;
active[i] = (int)(gparts[i].time_bin <= max_active_bin);
}
/* Distance to the CoM of the other cell. */
float dx = x[i] - CoM[0];
......@@ -294,8 +302,15 @@ gravity_cache_populate_no_mpole(const timebin_t max_active_bin,
y[i] = (float)(gparts[i].x[1] - shift[1]);
z[i] = (float)(gparts[i].x[2] - shift[2]);
epsilon[i] = gravity_get_softening(&gparts[i], grav_props);
m[i] = gparts[i].mass;
active[i] = (int)(gparts[i].time_bin <= max_active_bin);
/* Make a dummy particle out of the inhibted ones */
if (gparts[i].time_bin == time_bin_inhibited) {
m[i] = 0.f;
active[i] = 0;
} else {
m[i] = gparts[i].mass;
active[i] = (int)(gparts[i].time_bin <= max_active_bin);
}
}
#ifdef SWIFT_DEBUG_CHECKS
......
......@@ -121,6 +121,10 @@ static INLINE void runner_do_grav_down(struct runner *r, struct cell *c,
if (c->grav.multipole->pot.ti_init != e->ti_current)
error("c->field tensor not initialised");
/* Check that we are not updated an inhibited particle */
if (gp->time_bin == time_bin_inhibited)
error("Updating an inhibited particle!");
/* Check that the particle was initialised */
if (gp->initialised == 0)
error("Adding forces to an un-initialised gpart.");
......@@ -229,6 +233,14 @@ static INLINE void runner_dopair_grav_pp_full(
if (pjd < gcount_j && gparts_j[pjd].ti_drift != e->ti_current)
error("gpj not drifted to current time");
/* Check that we are not updated an inhibited particle */
if (gparts_i[pid].time_bin == time_bin_inhibited)
error("Updating an inhibited particle!");
/* Check that the particle we interact with was not inhibited */
if (gparts_j[pjd].time_bin == time_bin_inhibited && mass_j != 0.f)
error("Inhibited particle used as gravity source.");
/* Check that the particle was initialised */
if (gparts_i[pid].initialised == 0)
error("Adding forces to an un-initialised gpart.");
......@@ -359,6 +371,14 @@ static INLINE void runner_dopair_grav_pp_truncated(
if (pjd < gcount_j && gparts_j[pjd].ti_drift != e->ti_current)
error("gpj not drifted to current time");
/* Check that we are not updated an inhibited particle */
if (gparts_i[pid].time_bin == time_bin_inhibited)
error("Updating an inhibited particle!");
/* Check that the particle we interact with was not inhibited */
if (gparts_j[pjd].time_bin == time_bin_inhibited && mass_j != 0.f)
error("Inhibited particle used as gravity source.");
/* Check that the particle was initialised */
if (gparts_i[pid].initialised == 0)
error("Adding forces to an un-initialised gpart.");
......@@ -450,6 +470,10 @@ static INLINE void runner_dopair_grav_pm_full(
if (gparts_i[pid].ti_drift != e->ti_current)
error("gpi not drifted to current time");
/* Check that we are not updated an inhibited particle */
if (gparts_i[pid].time_bin == time_bin_inhibited)
error("Updating an inhibited particle!");
/* Check that the particle was initialised */
if (gparts_i[pid].initialised == 0)
error("Adding forces to an un-initialised gpart.");
......@@ -579,6 +603,10 @@ static INLINE void runner_dopair_grav_pm_truncated(
if (gparts_i[pid].ti_drift != e->ti_current)
error("gpi not drifted to current time");
/* Check that we are not updated an inhibited particle */
if (gparts_i[pid].time_bin == time_bin_inhibited)
error("Updating an inhibited particle!");
/* Check that the particle was initialised */
if (gparts_i[pid].initialised == 0)
error("Adding forces to an un-initialised gpart.");
......@@ -928,6 +956,14 @@ static INLINE void runner_doself_grav_pp_full(
if (pjd < gcount && gparts[pjd].ti_drift != e->ti_current)
error("gpj not drifted to current time");
/* Check that we are not updated an inhibited particle */
if (gparts[pid].time_bin == time_bin_inhibited)
error("Updating an inhibited particle!");
/* Check that the particle we interact with was not inhibited */
if (gparts[pjd].time_bin == time_bin_inhibited && mass_j != 0.f)
error("Inhibited particle used as gravity source.");
/* Check that the particle was initialised */
if (gparts[pid].initialised == 0)
error("Adding forces to an un-initialised gpart.");
......@@ -1042,6 +1078,14 @@ static INLINE void runner_doself_grav_pp_truncated(
if (pjd < gcount && gparts[pjd].ti_drift != e->ti_current)
error("gpj not drifted to current time");
/* Check that we are not updated an inhibited particle */
if (gparts[pid].time_bin == time_bin_inhibited)
error("Updating an inhibited particle!");
/* Check that the particle we interact with was not inhibited */
if (gparts[pjd].time_bin == time_bin_inhibited && mass_j != 0.f)
error("Inhibited particle used as gravity source.");
/* Check that the particle was initialised */
if (gparts[pid].initialised == 0)
error("Adding forces to an un-initialised gpart.");
......
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