Skip to content
Snippets Groups Projects
Commit 006ecb95 authored by Matthieu Schaller's avatar Matthieu Schaller
Browse files

Apply the same change to the timestep limiter loops

parent 1c275faf
No related branches found
No related tags found
1 merge request!1353Draft: Subtask speedup - Still requires work
Hide whitespace changes
Inline Side-by-side
src/runner_doiact_functions_limiter.h
+ 111
− 22
View file @ 006ecb95
...@@ -60,9 +60,15 @@ void DOPAIR1_NAIVE(struct runner *r, struct cell *restrict ci, ...@@ -60,9 +60,15 @@ void DOPAIR1_NAIVE(struct runner *r, struct cell *restrict ci,
if (ci->dmin != cj->dmin) error("Cells of different size!"); if (ci->dmin != cj->dmin) error("Cells of different size!");
#endif #endif
/* Get the depth limits (if any) */
const char min_depth = limit_max_h ? ci->depth : 0;
const char max_depth = limit_min_h ? ci->depth : CHAR_MAX;
#ifdef SWIFT_DEBUG_CHECKS
/* Get the limits in h (if any) */ /* Get the limits in h (if any) */
const float h_min = limit_min_h ? ci->dmin * 0.5 * (1. / kernel_gamma) : 0.; const float h_min = limit_min_h ? ci->h_min_allowed : 0.;
const float h_max = limit_max_h ? ci->dmin * (1. / kernel_gamma) : FLT_MAX; const float h_max = limit_max_h ? ci->h_max_allowed : FLT_MAX;
#endif
/* Get the relative distance between the pairs, wrapping. */ /* Get the relative distance between the pairs, wrapping. */
double shift[3] = {0.0, 0.0, 0.0}; double shift[3] = {0.0, 0.0, 0.0};
...@@ -83,6 +89,7 @@ void DOPAIR1_NAIVE(struct runner *r, struct cell *restrict ci, ...@@ -83,6 +89,7 @@ void DOPAIR1_NAIVE(struct runner *r, struct cell *restrict ci,
if (part_is_inhibited(pi, e)) continue; if (part_is_inhibited(pi, e)) continue;
const int pi_active = part_is_starting(pi, e); const int pi_active = part_is_starting(pi, e);
const char depth_i = pi->depth_h;
const float hi = pi->h; const float hi = pi->h;
const float hig2 = hi * hi * kernel_gamma2; const float hig2 = hi * hi * kernel_gamma2;
const float pix[3] = {(float)(pi->x[0] - (cj->loc[0] + shift[0])), const float pix[3] = {(float)(pi->x[0] - (cj->loc[0] + shift[0])),
...@@ -94,6 +101,7 @@ void DOPAIR1_NAIVE(struct runner *r, struct cell *restrict ci, ...@@ -94,6 +101,7 @@ void DOPAIR1_NAIVE(struct runner *r, struct cell *restrict ci,
/* Get a pointer to the jth particle. */ /* Get a pointer to the jth particle. */
struct part *restrict pj = &parts_j[pjd]; struct part *restrict pj = &parts_j[pjd];
const char depth_j = pj->depth_h;
/* Skip inhibited particles. */ /* Skip inhibited particles. */
if (part_is_inhibited(pj, e)) continue; if (part_is_inhibited(pj, e)) continue;
...@@ -117,16 +125,26 @@ void DOPAIR1_NAIVE(struct runner *r, struct cell *restrict ci, ...@@ -117,16 +125,26 @@ void DOPAIR1_NAIVE(struct runner *r, struct cell *restrict ci,
error("Particle pj not drifted to current time"); error("Particle pj not drifted to current time");
#endif #endif
const int doi = pi_active && (r2 < hig2) && (hi >= h_min) && (hi < h_max); const int doi = pi_active && (r2 < hig2) && (depth_i >= min_depth) &&
const int doj = pj_active && (r2 < hjg2) && (hj >= h_min) && (hj < h_max); (depth_i <= max_depth);
const int doj = pj_active && (r2 < hjg2) && (depth_j >= min_depth) &&
(depth_j <= max_depth);
/* Hit or miss? */ /* Hit or miss? */
if (doi) { if (doi) {
#ifdef SWIFT_DEBUG_CHECKS
if (hi < h_min || hi >= h_max) error("Inappropriate h for this level!");
#endif
IACT_NONSYM(r2, dx, hi, hj, pi, pj, a, H); IACT_NONSYM(r2, dx, hi, hj, pi, pj, a, H);
} }
if (doj) { if (doj) {
#ifdef SWIFT_DEBUG_CHECKS
if (hj < h_min || hj >= h_max) error("Inappropriate h for this level!");
#endif
dx[0] = -dx[0]; dx[0] = -dx[0];
dx[1] = -dx[1]; dx[1] = -dx[1];
dx[2] = -dx[2]; dx[2] = -dx[2];
...@@ -166,9 +184,15 @@ void DOSELF1_NAIVE(struct runner *r, const struct cell *c, ...@@ -166,9 +184,15 @@ void DOSELF1_NAIVE(struct runner *r, const struct cell *c,
const int count = c->hydro.count; const int count = c->hydro.count;
struct part *parts = c->hydro.parts; struct part *parts = c->hydro.parts;
/* Get the depth limits (if any) */
const char min_depth = limit_max_h ? c->depth : 0;
const char max_depth = limit_min_h ? c->depth : CHAR_MAX;
#ifdef SWIFT_DEBUG_CHECKS
/* Get the limits in h (if any) */ /* Get the limits in h (if any) */
const float h_min = limit_min_h ? c->dmin * 0.5 * (1. / kernel_gamma) : 0.; const float h_min = limit_min_h ? c->h_min_allowed : 0.;
const float h_max = limit_max_h ? c->dmin * (1. / kernel_gamma) : FLT_MAX; const float h_max = limit_max_h ? c->h_max_allowed : FLT_MAX;
#endif
/* Loop over the parts in ci. */ /* Loop over the parts in ci. */
for (int pid = 0; pid < count; pid++) { for (int pid = 0; pid < count; pid++) {
...@@ -180,6 +204,7 @@ void DOSELF1_NAIVE(struct runner *r, const struct cell *c, ...@@ -180,6 +204,7 @@ void DOSELF1_NAIVE(struct runner *r, const struct cell *c,
if (part_is_inhibited(pi, e)) continue; if (part_is_inhibited(pi, e)) continue;
const int pi_active = part_is_starting(pi, e); const int pi_active = part_is_starting(pi, e);
const char depth_i = pi->depth_h;
const float hi = pi->h; const float hi = pi->h;
const float hig2 = hi * hi * kernel_gamma2; const float hig2 = hi * hi * kernel_gamma2;
const float pix[3] = {(float)(pi->x[0] - c->loc[0]), const float pix[3] = {(float)(pi->x[0] - c->loc[0]),
...@@ -198,6 +223,7 @@ void DOSELF1_NAIVE(struct runner *r, const struct cell *c, ...@@ -198,6 +223,7 @@ void DOSELF1_NAIVE(struct runner *r, const struct cell *c,
const float hj = pj->h; const float hj = pj->h;
const float hjg2 = hj * hj * kernel_gamma2; const float hjg2 = hj * hj * kernel_gamma2;
const int pj_active = part_is_starting(pj, e); const int pj_active = part_is_starting(pj, e);
const char depth_j = pj->depth_h;
/* Compute the pairwise distance. */ /* Compute the pairwise distance. */
const float pjx[3] = {(float)(pj->x[0] - c->loc[0]), const float pjx[3] = {(float)(pj->x[0] - c->loc[0]),
...@@ -206,8 +232,10 @@ void DOSELF1_NAIVE(struct runner *r, const struct cell *c, ...@@ -206,8 +232,10 @@ void DOSELF1_NAIVE(struct runner *r, const struct cell *c,
float dx[3] = {pix[0] - pjx[0], pix[1] - pjx[1], pix[2] - pjx[2]}; float dx[3] = {pix[0] - pjx[0], pix[1] - pjx[1], pix[2] - pjx[2]};
const float r2 = dx[0] * dx[0] + dx[1] * dx[1] + dx[2] * dx[2]; const float r2 = dx[0] * dx[0] + dx[1] * dx[1] + dx[2] * dx[2];
const int doi = pi_active && (r2 < hig2) && (hi >= h_min) && (hi < h_max); const int doi = pi_active && (r2 < hig2) && (depth_i >= min_depth) &&
const int doj = pj_active && (r2 < hjg2) && (hj >= h_min) && (hj < h_max); (depth_i <= max_depth);
const int doj = pj_active && (r2 < hjg2) && (depth_j >= min_depth) &&
(depth_j <= max_depth);
#ifdef SWIFT_DEBUG_CHECKS #ifdef SWIFT_DEBUG_CHECKS
/* Check that particles have been drifted to the current time */ /* Check that particles have been drifted to the current time */
...@@ -220,12 +248,25 @@ void DOSELF1_NAIVE(struct runner *r, const struct cell *c, ...@@ -220,12 +248,25 @@ void DOSELF1_NAIVE(struct runner *r, const struct cell *c,
/* Hit or miss? */ /* Hit or miss? */
if (doi && doj) { if (doi && doj) {
#ifdef SWIFT_DEBUG_CHECKS
if (hi < h_min || hi >= h_max) error("Inappropriate h for this level!");
if (hj < h_min || hj >= h_max) error("Inappropriate h for this level!");
#endif
IACT(r2, dx, hi, hj, pi, pj, a, H); IACT(r2, dx, hi, hj, pi, pj, a, H);
} else if (doi) { } else if (doi) {
#ifdef SWIFT_DEBUG_CHECKS
if (hi < h_min || hi >= h_max) error("Inappropriate h for this level!");
#endif
IACT_NONSYM(r2, dx, hi, hj, pi, pj, a, H); IACT_NONSYM(r2, dx, hi, hj, pi, pj, a, H);
} else if (doj) { } else if (doj) {
#ifdef SWIFT_DEBUG_CHECKS
if (hj < h_min || hj >= h_max) error("Inappropriate h for this level!");
#endif
dx[0] = -dx[0]; dx[0] = -dx[0];
dx[1] = -dx[1]; dx[1] = -dx[1];
dx[2] = -dx[2]; dx[2] = -dx[2];
...@@ -277,9 +318,15 @@ void DOPAIR1(struct runner *r, const struct cell *restrict ci, ...@@ -277,9 +318,15 @@ void DOPAIR1(struct runner *r, const struct cell *restrict ci,
2. * max(ci->hydro.dx_max_part, cj->hydro.dx_max_part); 2. * max(ci->hydro.dx_max_part, cj->hydro.dx_max_part);
#endif /* SWIFT_DEBUG_CHECKS */ #endif /* SWIFT_DEBUG_CHECKS */
/* Get the depth limits (if any) */
const char min_depth = limit_max_h ? ci->depth : 0;
const char max_depth = limit_min_h ? ci->depth : CHAR_MAX;
#ifdef SWIFT_DEBUG_CHECKS
/* Get the limits in h (if any) */ /* Get the limits in h (if any) */
const float h_min = limit_min_h ? ci->dmin * 0.5 * (1. / kernel_gamma) : 0.; const float h_min = limit_min_h ? ci->h_min_allowed : 0.;
const float h_max = limit_max_h ? ci->dmin * (1. / kernel_gamma) : FLT_MAX; #endif
const float h_max = limit_max_h ? ci->h_max_allowed : FLT_MAX;
/* Get some other useful values. */ /* Get some other useful values. */
const double hi_max = const double hi_max =
...@@ -306,6 +353,7 @@ void DOPAIR1(struct runner *r, const struct cell *restrict ci, ...@@ -306,6 +353,7 @@ void DOPAIR1(struct runner *r, const struct cell *restrict ci,
/* Get a hold of the ith part in ci. */ /* Get a hold of the ith part in ci. */
struct part *restrict pi = &parts_i[sort_i[pid].i]; struct part *restrict pi = &parts_i[sort_i[pid].i];
const char depth_i = pi->depth_h;
const float hi = pi->h; const float hi = pi->h;
/* Skip inactive particles */ /* Skip inactive particles */
...@@ -317,8 +365,8 @@ void DOPAIR1(struct runner *r, const struct cell *restrict ci, ...@@ -317,8 +365,8 @@ void DOPAIR1(struct runner *r, const struct cell *restrict ci,
#endif #endif
/* Skip particles not in the range of h we care about */ /* Skip particles not in the range of h we care about */
if (hi >= h_max) continue; if (depth_i < min_depth) continue;
if (hi < h_min) continue; if (depth_i > max_depth) continue;
/* Is there anything we need to interact with ? */ /* Is there anything we need to interact with ? */
const double di = sort_i[pid].d + hi * kernel_gamma + dx_max - rshift; const double di = sort_i[pid].d + hi * kernel_gamma + dx_max - rshift;
...@@ -385,6 +433,11 @@ void DOPAIR1(struct runner *r, const struct cell *restrict ci, ...@@ -385,6 +433,11 @@ void DOPAIR1(struct runner *r, const struct cell *restrict ci,
/* Hit or miss? */ /* Hit or miss? */
if (r2 < hig2) { if (r2 < hig2) {
#ifdef SWIFT_DEBUG_CHECKS
if (hi < h_min || hi >= h_max)
error("Inappropriate h for this level!");
#endif
IACT_NONSYM(r2, dx, hi, hj, pi, pj, a, H); IACT_NONSYM(r2, dx, hi, hj, pi, pj, a, H);
} }
} /* loop over the parts in cj. */ } /* loop over the parts in cj. */
...@@ -399,6 +452,7 @@ void DOPAIR1(struct runner *r, const struct cell *restrict ci, ...@@ -399,6 +452,7 @@ void DOPAIR1(struct runner *r, const struct cell *restrict ci,
/* Get a hold of the jth part in cj. */ /* Get a hold of the jth part in cj. */
struct part *pj = &parts_j[sort_j[pjd].i]; struct part *pj = &parts_j[sort_j[pjd].i];
const char depth_j = pj->depth_h;
const float hj = pj->h; const float hj = pj->h;
/* Skip inactive particles */ /* Skip inactive particles */
...@@ -410,8 +464,8 @@ void DOPAIR1(struct runner *r, const struct cell *restrict ci, ...@@ -410,8 +464,8 @@ void DOPAIR1(struct runner *r, const struct cell *restrict ci,
#endif #endif
/* Skip particles not in the range of h we care about */ /* Skip particles not in the range of h we care about */
if (hj >= h_max) continue; if (depth_j < min_depth) continue;
if (hj < h_min) continue; if (depth_j > max_depth) continue;
/* Is there anything we need to interact with ? */ /* Is there anything we need to interact with ? */
const double dj = sort_j[pjd].d - hj * kernel_gamma - dx_max + rshift; const double dj = sort_j[pjd].d - hj * kernel_gamma - dx_max + rshift;
...@@ -478,6 +532,11 @@ void DOPAIR1(struct runner *r, const struct cell *restrict ci, ...@@ -478,6 +532,11 @@ void DOPAIR1(struct runner *r, const struct cell *restrict ci,
/* Hit or miss? */ /* Hit or miss? */
if (r2 < hjg2) { if (r2 < hjg2) {
#ifdef SWIFT_DEBUG_CHECKS
if (hj < h_min || hj >= h_max)
error("Inappropriate h for this level!");
#endif
IACT_NONSYM(r2, dx, hj, hi, pj, pi, a, H); IACT_NONSYM(r2, dx, hj, hi, pj, pi, a, H);
} }
} /* loop over the parts in ci. */ } /* loop over the parts in ci. */
...@@ -549,9 +608,15 @@ void DOSELF1(struct runner *r, const struct cell *c, const int limit_min_h, ...@@ -549,9 +608,15 @@ void DOSELF1(struct runner *r, const struct cell *c, const int limit_min_h,
struct part *parts = c->hydro.parts; struct part *parts = c->hydro.parts;
const int count = c->hydro.count; const int count = c->hydro.count;
/* Get the depth limits (if any) */
const char min_depth = limit_max_h ? c->depth : 0;
const char max_depth = limit_min_h ? c->depth : CHAR_MAX;
#ifdef SWIFT_DEBUG_CHECKS
/* Get the limits in h (if any) */ /* Get the limits in h (if any) */
const float h_min = limit_min_h ? c->dmin * 0.5 * (1. / kernel_gamma) : 0.; const float h_min = limit_min_h ? c->h_min_allowed : 0.;
const float h_max = limit_max_h ? c->dmin * (1. / kernel_gamma) : FLT_MAX; const float h_max = limit_max_h ? c->h_max_allowed : FLT_MAX;
#endif
/* Set up a list of the particles for which we want to compute interactions */ /* Set up a list of the particles for which we want to compute interactions */
int *indt = NULL; int *indt = NULL;
...@@ -561,8 +626,9 @@ void DOSELF1(struct runner *r, const struct cell *c, const int limit_min_h, ...@@ -561,8 +626,9 @@ void DOSELF1(struct runner *r, const struct cell *c, const int limit_min_h,
error("Failed to allocate indt."); error("Failed to allocate indt.");
for (int k = 0; k < count; k++) { for (int k = 0; k < count; k++) {
const struct part *p = &parts[k]; const struct part *p = &parts[k];
const float h = p->h; const char depth = p->depth_h;
if (part_is_starting(&parts[k], e) && (h >= h_min) && (h < h_max)) { if (part_is_starting(&parts[k], e) && (depth >= min_depth) &&
(depth <= max_depth)) {
indt[countdt] = k; indt[countdt] = k;
countdt += 1; countdt += 1;
} }
...@@ -580,6 +646,7 @@ void DOSELF1(struct runner *r, const struct cell *c, const int limit_min_h, ...@@ -580,6 +646,7 @@ void DOSELF1(struct runner *r, const struct cell *c, const int limit_min_h,
/* Get a pointer to the ith particle. */ /* Get a pointer to the ith particle. */
struct part *restrict pi = &parts[pid]; struct part *restrict pi = &parts[pid];
const char depth_i = pi->depth_h;
/* Skip inhibited particles. */ /* Skip inhibited particles. */
if (part_is_inhibited(pi, e)) continue; if (part_is_inhibited(pi, e)) continue;
...@@ -590,8 +657,8 @@ void DOSELF1(struct runner *r, const struct cell *c, const int limit_min_h, ...@@ -590,8 +657,8 @@ void DOSELF1(struct runner *r, const struct cell *c, const int limit_min_h,
const float hig2 = hi * hi * kernel_gamma2; const float hig2 = hi * hi * kernel_gamma2;
/* Is the ith particle active and in the range of h we care about? */ /* Is the ith particle active and in the range of h we care about? */
const int update_i = const int update_i = part_is_starting(pi, e) && (depth_i >= min_depth) &&
part_is_starting(pi, e) && (hi >= h_min) && (hi < h_max); (depth_i <= max_depth);
/* If false then it can only act as a neighbour of others */ /* If false then it can only act as a neighbour of others */
if (!update_i) { if (!update_i) {
...@@ -623,6 +690,11 @@ void DOSELF1(struct runner *r, const struct cell *c, const int limit_min_h, ...@@ -623,6 +690,11 @@ void DOSELF1(struct runner *r, const struct cell *c, const int limit_min_h,
/* Hit or miss? */ /* Hit or miss? */
if (r2 < hjg2) { if (r2 < hjg2) {
#ifdef SWIFT_DEBUG_CHECKS
if (hj < h_min || hj >= h_max)
error("Inappropriate h for this level!");
#endif
IACT_NONSYM(r2, dx, hj, hi, pj, pi, a, H); IACT_NONSYM(r2, dx, hj, hi, pj, pi, a, H);
} }
} /* loop over all the particles we want to update. */ } /* loop over all the particles we want to update. */
...@@ -643,6 +715,7 @@ void DOSELF1(struct runner *r, const struct cell *c, const int limit_min_h, ...@@ -643,6 +715,7 @@ void DOSELF1(struct runner *r, const struct cell *c, const int limit_min_h,
/* Get a pointer to the jth particle (by construction pi != pj). */ /* Get a pointer to the jth particle (by construction pi != pj). */
struct part *restrict pj = &parts[pjd]; struct part *restrict pj = &parts[pjd];
const char depth_j = pj->depth_h;
/* Skip inhibited particles. */ /* Skip inhibited particles. */
if (part_is_inhibited(pj, e)) continue; if (part_is_inhibited(pj, e)) continue;
...@@ -675,22 +748,38 @@ void DOSELF1(struct runner *r, const struct cell *c, const int limit_min_h, ...@@ -675,22 +748,38 @@ void DOSELF1(struct runner *r, const struct cell *c, const int limit_min_h,
* -> Check whether it is active * -> Check whether it is active
* -> Check whether it is in the right range of h * -> Check whether it is in the right range of h
* -> Check the distance to pi */ * -> Check the distance to pi */
const int doj = (part_is_starting(pj, e)) && (hj >= h_min) && const int doj = (part_is_starting(pj, e)) && (depth_j >= min_depth) &&
(hj < h_max) && (r2 < hjg2); (depth_j <= max_depth) && (r2 < hjg2);
/* Hit or miss? */ /* Hit or miss? */
if (doi && doj) { if (doi && doj) {
#ifdef SWIFT_DEBUG_CHECKS
if (hi < h_min || hi >= h_max)
error("Inappropriate h for this level!");
if (hj < h_min || hj >= h_max)
error("Inappropriate h for this level!");
#endif
/* Update both pi and pj */ /* Update both pi and pj */
IACT(r2, dx, hi, hj, pi, pj, a, H); IACT(r2, dx, hi, hj, pi, pj, a, H);
} else if (doi) { } else if (doi) {
#ifdef SWIFT_DEBUG_CHECKS
if (hi < h_min || hi >= h_max)
error("Inappropriate h for this level!");
#endif
/* Update only pi */ /* Update only pi */
IACT_NONSYM(r2, dx, hi, hj, pi, pj, a, H); IACT_NONSYM(r2, dx, hi, hj, pi, pj, a, H);
} else if (doj) { } else if (doj) {
#ifdef SWIFT_DEBUG_CHECKS
if (hj < h_min || hj >= h_max)
error("Inappropriate h for this level!");
#endif
/* Update only pj */ /* Update only pj */
dx[0] = -dx[0]; dx[0] = -dx[0];
......
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment