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

Added an extra loop over neighbours between density and force to compute the GIZMO gradients.

parent 3265764b
......@@ -704,6 +704,9 @@ void cell_clean_links(struct cell *c, void *data) {
c->density = NULL;
c->nr_density = 0;
c->gradient = NULL;
c->nr_gradient = 0;
c->force = NULL;
c->nr_force = 0;
}
......
......@@ -118,11 +118,11 @@ struct cell {
int sortsize, gsortsize;
/* The tasks computing this cell's density. */
struct link *density, *force, *grav;
int nr_density, nr_force, nr_grav;
struct link *density, *gradient, *force, *grav;
int nr_density, nr_gradient, nr_force, nr_grav;
/* The hierarchical tasks. */
struct task *ghost, *init, *drift, *kick;
struct task *extra_ghost, *ghost, *init, *drift, *kick;
/* Task receiving data. */
struct task *recv_xv, *recv_rho, *recv_ti;
......
......@@ -225,6 +225,12 @@ void engine_make_hydro_hierarchical_tasks(struct engine *e, struct cell *c,
/* Generate the ghost task. */
c->ghost = scheduler_addtask(s, task_type_ghost, task_subtype_none, 0, 0,
c, NULL, 0);
#ifdef EXTRA_HYDRO_LOOP
/* Generate the extra ghost task. */
c->extra_ghost = scheduler_addtask(s, task_type_extra_ghost,
task_subtype_none, 0, 0, c, NULL, 0);
#endif
}
}
......@@ -1477,19 +1483,46 @@ void engine_link_gravity_tasks(struct engine *e) {
}
}
#ifdef EXTRA_HYDRO_LOOP
/**
* @brief Creates the dependency network for the hydro tasks of a given cell.
*
* @param sched The #scheduler.
* @param density The density task to link.
* @param gradient The gradient task to link.
* @param force The force task to link.
* @param c The cell.
*/
static inline void engine_make_hydro_loops_dependencies(struct scheduler *sched,
struct task *density,
struct task *gradient,
struct task *force,
struct cell *c) {
/* init --> density loop --> ghost --> gradient loop --> extra_ghost */
/* extra_ghost --> force loop --> kick */
scheduler_addunlock(sched, c->super->init, density);
scheduler_addunlock(sched, density, c->super->ghost);
scheduler_addunlock(sched, c->super->ghost, gradient);
scheduler_addunlock(sched, gradient, c->super->extra_ghost);
scheduler_addunlock(sched, c->super->extra_ghost, force);
scheduler_addunlock(sched, force, c->super->kick);
}
#else
/**
* @brief Creates the dependency network for the hydro tasks of a given cell.
*
* @param sched The #scheduler.
* @param density The density task to link.
* @param force The force task to link.
* @param c The cell.
*/
static inline void engine_make_hydro_loops_dependencies(struct scheduler *sched,
struct task *density,
struct task *force,
struct cell *c) {
/* init --> density loop --> ghost --> force loop --> kick */
scheduler_addunlock(sched, c->super->init, density);
scheduler_addunlock(sched, density, c->super->ghost);
......@@ -1497,6 +1530,7 @@ static inline void engine_make_hydro_loops_dependencies(struct scheduler *sched,
scheduler_addunlock(sched, force, c->super->kick);
}
#endif
/**
* @brief Duplicates the first hydro loop and construct all the
* dependencies for the hydro part
......@@ -1526,6 +1560,24 @@ void engine_make_extra_hydroloop_tasks(struct engine *e) {
/* Self-interaction? */
if (t->type == task_type_self && t->subtype == task_subtype_density) {
#ifdef EXTRA_HYDRO_LOOP
/* Start by constructing the task for the second and third hydro loop */
struct task *t2 = scheduler_addtask(
sched, task_type_self, task_subtype_gradient, 0, 0, t->ci, NULL, 0);
struct task *t3 = scheduler_addtask(
sched, task_type_self, task_subtype_force, 0, 0, t->ci, NULL, 0);
/* Add the link between the new loops and the cell */
t->ci->gradient = engine_addlink(e, t->ci->gradient, t2);
atomic_inc(&t->ci->nr_gradient);
t->ci->force = engine_addlink(e, t->ci->force, t3);
atomic_inc(&t->ci->nr_force);
/* Now, build all the dependencies for the hydro */
engine_make_hydro_loops_dependencies(sched, t, t2, t3, t->ci);
#else
/* Start by constructing the task for the second hydro loop */
struct task *t2 = scheduler_addtask(
sched, task_type_self, task_subtype_force, 0, 0, t->ci, NULL, 0);
......@@ -1536,11 +1588,40 @@ void engine_make_extra_hydroloop_tasks(struct engine *e) {
/* Now, build all the dependencies for the hydro */
engine_make_hydro_loops_dependencies(sched, t, t2, t->ci);
#endif
}
/* Otherwise, pair interaction? */
else if (t->type == task_type_pair && t->subtype == task_subtype_density) {
#ifdef EXTRA_HYDRO_LOOP
/* Start by constructing the task for the second and third hydro loop */
struct task *t2 = scheduler_addtask(
sched, task_type_pair, task_subtype_gradient, 0, 0, t->ci, t->cj, 0);
struct task *t3 = scheduler_addtask(
sched, task_type_pair, task_subtype_force, 0, 0, t->ci, t->cj, 0);
/* Add the link between the new loop and both cells */
t->ci->gradient = engine_addlink(e, t->ci->gradient, t2);
atomic_inc(&t->ci->nr_gradient);
t->cj->gradient = engine_addlink(e, t->cj->gradient, t2);
atomic_inc(&t->cj->nr_gradient);
t->ci->force = engine_addlink(e, t->ci->force, t3);
atomic_inc(&t->ci->nr_force);
t->cj->force = engine_addlink(e, t->cj->force, t3);
atomic_inc(&t->cj->nr_force);
/* Now, build all the dependencies for the hydro for the cells */
/* that are local and are not descendant of the same super-cells */
if (t->ci->nodeID == nodeID) {
engine_make_hydro_loops_dependencies(sched, t, t2, t3, t->ci);
}
if (t->cj->nodeID == nodeID && t->ci->super != t->cj->super) {
engine_make_hydro_loops_dependencies(sched, t, t2, t3, t->cj);
}
#else
/* Start by constructing the task for the second hydro loop */
struct task *t2 = scheduler_addtask(
sched, task_type_pair, task_subtype_force, 0, 0, t->ci, t->cj, 0);
......@@ -1559,12 +1640,38 @@ void engine_make_extra_hydroloop_tasks(struct engine *e) {
if (t->cj->nodeID == nodeID && t->ci->super != t->cj->super) {
engine_make_hydro_loops_dependencies(sched, t, t2, t->cj);
}
#endif
}
/* Otherwise, sub-self interaction? */
else if (t->type == task_type_sub_self &&
t->subtype == task_subtype_density) {
#ifdef EXTRA_HYDRO_LOOP
/* Start by constructing the task for the second and third hydro loop */
struct task *t2 =
scheduler_addtask(sched, task_type_sub_self, task_subtype_gradient,
t->flags, 0, t->ci, t->cj, 0);
struct task *t3 =
scheduler_addtask(sched, task_type_sub_self, task_subtype_force,
t->flags, 0, t->ci, t->cj, 0);
/* Add the link between the new loop and the cell */
t->ci->gradient = engine_addlink(e, t->ci->gradient, t2);
atomic_inc(&t->ci->nr_gradient);
t->ci->force = engine_addlink(e, t->ci->force, t3);
atomic_inc(&t->ci->nr_force);
/* Now, build all the dependencies for the hydro for the cells */
/* that are local and are not descendant of the same super-cells */
if (t->ci->nodeID == nodeID) {
engine_make_hydro_loops_dependencies(sched, t, t2, t3, t->ci);
}
#else
/* Start by constructing the task for the second hydro loop */
struct task *t2 =
scheduler_addtask(sched, task_type_sub_self, task_subtype_force,
......@@ -1579,12 +1686,43 @@ void engine_make_extra_hydroloop_tasks(struct engine *e) {
if (t->ci->nodeID == nodeID) {
engine_make_hydro_loops_dependencies(sched, t, t2, t->ci);
}
#endif
}
/* Otherwise, sub-pair interaction? */
else if (t->type == task_type_sub_pair &&
t->subtype == task_subtype_density) {
#ifdef EXTRA_HYDRO_LOOP
/* Start by constructing the task for the second and third hydro loop */
struct task *t2 =
scheduler_addtask(sched, task_type_sub_pair, task_subtype_gradient,
t->flags, 0, t->ci, t->cj, 0);
struct task *t3 =
scheduler_addtask(sched, task_type_sub_pair, task_subtype_force,
t->flags, 0, t->ci, t->cj, 0);
/* Add the link between the new loop and both cells */
t->ci->gradient = engine_addlink(e, t->ci->gradient, t2);
atomic_inc(&t->ci->nr_gradient);
t->cj->gradient = engine_addlink(e, t->cj->gradient, t2);
atomic_inc(&t->cj->nr_gradient);
t->ci->force = engine_addlink(e, t->ci->force, t3);
atomic_inc(&t->ci->nr_force);
t->cj->force = engine_addlink(e, t->cj->force, t3);
atomic_inc(&t->cj->nr_force);
/* Now, build all the dependencies for the hydro for the cells */
/* that are local and are not descendant of the same super-cells */
if (t->ci->nodeID == nodeID) {
engine_make_hydro_loops_dependencies(sched, t, t2, t3, t->ci);
}
if (t->cj->nodeID == nodeID && t->ci->super != t->cj->super) {
engine_make_hydro_loops_dependencies(sched, t, t2, t3, t->cj);
}
#else
/* Start by constructing the task for the second hydro loop */
struct task *t2 =
scheduler_addtask(sched, task_type_sub_pair, task_subtype_force,
......@@ -1604,8 +1742,8 @@ void engine_make_extra_hydroloop_tasks(struct engine *e) {
if (t->cj->nodeID == nodeID && t->ci->super != t->cj->super) {
engine_make_hydro_loops_dependencies(sched, t, t2, t->cj);
}
#endif
}
/* External gravity tasks should depend on init and unlock the kick */
else if (t->type == task_type_grav_external) {
scheduler_addunlock(sched, t->ci->init, t);
......@@ -2577,8 +2715,12 @@ void engine_step(struct engine *e) {
mask |= 1 << task_type_sub_self;
mask |= 1 << task_type_sub_pair;
mask |= 1 << task_type_ghost;
mask |=
1 << task_type_extra_ghost; /* Adding unnecessary things to the mask
does not harm */
submask |= 1 << task_subtype_density;
submask |= 1 << task_subtype_gradient;
submask |= 1 << task_subtype_force;
}
......
......@@ -94,15 +94,14 @@ __attribute__((always_inline)) INLINE static void runner_iact_nonsym_density(
hydro_gradients_density_loop(pi, pj, wi_dx, 0.0f, dx, r, 0);
}
__attribute__((always_inline)) INLINE static void runner_iact_hydro_loop2(
__attribute__((always_inline)) INLINE static void runner_iact_gradient(
float r2, float *dx, float hi, float hj, struct part *pi, struct part *pj) {
hydro_gradients_gradient_loop(r2, dx, hi, hj, pi, pj, 1);
}
__attribute__((always_inline)) INLINE static void
runner_iact_nonsym_hydro_loop2(float r2, float *dx, float hi, float hj,
struct part *pi, struct part *pj) {
__attribute__((always_inline)) INLINE static void runner_iact_nonsym_gradient(
float r2, float *dx, float hi, float hj, struct part *pi, struct part *pj) {
hydro_gradients_gradient_loop(r2, dx, hi, hj, pi, pj, 0);
}
......
......@@ -47,6 +47,7 @@
#include "./hydro/Default/hydro_part.h"
#elif defined(GIZMO_SPH)
#include "./hydro/Gizmo/hydro_part.h"
#define EXTRA_HYDRO_LOOP
#else
#error "Invalid choice of SPH variant"
#endif
......
......@@ -82,6 +82,13 @@ const char runner_flip[27] = {1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0,
#define FUNCTION density
#include "runner_doiact.h"
#ifdef EXTRA_HYDRO_LOOP
/* Import the gradient loop functions. */
#undef FUNCTION
#define FUNCTION gradient
#include "runner_doiact.h"
#endif
/* Import the force loop functions. */
#undef FUNCTION
#define FUNCTION force
......@@ -421,8 +428,11 @@ void runner_do_init(struct runner *r, struct cell *c, int timer) {
if (timer) TIMER_TOC(timer_init);
}
void runner_do_extra_ghost(struct runner *r, struct cell *c) {}
/**
* @brief Intermediate task between density and force
* @brief Intermediate task after the density to check that the smoothing
* lengths are correct.
*
* @param r The runner thread.
* @param c The cell.
......@@ -1074,8 +1084,11 @@ void *runner_main(void *data) {
/* Different types of tasks... */
switch (t->type) {
case task_type_self:
if (t->subtype == task_subtype_density)
runner_doself1_density(r, ci);
if (t->subtype == task_subtype_density) runner_doself1_density(r, ci);
#ifdef EXTRA_HYDRO_LOOP
else if (t->subtype == task_subtype_gradient)
runner_doself1_gradient(r, ci);
#endif
else if (t->subtype == task_subtype_force)
runner_doself2_force(r, ci);
else if (t->subtype == task_subtype_grav)
......@@ -1086,6 +1099,10 @@ void *runner_main(void *data) {
case task_type_pair:
if (t->subtype == task_subtype_density)
runner_dopair1_density(r, ci, cj);
#ifdef EXTRA_HYDRO_LOOP
else if (t->subtype == task_subtype_gradient)
runner_dopair1_gradient(r, ci, cj);
#endif
else if (t->subtype == task_subtype_force)
runner_dopair2_force(r, ci, cj);
else if (t->subtype == task_subtype_grav)
......@@ -1099,6 +1116,10 @@ void *runner_main(void *data) {
case task_type_sub_self:
if (t->subtype == task_subtype_density)
runner_dosub_self1_density(r, ci, 1);
#ifdef EXTRA_HYDRO_LOOP
else if (t->subtype == task_subtype_gradient)
runner_dosub_self1_gradient(r, ci, 1);
#endif
else if (t->subtype == task_subtype_force)
runner_dosub_self2_force(r, ci, 1);
else if (t->subtype == task_subtype_grav)
......@@ -1109,6 +1130,10 @@ void *runner_main(void *data) {
case task_type_sub_pair:
if (t->subtype == task_subtype_density)
runner_dosub_pair1_density(r, ci, cj, t->flags, 1);
#ifdef EXTRA_HYDRO_LOOP
else if (t->subtype == task_subtype_gradient)
runner_dosub_pair1_gradient(r, ci, cj, t->flags, 1);
#endif
else if (t->subtype == task_subtype_force)
runner_dosub_pair2_force(r, ci, cj, t->flags, 1);
else if (t->subtype == task_subtype_grav)
......@@ -1122,6 +1147,9 @@ void *runner_main(void *data) {
case task_type_ghost:
runner_do_ghost(r, ci);
break;
case task_type_extra_ghost:
runner_do_extra_ghost(r, ci);
break;
case task_type_drift:
runner_do_drift(r, ci, 1);
break;
......
......@@ -351,14 +351,17 @@ void space_regrid(struct space *s, double cell_max, int verbose) {
s->cells[k].sorts = NULL;
s->cells[k].nr_tasks = 0;
s->cells[k].nr_density = 0;
s->cells[k].nr_gradient = 0;
s->cells[k].nr_force = 0;
s->cells[k].density = NULL;
s->cells[k].gradient = NULL;
s->cells[k].force = NULL;
s->cells[k].dx_max = 0.0f;
s->cells[k].sorted = 0;
s->cells[k].count = 0;
s->cells[k].gcount = 0;
s->cells[k].init = NULL;
s->cells[k].extra_ghost = NULL;
s->cells[k].ghost = NULL;
s->cells[k].drift = NULL;
s->cells[k].kick = NULL;
......
......@@ -48,14 +48,14 @@
/* Task type names. */
const char *taskID_names[task_type_count] = {
"none", "sort", "self", "pair", "sub_self",
"sub_pair", "init", "ghost", "drift", "kick",
"kick_fixdt", "send", "recv", "grav_gather_m", "grav_fft",
"grav_mm", "grav_up", "grav_external", "part_sort", "gpart_sort",
"split_cell", "rewait"};
const char *subtaskID_names[task_subtype_count] = {"none", "density", "force",
"grav", "tend"};
"none", "sort", "self", "pair", "sub_self",
"sub_pair", "init", "ghost", "extra_ghost", "drift",
"kick", "kick_fixdt", "send", "recv", "grav_gather_m",
"grav_fft", "grav_mm", "grav_up", "grav_external", "part_sort",
"gpart_sort", "split_cell", "rewait"};
const char *subtaskID_names[task_subtype_count] = {
"none", "density", "gradient", "force", "grav", "tend"};
/**
* @brief Computes the overlap between the parts array of two given cells.
......
......@@ -41,6 +41,7 @@ enum task_types {
task_type_sub_pair,
task_type_init,
task_type_ghost,
task_type_extra_ghost,
task_type_drift,
task_type_kick,
task_type_kick_fixdt,
......@@ -64,6 +65,7 @@ extern const char *taskID_names[];
enum task_subtypes {
task_subtype_none = 0,
task_subtype_density,
task_subtype_gradient,
task_subtype_force,
task_subtype_grav,
task_subtype_tend,
......
......@@ -37,11 +37,11 @@ ticks timers[timer_count];
* To reset all timers, use the mask #timers_mask_all.
*/
void timers_reset(unsigned int mask) {
void timers_reset(unsigned long long mask) {
int k;
/* Loop over the timers and set the masked ones to zero. */
for (k = 0; k < timer_count; k++)
if (mask & (1 << k)) timers[k] = 0;
if (mask & (1ull << k)) timers[k] = 0;
}
......@@ -36,17 +36,21 @@ enum {
timer_kick,
timer_dosort,
timer_doself_density,
timer_doself_gradient,
timer_doself_force,
timer_doself_grav_pp,
timer_dopair_density,
timer_dopair_gradient,
timer_dopair_force,
timer_dopair_grav_pm,
timer_dopair_grav_pp,
timer_dograv_external,
timer_dosub_self_density,
timer_dosub_self_gradient,
timer_dosub_self_force,
timer_dosub_self_grav,
timer_dosub_pair_density,
timer_dosub_pair_gradient,
timer_dosub_pair_force,
timer_dosub_pair_grav,
timer_dopair_subset,
......@@ -64,7 +68,7 @@ enum {
extern ticks timers[timer_count];
/* Mask for all timers. */
#define timers_mask_all ((1 << timer_count) - 1)
#define timers_mask_all ((1ull << timer_count) - 1)
/* Define the timer macros. */
#ifdef TIMER
......@@ -74,7 +78,7 @@ extern ticks timers[timer_count];
#define TIMER_TOC(t) timers_toc(t, tic)
#define TIMER_TIC2 ticks tic2 = getticks();
#define TIMER_TOC2(t) timers_toc(t, tic2)
INLINE static ticks timers_toc(int t, ticks tic) {
INLINE static ticks timers_toc(unsigned int t, ticks tic) {
ticks d = (getticks() - tic);
atomic_add(&timers[t], d);
return d;
......@@ -87,6 +91,6 @@ INLINE static ticks timers_toc(int t, ticks tic) {
#endif
/* Function prototypes. */
void timers_reset(unsigned int mask);
void timers_reset(unsigned long long mask);
#endif /* SWIFT_TIMERS_H */
Markdown is supported
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