diff --git a/src/engine.c b/src/engine.c index f152177338e6e2a43c5b9fe8351ff372d0aadbe6..b577ae3bf92234749802571f7b1d1c2462cb7d07 100644 --- a/src/engine.c +++ b/src/engine.c @@ -2110,11 +2110,11 @@ void engine_init_particles(struct engine *e) { if (e->policy & engine_policy_self_gravity) { mask |= 1 << task_type_grav_up; - mask |= 1 << task_type_self; - mask |= 1 << task_type_pair; - mask |= 1 << task_type_sub; + /* mask |= 1 << task_type_self; */ + /* mask |= 1 << task_type_pair; */ + /* mask |= 1 << task_type_sub; */ - submask |= 1 << task_subtype_grav; + /* submask |= 1 << task_subtype_grav; */ } /* Add the tasks corresponding to external gravity to the masks */ @@ -2335,20 +2335,20 @@ void engine_step(struct engine *e) { fclose(file); /* Check the gravity accelerations */ - /* struct gpart *temp = malloc(s->nr_gparts * sizeof(struct gpart)); */ - /* memcpy(temp, s->gparts, s->nr_gparts * sizeof(struct gpart)); */ - /* gravity_n2(temp, s->nr_gparts); */ - /* file = fopen("grav_brute.dat", "w"); */ - /* for(size_t k = 0; k < s->nr_gparts; ++k) { */ - /* fprintf(file, "%lld %f %f %f %f %f %f\n", temp[k].id, */ - /* temp[k].x[0], temp[k].x[1], temp[k].x[2], */ - /* temp[k].a_grav[0], temp[k].a_grav[1], temp[k].a_grav[2]); */ - /* } */ - /* fclose(file); */ + struct gpart *temp = malloc(s->nr_gparts * sizeof(struct gpart)); + memcpy(temp, s->gparts, s->nr_gparts * sizeof(struct gpart)); + gravity_n2(temp, s->nr_gparts); + file = fopen("grav_brute.dat", "w"); + for(size_t k = 0; k < s->nr_gparts; ++k) { + fprintf(file, "%lld %f %f %f %f %f %f\n", temp[k].id, + temp[k].x[0], temp[k].x[1], temp[k].x[2], + temp[k].a_grav[0], temp[k].a_grav[1], temp[k].a_grav[2]); + } + fclose(file); - /* free(temp); */ + free(temp); - float mass = 0.f; + double mass = 0.f; for (int i = 0; i < s->nr_cells; ++i) mass += s->cells[i].multipole.mass; message("Total mass: %f", mass); diff --git a/src/gravity/Default/gravity_part.h b/src/gravity/Default/gravity_part.h index 71d8f34aa8d38ebd5d0f2fd681e42c0662dd166e..2d3bc69ba6f9414a8e9ccaf1b5662227184f7aa1 100644 --- a/src/gravity/Default/gravity_part.h +++ b/src/gravity/Default/gravity_part.h @@ -47,7 +47,7 @@ struct gpart { /* float tx; */ /* float tv; */ - float mass_interacted; + double mass_interacted; /* Anonymous union for id/part. */ union { diff --git a/src/runner_doiact_grav.h b/src/runner_doiact_grav.h index 8ac9e5b65b5a00ebb9caf4bb67aa8943a77f4176..9ec64f74f1f41e8a68a36e1efa42f7575867a4c5 100644 --- a/src/runner_doiact_grav.h +++ b/src/runner_doiact_grav.h @@ -25,7 +25,7 @@ #include "clocks.h" #include "part.h" -#define ICHECK 4934 +#define ICHECK 1 /** * @brief Compute the recursive upward sweep, i.e. construct the @@ -125,8 +125,9 @@ __attribute__((always_inline)) INLINE static void runner_dopair_grav_pm( if (gp->ti_end > ti_current) continue; if (gp->id == -ICHECK) - message("id=%lld mass=%f CoM=[%f %f %f]\n", gp->id, multi.mass, - multi.CoM[0], multi.CoM[1], multi.CoM[2]); + message("id=%lld x=[%f %f %f] mass=%f CoM=[ %f %f %f ] size= %f\n", gp->id, + gp->x[0], gp->x[1], gp->x[2], multi.mass, + multi.CoM[0], multi.CoM[1], multi.CoM[2], cj->h[0]); /* Compute the pairwise distance. */ const float dx[3] = {multi.CoM[0] - gp->x[0], // x @@ -210,6 +211,23 @@ __attribute__((always_inline)) INLINE static void runner_dopair_grav_pp( /* Anything to do here? */ // if (ci->ti_end_min > ti_current && cj->ti_end_min > ti_current) return; + for (int pid = 0; pid < gcount_i; pid++) { + + struct gpart *restrict gp = &gparts_i[pid]; + if (gp->id == -ICHECK) + message("id=%lld size=%f\n", gp->id, cj->h[0]); + + } + + for (int pid = 0; pid < gcount_j; pid++) { + + struct gpart *restrict gp = &gparts_j[pid]; + if (gp->id == -ICHECK) + message("id=%lld size=%f\n", gp->id, ci->h[0]); + + } + + /* Loop over all particles in ci... */ for (int pid = 0; pid < gcount_i; pid++) { @@ -371,23 +389,23 @@ static void runner_dopair_grav(struct runner *r, struct cell *ci, /* cj->loc[2], cj->h[0]); */ } - for (int i = 0; i < gcount_i; i++) { + /* for (int i = 0; i < gcount_i; i++) { */ - struct gpart *const gp = &ci->gparts[i]; + /* struct gpart *const gp = &ci->gparts[i]; */ - if (gp->id == -ICHECK) - message("id=%lld a=[%f %f %f]\n", gp->id, gp->a_grav[0], gp->a_grav[1], - gp->a_grav[2]); - } + /* if (gp->id == -ICHECK) */ + /* message("id=%lld a=[%f %f %f]\n", gp->id, gp->a_grav[0], gp->a_grav[1], */ + /* gp->a_grav[2]); */ + /* } */ - for (int i = 0; i < gcount_j; i++) { + /* for (int i = 0; i < gcount_j; i++) { */ - struct gpart *const gp = &cj->gparts[i]; + /* struct gpart *const gp = &cj->gparts[i]; */ - if (gp->id == -ICHECK) - message("id=%lld a=[%f %f %f]\n", gp->id, gp->a_grav[0], gp->a_grav[1], - gp->a_grav[2]); - } + /* if (gp->id == -ICHECK) */ + /* message("id=%lld a=[%f %f %f]\n", gp->id, gp->a_grav[0], gp->a_grav[1], */ + /* gp->a_grav[2]); */ + /* } */ /* Are both cells split ? */ if (ci->split && cj->split) { diff --git a/src/scheduler.c b/src/scheduler.c index b9aac5a160738510a06dc5a2e1538c8d183c78f8..6271f80ca306f3546f12360f18eb7db038d97fa3 100644 --- a/src/scheduler.c +++ b/src/scheduler.c @@ -172,8 +172,8 @@ void scheduler_splittasks(struct scheduler *s) { if (ci->split) { /* Make a sub? */ - if (scheduler_dosub && (ci->count * ci->count < space_subsize || - ci->gcount * ci->gcount < space_subsize)) { + if (scheduler_dosub && (ci->count * ci->count < space_subsize)) { + // ci->gcount * ci->gcount < space_subsize)) { /* convert to a self-subtask. */ t->type = task_type_sub; @@ -238,8 +238,8 @@ void scheduler_splittasks(struct scheduler *s) { /* Replace by a single sub-task? */ if (scheduler_dosub && - (ci->count * cj->count * sid_scale[sid] < space_subsize || - ci->gcount * cj->gcount * sid_scale[sid] < space_subsize) && + (ci->count * cj->count * sid_scale[sid] < space_subsize) &&// || + //ci->gcount * cj->gcount * sid_scale[sid] < space_subsize) && sid != 0 && sid != 2 && sid != 6 && sid != 8) { /* Make this task a sub task. */