diff --git a/src/cell.c b/src/cell.c index 5b7941634bd0444a178b9c2e0ad2ff825e3ef9b9..2b34bfcada6fbe06672c0080a3e5b383aae2db8e 100644 --- a/src/cell.c +++ b/src/cell.c @@ -721,7 +721,8 @@ void cell_check_drift_point(struct cell *c, void *data) { /* for (int i = 0; i < c->count; ++i) */ /* if (c->parts[i].ti_old != ti_current) */ /* error( */ - /* "Particle in an incorrect time-zone! part->ti_old=%d c->ti_old=%d " */ + /* "Particle in an incorrect time-zone! part->ti_old=%d c->ti_old=%d " + */ /* "ti_current=%d", */ /* c->parts[i].ti_old, c->ti_old, ti_current); */ } @@ -1009,11 +1010,11 @@ void cell_drift(struct cell *c, const struct engine *e) { const double dt = (ti_current - ti_old) * timeBase; float dx_max = 0.f, dx2_max = 0.f, h_max = 0.f; - //message("DRFIT ! ti_old=%d ti_current=%d", ti_old, ti_current); + // message("DRFIT ! ti_old=%d ti_current=%d", ti_old, ti_current); /* Check that we are actually going to move forward. */ if (ti_current < ti_old) error("Attempt to drift to the past"); - //if (ti_current == ti_old) return; + // if (ti_current == ti_old) return; /* Are we not in a leaf ? */ if (c->split) { @@ -1057,7 +1058,7 @@ void cell_drift(struct cell *c, const struct engine *e) { /* Drift... */ drift_part(p, xp, dt, timeBase, ti_old, ti_current); - //p->ti_old = ti_current; + // p->ti_old = ti_current; /* Compute (square of) motion since last cell construction */ const float dx2 = xp->x_diff[0] * xp->x_diff[0] + diff --git a/src/engine.c b/src/engine.c index 0650547117e5dd8239ede4f654226813663ea0b4..08811e3502df023b2935a2345dfef350f9d4df65 100644 --- a/src/engine.c +++ b/src/engine.c @@ -142,6 +142,11 @@ void engine_make_hierarchical_tasks(struct engine *e, struct cell *c) { c->kick = scheduler_addtask(s, task_type_kick, task_subtype_none, 0, 0, c, NULL, 0); + c->drift = scheduler_addtask(s, task_type_drift, task_subtype_none, 0, 0, + c, NULL, 0); + + scheduler_addunlock(s, c->drift, c->init); + /* Generate the ghost task. */ if (is_hydro) c->ghost = scheduler_addtask(s, task_type_ghost, task_subtype_none, 0, @@ -1987,10 +1992,10 @@ void engine_marktasks_mapper(void *map_data, int num_elements, else if (t->type == task_type_self || t->type == task_type_sub_self) { /* Local pointers. */ - const struct cell *ci = t->ci; + // const struct cell *ci = t->ci; /* Activate the drift */ - if (ci->drift) scheduler_activate(s, ci->drift); + // if (ci->drift) scheduler_activate(s, ci->drift); } /* Pair? */ @@ -2001,8 +2006,8 @@ void engine_marktasks_mapper(void *map_data, int num_elements, const struct cell *cj = t->cj; /* Activate the drift on both sides */ - if (ci->drift) scheduler_activate(s, ci->drift); - if (cj->drift) scheduler_activate(s, cj->drift); + // if (ci->drift) scheduler_activate(s, ci->drift); + // if (cj->drift) scheduler_activate(s, cj->drift); /* Too much particle movement? */ if (t->tight && @@ -2100,6 +2105,11 @@ void engine_marktasks_mapper(void *map_data, int num_elements, if (t->ci->ti_end_min <= ti_end) scheduler_activate(s, t); } + /* Drift? */ + else if (t->type == task_type_drift) { + if (t->ci->ti_end_min <= ti_end) scheduler_activate(s, t); + } + /* Init? */ else if (t->type == task_type_init) { if (t->ci->ti_end_min <= ti_end) scheduler_activate(s, t); @@ -2482,7 +2492,8 @@ void engine_skip_force_and_kick(struct engine *e) { /* Skip everything that updates the particles */ if (t->subtype == task_subtype_force || t->type == task_type_kick || - t->type == task_type_cooling || t->type == task_type_sourceterms) + t->type == task_type_cooling || t->type == task_type_sourceterms || + t->type == task_type_drift) t->skip = 1; } } diff --git a/src/hydro/Gadget2/hydro.h b/src/hydro/Gadget2/hydro.h index a1321da89a0cb24591b0c489bf5afe38b4745fe8..c681db18a7b2d4fdafe7e46c784450b5ce256eca 100644 --- a/src/hydro/Gadget2/hydro.h +++ b/src/hydro/Gadget2/hydro.h @@ -212,7 +212,7 @@ __attribute__((always_inline)) INLINE static float hydro_compute_timestep( __attribute__((always_inline)) INLINE static void hydro_first_init_part( struct part *restrict p, struct xpart *restrict xp) { - //p->ti_old = 0; + // p->ti_old = 0; p->ti_begin = 0; p->ti_end = 0; xp->v_full[0] = p->v[0]; diff --git a/src/hydro/Gadget2/hydro_part.h b/src/hydro/Gadget2/hydro_part.h index ee8114b8eb7aa62123bcb11627965a4862da9511..61158a71d518d7843dd6ea04e21385b1f3d2dc29 100644 --- a/src/hydro/Gadget2/hydro_part.h +++ b/src/hydro/Gadget2/hydro_part.h @@ -55,7 +55,7 @@ struct part { /* Pointer to corresponding gravity part. */ struct gpart* gpart; - + /* Particle position. */ double x[3]; @@ -77,8 +77,8 @@ struct part { /* Particle time of end of time-step. */ int ti_end; - //int ti_old; - + // int ti_old; + /* Particle density. */ float rho; @@ -88,7 +88,7 @@ struct part { /* Entropy time derivative */ float entropy_dt; - //union { + union { struct { @@ -130,7 +130,7 @@ struct part { float h_dt; } force; - //}; + }; } SWIFT_STRUCT_ALIGN; diff --git a/src/runner.c b/src/runner.c index 975068f3486cf1b7211a2f1fa04e366042d4b68e..35148d506f615e0cd67673ba0e87c5d7ab09e979 100644 --- a/src/runner.c +++ b/src/runner.c @@ -736,11 +736,11 @@ void runner_do_ghost(struct runner *r, struct cell *c, int timer) { if (count) { message("Smoothing length failed to converge on %i particles.", count); - for(int i=0; i<count; ++i) { - struct part *restrict p = &parts[pid[i]]; + for (int i = 0; i < count; ++i) { + struct part *restrict p = &parts[pid[i]]; struct xpart *restrict xp = &xparts[pid[i]]; - printParticle_single(p, xp); + printParticle_single(p, xp); } } diff --git a/src/runner_doiact.h b/src/runner_doiact.h index bfc7bcf40465a9793e55579efee3479e30460d35..8b860a013e7f31d24521fe5edf2cbcd5f0660088 100644 --- a/src/runner_doiact.h +++ b/src/runner_doiact.h @@ -721,8 +721,8 @@ void DOPAIR1(struct runner *r, struct cell *ci, struct cell *cj) { /* Anything to do here? */ if (!cell_is_active(ci, e) && !cell_is_active(cj, e)) return; - if(!cell_is_drifted(ci, e)) cell_drift(ci, e); - if(!cell_is_drifted(cj, e)) cell_drift(cj, e); + if (!cell_is_drifted(ci, e)) cell_drift(ci, e); + if (!cell_is_drifted(cj, e)) cell_drift(cj, e); /* Get the sort ID. */ double shift[3] = {0.0, 0.0, 0.0}; @@ -762,8 +762,6 @@ void DOPAIR1(struct runner *r, struct cell *ci, struct cell *cj) { const double di = sort_i[pid].d + hi * kernel_gamma + dx_max - rshift; if (di < dj_min) continue; - - double pix[3]; for (int k = 0; k < 3; k++) pix[k] = pi->x[k] - shift[k]; const float hig2 = hi * hi * kernel_gamma2; @@ -774,7 +772,6 @@ void DOPAIR1(struct runner *r, struct cell *ci, struct cell *cj) { /* Get a pointer to the jth particle. */ struct part *restrict pj = &parts_j[sort_j[pjd].i]; - /* Compute the pairwise distance. */ float r2 = 0.0f; float dx[3]; @@ -827,8 +824,6 @@ void DOPAIR1(struct runner *r, struct cell *ci, struct cell *cj) { const double dj = sort_j[pjd].d - hj * kernel_gamma - dx_max - rshift; if (dj > di_max) continue; - - double pjx[3]; for (int k = 0; k < 3; k++) pjx[k] = pj->x[k] + shift[k]; const float hjg2 = hj * hj * kernel_gamma2; @@ -839,8 +834,6 @@ void DOPAIR1(struct runner *r, struct cell *ci, struct cell *cj) { /* Get a pointer to the jth particle. */ struct part *restrict pi = &parts_i[sort_i[pid].i]; - - /* Compute the pairwise distance. */ float r2 = 0.0f; float dx[3]; @@ -923,9 +916,9 @@ void DOPAIR2(struct runner *r, struct cell *ci, struct cell *cj) { /* Anything to do here? */ if (!cell_is_active(ci, e) && !cell_is_active(cj, e)) return; - if(!cell_is_drifted(ci, e)) error("Cell ci not drifted"); - if(!cell_is_drifted(cj, e)) error("Cell cj not drifted"); - + if (!cell_is_drifted(ci, e)) error("Cell ci not drifted"); + if (!cell_is_drifted(cj, e)) error("Cell cj not drifted"); + /* Get the shift ID. */ double shift[3] = {0.0, 0.0, 0.0}; const int sid = space_getsid(e->s, &ci, &cj, shift); @@ -993,8 +986,6 @@ void DOPAIR2(struct runner *r, struct cell *ci, struct cell *cj) { const double di = sort_i[pid].d + hi * kernel_gamma + dx_max - rshift; if (di < dj_min) continue; - - double pix[3]; for (int k = 0; k < 3; k++) pix[k] = pi->x[k] - shift[k]; const float hig2 = hi * hi * kernel_gamma2; @@ -1009,8 +1000,6 @@ void DOPAIR2(struct runner *r, struct cell *ci, struct cell *cj) { struct part *restrict pj = &parts_j[sortdt_j[pjd].i]; const float hj = pj->h; - - /* Compute the pairwise distance. */ float r2 = 0.0f; float dx[3]; @@ -1062,8 +1051,6 @@ void DOPAIR2(struct runner *r, struct cell *ci, struct cell *cj) { struct part *restrict pj = &parts_j[sort_j[pjd].i]; const float hj = pj->h; - - /* Compute the pairwise distance. */ float r2 = 0.0f; float dx[3]; @@ -1143,8 +1130,6 @@ void DOPAIR2(struct runner *r, struct cell *ci, struct cell *cj) { const double dj = sort_j[pjd].d - hj * kernel_gamma - dx_max - rshift; if (dj > di_max) continue; - - double pjx[3]; for (int k = 0; k < 3; k++) pjx[k] = pj->x[k] + shift[k]; const float hjg2 = hj * hj * kernel_gamma2; @@ -1159,8 +1144,6 @@ void DOPAIR2(struct runner *r, struct cell *ci, struct cell *cj) { struct part *restrict pi = &parts_i[sortdt_i[pid].i]; const float hi = pi->h; - - /* Compute the pairwise distance. */ float r2 = 0.0f; float dx[3]; @@ -1211,8 +1194,6 @@ void DOPAIR2(struct runner *r, struct cell *ci, struct cell *cj) { struct part *restrict pi = &parts_i[sort_i[pid].i]; const float hi = pi->h; - - /* Compute the pairwise distance. */ float r2 = 0.0f; float dx[3]; @@ -1328,7 +1309,7 @@ void DOSELF1(struct runner *r, struct cell *restrict c) { if (!cell_is_active(c, e)) return; - if(!cell_is_drifted(c, e)) cell_drift(c, e); + if (!cell_is_drifted(c, e)) cell_drift(c, e); struct part *restrict parts = c->parts; const int count = c->count; @@ -1357,8 +1338,6 @@ void DOSELF1(struct runner *r, struct cell *restrict c) { const float hi = pi->h; const float hig2 = hi * hi * kernel_gamma2; - - /* Is the ith particle inactive? */ if (!part_is_active(pi, e)) { @@ -1369,8 +1348,6 @@ void DOSELF1(struct runner *r, struct cell *restrict c) { struct part *restrict pj = &parts[indt[pjd]]; const float hj = pj->h; - - /* Compute the pairwise distance. */ float r2 = 0.0f; float dx[3]; @@ -1425,8 +1402,6 @@ void DOSELF1(struct runner *r, struct cell *restrict c) { struct part *restrict pj = &parts[pjd]; const float hj = pj->h; - - /* Compute the pairwise distance. */ float r2 = 0.0f; float dx[3]; @@ -1567,7 +1542,7 @@ void DOSELF2(struct runner *r, struct cell *restrict c) { if (!cell_is_active(c, e)) return; - if(!cell_is_drifted(c, e)) error("Cell is not drifted"); + if (!cell_is_drifted(c, e)) error("Cell is not drifted"); struct part *restrict parts = c->parts; const int count = c->count; @@ -1596,8 +1571,6 @@ void DOSELF2(struct runner *r, struct cell *restrict c) { const float hi = pi->h; const float hig2 = hi * hi * kernel_gamma2; - - /* Is the ith particle not active? */ if (!part_is_active(pi, e)) { @@ -1608,8 +1581,6 @@ void DOSELF2(struct runner *r, struct cell *restrict c) { struct part *restrict pj = &parts[indt[pjd]]; const float hj = pj->h; - - /* Compute the pairwise distance. */ float r2 = 0.0f; float dx[3]; @@ -1664,8 +1635,6 @@ void DOSELF2(struct runner *r, struct cell *restrict c) { struct part *restrict pj = &parts[pjd]; const float hj = pj->h; - - /* Compute the pairwise distance. */ float r2 = 0.0f; float dx[3]; diff --git a/src/scheduler.c b/src/scheduler.c index b530b24ae864b65590bbd0c8ffeffa738ca5fa3b..82ea453b6369259b077eddd31db1d00acf7da37f 100644 --- a/src/scheduler.c +++ b/src/scheduler.c @@ -216,7 +216,6 @@ static void scheduler_splittask(struct task *t, struct scheduler *s) { to make sure we get ci and cj swapped if needed. */ double shift[3]; const int sid = space_getsid(s->space, &ci, &cj, shift); - const int did = space_getdid(s->space, ci, cj); /* Should this task be split-up? */ if (ci->split && cj->split && @@ -617,28 +616,10 @@ static void scheduler_splittask(struct task *t, struct scheduler *s) { 1 << sid, 0, ci, NULL, 0); else ci->sorts->flags |= (1 << sid); + lock_unlock_blind(&ci->lock); scheduler_addunlock(s, ci->sorts, t); - /* Create the drift for ci. */ - if (ci->drift == NULL) { - // ci->drift = scheduler_addtask(s, task_type_drift, - // task_subtype_none, - // 1 << sid, 0, ci, NULL, 0); - // scheduler_addunlock(s, ci->drift, ci->sorts); - - // scheduler_addunlock(s, ci->drift, t); - } - - /* if(!(ci->drift->flags & (1 << sid))) { */ - /* scheduler_addunlock(s, ci->drift, ci->sorts); */ - /* ci->drift->flags |= (1 << sid); */ - /* } */ - - if (did == 0 || did > 31) message("did=%d 1<<did=%d", did, 1 << did); - - lock_unlock_blind(&ci->lock); - /* Create the sort for cj. */ lock_lock(&cj->lock); if (cj->sorts == NULL) @@ -646,26 +627,9 @@ static void scheduler_splittask(struct task *t, struct scheduler *s) { 1 << sid, 0, cj, NULL, 0); else cj->sorts->flags |= (1 << sid); + lock_unlock_blind(&cj->lock); scheduler_addunlock(s, cj->sorts, t); - - /* Create the drift for cj. */ - if (cj->drift == NULL) { - // cj->drift = scheduler_addtask(s, task_type_drift, - // task_subtype_none, - // 1 << sid, 0, cj, NULL, 0); - // scheduler_addunlock(s, cj->drift, cj->sorts); - - // scheduler_addunlock(s, cj->drift, t); - // scheduler_addunlock(s, cj->drift, cj->init); - } - - /* if(!(cj->drift->flags & (1 << sid))) { */ - /* scheduler_addunlock(s, cj->drift, cj->sorts); */ - /* cj->drift->flags |= (1 << sid); */ - /* } */ - - lock_unlock_blind(&cj->lock); } } /* pair interaction? */ diff --git a/src/space.c b/src/space.c index 2c67cfce2e1f306f2a27c4b91f37a36e3cafe9a2..145b0bf8b84780ca543d0b27667384a362c04bd9 100644 --- a/src/space.c +++ b/src/space.c @@ -168,32 +168,6 @@ int space_getsid(struct space *s, struct cell **ci, struct cell **cj, return sid; } -int space_getdid(struct space *s, struct cell *ci, struct cell *cj) { - - /* Get the relative distance between the pairs, wrapping. */ - const int periodic = s->periodic; - double dx[3]; - double shift[3]; - for (int k = 0; k < 3; k++) { - dx[k] = cj->loc[k] - ci->loc[k]; - if (periodic && dx[k] < -s->dim[k] / 2) - shift[k] = s->dim[k]; - else if (periodic && dx[k] > s->dim[k] / 2) - shift[k] = -s->dim[k]; - else - shift[k] = 0.0; - dx[k] += shift[k]; - } - - /* Get the drift index. */ - int did = 0; - for (int k = 0; k < 3; k++) - did = 3 * did + ((dx[k] < 0.0) ? 0 : ((dx[k] > 0.0) ? 2 : 1)); - - /* Return the drift ID. */ - return did; -} - /** * @brief Recursively dismantle a cell tree. * diff --git a/src/space.h b/src/space.h index 2f53e24ad4249548e3468358bece3ae9deee4ed2..53cf2d0c8fa548ae19aa7452abb38c3e3e028165 100644 --- a/src/space.h +++ b/src/space.h @@ -150,7 +150,6 @@ void space_gparts_sort(struct space *s, int *ind, size_t N, int min, int max, struct cell *space_getcell(struct space *s); int space_getsid(struct space *s, struct cell **ci, struct cell **cj, double *shift); -int space_getdid(struct space *s, struct cell *ci, struct cell *cj); void space_init(struct space *s, const struct swift_params *params, double dim[3], struct part *parts, struct gpart *gparts, size_t Npart, size_t Ngpart, int periodic, int gravity,