diff --git a/configure.ac b/configure.ac index 6c0a06005778c499c53ae8e596a0685a78886542..e979580470671143348d986e82ad1d1861eb3129 100644 --- a/configure.ac +++ b/configure.ac @@ -215,6 +215,21 @@ elif test "$gravity_force_checks" != "no"; then AC_DEFINE_UNQUOTED([SWIFT_GRAVITY_FORCE_CHECKS], [$enableval] ,[Enable gravity brute-force checks]) fi +# Check if we want to zero the gravity forces for all particles below some ID. +AC_ARG_ENABLE([no-gravity-below-id], + [AS_HELP_STRING([--enable-no-gravity-below-id], + [Zeros the gravitational acceleration of all particles with an ID smaller than @<:@N@:>@] + )], + [no_gravity_below_id="$enableval"], + [no_gravity_below_id="no"] +) +if test "$no_gravity_below_id" == "yes"; then + AC_MSG_ERROR(Need to specify the ID below which particles get zero forces when using --enable-no-gravity-below-id!) +elif test "$no_gravity_below_id" != "no"; then + AC_DEFINE_UNQUOTED([SWIFT_NO_GRAVITY_BELOW_ID], [$enableval] ,[Particles with smaller ID than this will have zero gravity forces]) +fi + + # Define HAVE_POSIX_MEMALIGN if it works. AX_FUNC_POSIX_MEMALIGN @@ -854,12 +869,14 @@ AC_MSG_RESULT([ Adiabatic index : $with_gamma Riemann solver : $with_riemann Cooling function : $with_cooling - External potential : $with_potential - Multipole order : $with_multipole_order - - Task debugging : $enable_task_debugging - Debugging checks : $enable_debugging_checks - Gravity checks : $gravity_force_checks + + External potential : $with_potential + Multipole order : $with_multipole_order + No gravity below ID : $no_gravity_below_id + + Task debugging : $enable_task_debugging + Debugging checks : $enable_debugging_checks + Gravity checks : $gravity_force_checks ]) # Make sure the latest git revision string gets included diff --git a/src/runner.c b/src/runner.c index 29d715a329940313041f1d8e15c8f067c5c24408..18e48c945db389e935924b3356d68d63349f45cf 100644 --- a/src/runner.c +++ b/src/runner.c @@ -1356,9 +1356,8 @@ void runner_do_end_force(struct runner *r, struct cell *c, int timer) { if (part_is_active(p, e)) { - /* First, finish the force loop */ + /* Finish the force loop */ hydro_end_force(p); - if (p->gpart != NULL) gravity_end_force(p->gpart, const_G); } } @@ -1368,25 +1367,45 @@ void runner_do_end_force(struct runner *r, struct cell *c, int timer) { /* Get a handle on the gpart. */ struct gpart *restrict gp = &gparts[k]; - if (gp->type == swift_type_dark_matter) { + if (gpart_is_active(gp, e)) { - if (gpart_is_active(gp, e)) gravity_end_force(gp, const_G); - } + /* Finish the force calculation */ + gravity_end_force(gp, const_G); + +#ifdef SWIFT_NO_GRAVITY_BELOW_ID + /* Cancel gravity forces of these particles */ + if ((gp->type == swift_type_dark_matter && + gp->id_or_neg_offset < SWIFT_NO_GRAVITY_BELOW_ID) || + (gp->type == swift_type_gas && + parts[-gp->id_or_neg_offset].id < SWIFT_NO_GRAVITY_BELOW_ID) || + (gp->type == swift_type_star && + sparts[-gp->id_or_neg_offset].id < SWIFT_NO_GRAVITY_BELOW_ID)) { + + /* Don't move ! */ + gp->a_grav[0] = 0.f; + gp->a_grav[1] = 0.f; + gp->a_grav[2] = 0.f; + } +#endif #ifdef SWIFT_DEBUG_CHECKS - if (e->policy & engine_policy_self_gravity && gpart_is_active(gp, e)) { + if (e->policy & engine_policy_self_gravity) { - /* Check that this gpart has interacted with all the other particles - * (via direct or multipoles) in the box */ - gp->num_interacted++; - if (gp->num_interacted != (long long)e->s->nr_gparts) - error( - "g-particle (id=%lld, type=%d) did not interact gravitationally " - "with all other gparts gp->num_interacted=%lld, total_gparts=%zd", - gp->id_or_neg_offset, gp->type, gp->num_interacted, - e->s->nr_gparts); - } + /* Check that this gpart has interacted with all the other + * particles + * (via direct or multipoles) in the box */ + gp->num_interacted++; + if (gp->num_interacted != (long long)e->s->nr_gparts) + error( + "g-particle (id=%lld, type=%d) did not interact " + "gravitationally " + "with all other gparts gp->num_interacted=%lld, " + "total_gparts=%zd", + gp->id_or_neg_offset, gp->type, gp->num_interacted, + e->s->nr_gparts); + } #endif + } } /* Loop over the star particles in this cell. */ @@ -1396,9 +1415,8 @@ void runner_do_end_force(struct runner *r, struct cell *c, int timer) { struct spart *restrict sp = &sparts[k]; if (spart_is_active(sp, e)) { - /* First, finish the force loop */ + /* Finish the force loop */ star_end_force(sp); - gravity_end_force(sp->gpart, const_G); } } }