diff --git a/examples/ExternalPointMass/externalPointMass.yml b/examples/ExternalPointMass/externalPointMass.yml
index 9bc67721b7ef99fb8097d21adfc6532c6c32a865..b001c42624505c450d27914873a712fbeac6b4c3 100644
--- a/examples/ExternalPointMass/externalPointMass.yml
+++ b/examples/ExternalPointMass/externalPointMass.yml
@@ -28,7 +28,7 @@ SPH:
delta_neighbours: 1. # The tolerance for the targetted number of neighbours.
CFL_condition: 0.1 # Courant-Friedrich-Levy condition for time integration.
max_ghost_iterations: 30 # Maximal number of iterations allowed to converge towards the smoothing length.
- max_smoothing_length: 3. # Maximal smoothing length allowed (in internal units).
+ max_smoothing_length: 10. # Maximal smoothing length allowed (in internal units).
# Parameters related to the initial conditions
InitialConditions:
diff --git a/src/gravity/Default/gravity.h b/src/gravity/Default/gravity.h
index 2040c73602a93f010e623bab2da9b93dc8ac0b2a..4dd8d5b7c7194ae072f90e7703907500c208b9b9 100644
--- a/src/gravity/Default/gravity.h
+++ b/src/gravity/Default/gravity.h
@@ -22,7 +22,8 @@
#include "potentials.h"
/**
- * @brief Computes the gravity time-step of a given particle.
+ * @brief Computes the gravity time-step of a given particle due to an external
+ *potential.
*
* This function only branches towards the potential chosen by the user.
*
@@ -30,20 +31,40 @@
* @param phys_const The physical constants in internal units.
* @param gp Pointer to the g-particle data.
*/
-__attribute__((always_inline)) INLINE static float gravity_compute_timestep(
- const struct external_potential* potential,
- const struct phys_const* const phys_const, const struct gpart* const gp) {
+__attribute__((always_inline))
+ INLINE static float gravity_compute_timestep_external(
+ const struct external_potential* potential,
+ const struct phys_const* const phys_const,
+ const struct gpart* const gp) {
float dt = FLT_MAX;
#ifdef EXTERNAL_POTENTIAL_POINTMASS
- dt =
- fminf(dt, external_gravity_pointmass_timestep(potential, phys_const, gp));
+ dt =
+ fminf(dt, external_gravity_pointmass_timestep(potential, phys_const, gp));
#endif
return dt;
}
+/**
+ * @brief Computes the gravity time-step of a given particle due to self-gravity
+ *
+ * This function only branches towards the potential chosen by the user.
+ *
+ * @param phys_const The physical constants in internal units.
+ * @param gp Pointer to the g-particle data.
+ */
+__attribute__((always_inline))
+ INLINE static float gravity_compute_timestep_self(
+ const struct phys_const* const phys_const,
+ const struct gpart* const gp) {
+
+ float dt = FLT_MAX;
+
+ return dt;
+}
+
/**
* @brief Initialises the g-particles for the first time
*
diff --git a/src/runner.c b/src/runner.c
index bb24002e14442c30be6e10302a4dca98e02643b4..fddbe5dd457b9d40e3222e43b07bc5d12e812570 100644
--- a/src/runner.c
+++ b/src/runner.c
@@ -111,6 +111,8 @@ void runner_dograv_external(struct runner *r, struct cell *c, int timer) {
struct gpart *g, *gparts = c->gparts;
int i, k, gcount = c->gcount;
const int ti_current = r->e->ti_current;
+ const struct external_potential *potential = r->e->potential;
+ const struct phys_const *constants = r->e->physical_constants;
TIMER_TIC;
@@ -134,7 +136,7 @@ void runner_dograv_external(struct runner *r, struct cell *c, int timer) {
/* Is this part within the time step? */
if (g->ti_end <= ti_current) {
- external_gravity(r->e->potential, r->e->physical_constants, g);
+ external_gravity(potential, constants, g);
/* /\* check for energy and angular momentum conservation - begin by */
/* * synchronizing velocity*\/ */
@@ -900,6 +902,8 @@ void runner_dokick(struct runner *r, struct cell *c, int timer) {
struct part *const parts = c->parts;
struct xpart *const xparts = c->xparts;
struct gpart *const gparts = c->gparts;
+ const struct external_potential *potential = r->e->potential;
+ const struct phys_const *constants = r->e->physical_constants;
const int is_fixdt =
(r->e->policy & engine_policy_fixdt) == engine_policy_fixdt;
@@ -941,8 +945,12 @@ void runner_dokick(struct runner *r, struct cell *c, int timer) {
} else {
/* Compute the next timestep (gravity condition) */
- float new_dt = gravity_compute_timestep(r->e->potential,
- r->e->physical_constants, gp);
+ const float new_dt_external =
+ gravity_compute_timestep_external(potential, constants, gp);
+ const float new_dt_self =
+ gravity_compute_timestep_self(constants, gp);
+
+ float new_dt = fminf(new_dt_external, new_dt_self);
/* Limit timestep within the allowed range */
new_dt = fminf(new_dt, global_dt_max);
@@ -1026,10 +1034,17 @@ void runner_dokick(struct runner *r, struct cell *c, int timer) {
/* Compute the next timestep (gravity condition) */
float new_dt_grav = FLT_MAX;
- if (p->gpart != NULL)
- new_dt_grav = gravity_compute_timestep(
- r->e->potential, r->e->physical_constants, p->gpart);
+ if (p->gpart != NULL) {
+
+ const float new_dt_external = gravity_compute_timestep_external(
+ potential, constants, p->gpart);
+ const float new_dt_self =
+ gravity_compute_timestep_self(constants, p->gpart);
+
+ new_dt_grav = fminf(new_dt_external, new_dt_self);
+ }
+ /* Final time-step is minimum of hydro and gravity */
float new_dt = fminf(new_dt_hydro, new_dt_grav);
/* Limit change in h */
diff --git a/src/space.c b/src/space.c
index 2e4fe71bfd248948ffd54eb17782269bfba178fb..9727ad664ad580b55fdc18d126f862eaa49428f0 100644
--- a/src/space.c
+++ b/src/space.c
@@ -177,10 +177,6 @@ void space_regrid(struct space *s, double cell_max, int verbose) {
}
s->h_max = h_max;
}
- } else {
- /* It would be nice to replace this with something more physical or
- * meaningful */
- h_max = s->dim[0] / 16.f;
}
/* If we are running in parallel, make sure everybody agrees on