Commit d13e9833 authored by Tom Theuns's avatar Tom Theuns
Browse files

added external_gravity to task mask. Need to replace part by gpart in dograv_external

parent c6cb767c
......@@ -62,6 +62,8 @@
#define ENGINE_HYDRO 0
#endif
/**
* @brief Main routine that loads a few particles and generates some output.
*
......@@ -433,7 +435,7 @@ int main(int argc, char *argv[]) {
tic = getticks();
if (myrank == 0) message("nr_nodes is %i.", nr_nodes);
engine_init(&e, &s, dt_max, nr_threads, nr_queues, nr_nodes, myrank,
ENGINE_POLICY | engine_policy_steal | ENGINE_HYDRO, 0,
ENGINE_POLICY | engine_policy_steal | ENGINE_HYDRO | engine_policy_external_gravity , 0,
time_end, dt_min, dt_max);
if (myrank == 0)
message("engine_init took %.3f ms.",
......
......@@ -122,6 +122,9 @@ struct cell {
/* Tasks for gravity tree. */
struct task *grav_up, *grav_down;
/* Task for external gravity */
struct task *grav_external;
/* Number of tasks that are associated with this cell. */
int nr_tasks;
......
......@@ -127,12 +127,12 @@ void engine_mkghosts(struct engine *e, struct cell *c, struct cell *super) {
if(c->gcount > 0)
{
/* /\* Add the gravity tasks *\/ */
/* c->grav_external = scheduler_addtask(s, task_type_grav_external, task_subtype_none, 0, 0, */
/* c, NULL, 0); */
/* Add the gravity tasks */
c->grav_external = scheduler_addtask(s, task_type_grav_external, task_subtype_none, 0, 0,
c, NULL, 0);
/* /\* Enforce gravity calculated before kick *\/ */
/* scheduler_addunlock(s, c->grav_external, c->kick); */
/* Enforce gravity calculated before kick */
scheduler_addunlock(s, c->grav_external, c->kick);
}
}
}
......@@ -791,20 +791,8 @@ void engine_maketasks(struct engine *e) {
}
#ifdef GRAVITY
#ifdef DEFAULT_GRAVITY
/* Add the gravity mm tasks. */
for (i = 0; i < nr_cells; i++)
if (cells[i].gcount > 0) {
scheduler_addtask(sched, task_type_grav_mm, task_subtype_none, -1, 0,
&cells[i], NULL, 0);
for (j = i + 1; j < nr_cells; j++)
if (cells[j].gcount > 0)
scheduler_addtask(sched, task_type_grav_mm, task_subtype_none, -1, 0,
&cells[i], &cells[j], 0);
}
#endif
/* #ifdef EXTERNAL_POTENTIAL */
/* /\* Add the external potential gravity *\/ */
/* #ifdef DEFAULT_GRAVITY */
/* /\* Add the gravity mm tasks. *\/ */
/* for (i = 0; i < nr_cells; i++) */
/* if (cells[i].gcount > 0) { */
/* scheduler_addtask(sched, task_type_grav_mm, task_subtype_none, -1, 0, */
......@@ -814,10 +802,7 @@ void engine_maketasks(struct engine *e) {
/* scheduler_addtask(sched, task_type_grav_mm, task_subtype_none, -1, 0, */
/* &cells[i], &cells[j], 0); */
/* } */
/* #endif */
#endif
scheduler_print_tasks(sched, "sched.txt");
/* Split the tasks. */
scheduler_splittasks(sched);
......@@ -1017,6 +1002,10 @@ void engine_maketasks(struct engine *e) {
/* Set the tasks age. */
e->tasks_age = 0;
/* print all possibly scheduled tasks */
scheduler_print_tasks(sched, "sched.txt");
}
/**
......@@ -1463,15 +1452,15 @@ void engine_init_particles(struct engine *e) {
/* Add the tasks corresponding to self-gravity to the masks */
if ((e->policy & engine_policy_self_gravity) == engine_policy_self_gravity) {
/* Nothing here for now */
}
/* Add the tasks corresponding to self-gravity to the masks */
/* Add the tasks corresponding to external gravity to the masks */
if ((e->policy & engine_policy_external_gravity) ==
engine_policy_external_gravity) {
mask |= 1 << task_type_grav_external;
/* Nothing here for now */
}
/* Add MPI tasks if need be */
......@@ -1634,8 +1623,7 @@ void engine_step(struct engine *e) {
/* Add the tasks corresponding to self-gravity to the masks */
if ((e->policy & engine_policy_external_gravity) ==
engine_policy_external_gravity) {
/* Nothing here for now */
mask |= 1 << task_type_grav_external;
}
/* Add MPI tasks if need be */
......
......@@ -121,33 +121,33 @@ void runner_dograv_external(struct runner *r, struct cell *c) {
#ifdef TASK_VERBOSE
OUT
#endif
/* /\* Loop over the parts in this cell. *\/ */
/* for (i = 0; i < count; i++) { */
/* /\* Get a direct pointer on the part. *\/ */
/* p = &parts[i]; */
/* /\* Is this part within the time step? *\/ */
/* if (p->t_end <= t_end) { */
/* rinv = 1 / sqrtf((p->x[0]-External_Potential_X)*(p->x[0]-External_Potential_X) + (p->x[1]-External_Potential_Y)*(p->x[1]-External_Potential_Y) + (p->x[2]-External_Potential_Z)*(p->x[2]-External_Potential_Z)); */
/* /\* check for energy and angular momentum conservation *\/ */
/* E = 0.5 * ((p->v[0]*p->v[0]) + (p->v[1]*p->v[1]) + (p->v[2]*p->v[2])) - const_G * External_Potential_Mass * rinv; */
/* L[0] = (p->x[1] - External_Potential_X)*p->v[2] - (p->x[2] - External_Potential_X)*p->v[1]; */
/* L[1] = (p->x[2] - External_Potential_Y)*p->v[0] - (p->x[0] - External_Potential_Y)*p->v[2]; */
/* L[2] = (p->x[0] - External_Potential_Z)*p->v[1] - (p->x[1] - External_Potential_Z)*p->v[0]; */
/* if(p->id == 0) { */
/* message("update %f\t %f\t %f\t %f\t %f\t %f\t %f\t %f\t %f\n", r->e->time, 1./rinv, p->x[0], p->x[1], p->x[2], E, L[0], L[1], L[2]); */
/* message(" ... %f %f %f\n", p->v[0], p->v[1], p->v[2]); */
/* message(" ... %e %e\n", const_G, External_Potential_Mass); */
/* } */
/* p->grav_accel[0] = - const_G * External_Potential_Mass * (p->x[0] - External_Potential_X) * rinv * rinv * rinv; */
/* p->grav_accel[1] = - const_G * External_Potential_Mass * (p->x[1] - External_Potential_Y) * rinv * rinv * rinv; */
/* p->grav_accel[2] = - const_G * External_Potential_Mass * (p->x[2] - External_Potential_Z) * rinv * rinv * rinv; */
/* Loop over the parts in this cell. */
for (i = 0; i < count; i++) {
/* Get a direct pointer on the part. */
p = &parts[i];
/* Is this part within the time step? */
if (p->t_end <= t_end) {
rinv = 1 / sqrtf((p->x[0]-External_Potential_X)*(p->x[0]-External_Potential_X) + (p->x[1]-External_Potential_Y)*(p->x[1]-External_Potential_Y) + (p->x[2]-External_Potential_Z)*(p->x[2]-External_Potential_Z));
/* check for energy and angular momentum conservation */
E = 0.5 * ((p->v[0]*p->v[0]) + (p->v[1]*p->v[1]) + (p->v[2]*p->v[2])) - const_G * External_Potential_Mass * rinv;
L[0] = (p->x[1] - External_Potential_X)*p->v[2] - (p->x[2] - External_Potential_X)*p->v[1];
L[1] = (p->x[2] - External_Potential_Y)*p->v[0] - (p->x[0] - External_Potential_Y)*p->v[2];
L[2] = (p->x[0] - External_Potential_Z)*p->v[1] - (p->x[1] - External_Potential_Z)*p->v[0];
if(p->id == 0) {
message("update %f\t %f\t %f\t %f\t %f\t %f\t %f\t %f\t %f\n", r->e->time, 1./rinv, p->x[0], p->x[1], p->x[2], E, L[0], L[1], L[2]);
message(" ... %f %f %f\n", p->v[0], p->v[1], p->v[2]);
message(" ... %e %e\n", const_G, External_Potential_Mass);
}
p->grav_accel[0] = - const_G * External_Potential_Mass * (p->x[0] - External_Potential_X) * rinv * rinv * rinv;
p->grav_accel[1] = - const_G * External_Potential_Mass * (p->x[1] - External_Potential_Y) * rinv * rinv * rinv;
p->grav_accel[2] = - const_G * External_Potential_Mass * (p->x[2] - External_Potential_Z) * rinv * rinv * rinv;
/* } */
/* } */
/* TIMER_TOC(timer_dograv_external); */
}
}
TIMER_TOC(timer_dograv_external);
}
......@@ -1077,6 +1077,10 @@ void runner_dokick(struct runner *r, struct cell *c, int timer) {
g->v[1] += g->a_grav_external[1] * dt;
g->v[2] += g->a_grav_external[2] * dt;
/* Minimal time for next end of time-step */
ti_end_min = min(g->ti_end, ti_end_min);
ti_end_max = max(g->ti_end, ti_end_max);
/* Number of updated particles */
updated++;
......@@ -1250,6 +1254,9 @@ void *runner_main(void *data) {
case task_type_grav_down:
runner_dograv_down(r, t->ci);
break;
case task_type_grav_external:
runner_dograv_external(r, t->ci);
break;
case task_type_psort:
space_do_parts_sort();
break;
......
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