diff --git a/examples/EAGLE_6/eagle_6.yml b/examples/EAGLE_6/eagle_6.yml index de1028ca6691223123e677c873ad3ee2f2ddfe92..0be082dd758a974ff9c9b1c0284e3a2a3efa6d00 100644 --- a/examples/EAGLE_6/eagle_6.yml +++ b/examples/EAGLE_6/eagle_6.yml @@ -24,7 +24,8 @@ TimeIntegration: time_end: 1e-2 # The end time of the simulation (in internal units). dt_min: 1e-10 # The minimal time-step size of the simulation (in internal units). dt_max: 1e-4 # The maximal time-step size of the simulation (in internal units). - + max_dt_RMS_factor: 0.25 + # Parameters governing the snapshots Snapshots: basename: eagle # Common part of the name of output files diff --git a/src/atomic.h b/src/atomic.h index b09ed3dd22001586cfde8545e636de67a819c003..4e3407cc3344445ba4d7ec7a7c0715ab0c591859 100644 --- a/src/atomic.h +++ b/src/atomic.h @@ -24,6 +24,7 @@ /* Includes. */ #include "inline.h" +#include "minmax.h" #define atomic_add(v, i) __sync_fetch_and_add(v, i) #define atomic_sub(v, i) __sync_fetch_and_sub(v, i) @@ -33,4 +34,40 @@ #define atomic_cas(v, o, n) __sync_val_compare_and_swap(v, o, n) #define atomic_swap(v, n) __sync_lock_test_and_set(v, n) +/** + * @param Atomic min operation on floats. + */ +__attribute__((always_inline)) INLINE void atomic_min_f(float const* x, + float y) { + int done = 0; + while (!done) { + const float val = *x; + done = __sync_bool_compare_and_swap((int*)x, val, min(val, y)); + } +} + +/** + * @param Atomic max operation on floats. + */ +__attribute__((always_inline)) INLINE void atomic_max_f(float const* x, + float y) { + int done = 0; + while (!done) { + const float val = *x; + done = __sync_bool_compare_and_swap((int*)x, val, max(val, y)); + } +} + +/** + * @param Atomic add operation on floats. + */ +__attribute__((always_inline)) INLINE void atomic_add_f(float const* x, + float y) { + int done = 0; + while (!done) { + const float val = *x; + done = __sync_bool_compare_and_swap((int*)x, val, val + y); + } +} + #endif /* SWIFT_ATOMIC_H */ diff --git a/src/engine.c b/src/engine.c index af6150852e273fbc4dff66574a311f3e6b4f7a7f..bb9a478b4933106c3753a31388bdf78117366147 100644 --- a/src/engine.c +++ b/src/engine.c @@ -3804,6 +3804,10 @@ void engine_rebuild(struct engine *e, int clean_smoothing_length_values) { /* Re-build the space. */ space_rebuild(e->s, e->verbose); + /* Re-compute the maximal RMS displacement constraint */ + if (e->policy & engine_policy_cosmology) + engine_recompute_displacement_constraint(e); + #ifdef SWIFT_DEBUG_CHECKS part_verify_links(e->s->parts, e->s->gparts, e->s->sparts, e->s->nr_parts, e->s->nr_gparts, e->s->nr_sparts, e->verbose); @@ -5422,8 +5426,8 @@ void engine_init(struct engine *e, struct space *s, e->dt_min = parser_get_param_double(params, "TimeIntegration:dt_min"); e->dt_max = parser_get_param_double(params, "TimeIntegration:dt_max"); e->dt_max_RMS_displacement = FLT_MAX; - e->max_RMS_displacement_factor = - parser_get_param_double(params, "TimeIntegration:max_RMS_factor"); + e->max_RMS_displacement_factor = parser_get_opt_param_double( + params, "TimeIntegration:max_dt_RMS_factor", 1.); e->a_first_statistics = parser_get_opt_param_double(params, "Statistics:scale_factor_first", 0.1); e->time_first_statistics = @@ -6124,6 +6128,17 @@ void engine_compute_next_statistics_time(struct engine *e) { } } +/** + * @brief Computes the maximal time-step allowed by the max RMS displacement + * condition. + * + * @param e The #engine. + */ +void engine_recompute_displacement_constraint(struct engine *e) { + + message("max_dt_RMS_displacement = %e", e->dt_max_RMS_displacement); +} + /** * @brief Frees up the memory allocated for this #engine */ diff --git a/src/engine.h b/src/engine.h index e137595390d1ba75f0d318485f950780a57927a2..a6707e0426d485a61464c440915c94d6e74a580b 100644 --- a/src/engine.h +++ b/src/engine.h @@ -341,6 +341,7 @@ struct engine { void engine_barrier(struct engine *e); void engine_compute_next_snapshot_time(struct engine *e); void engine_compute_next_statistics_time(struct engine *e); +void engine_recompute_displacement_constraint(struct engine *e); void engine_unskip(struct engine *e); void engine_drift_all(struct engine *e); void engine_drift_top_multipoles(struct engine *e); diff --git a/src/space.c b/src/space.c index 9dc6037820063847d47624bad1f3389f5e113975..f0c5d4c8b91e329e9b593266a4cec9e678798af5 100644 --- a/src/space.c +++ b/src/space.c @@ -1087,6 +1087,9 @@ void space_parts_get_cell_index_mapper(void *map_data, int nr_parts, if (cell_counts == NULL) error("Failed to allocate temporary cell count buffer."); + /* Init the local minimal part mass */ + float min_mass = FLT_MAX; + /* Loop over the parts. */ for (int k = 0; k < nr_parts; k++) { @@ -1119,6 +1122,9 @@ void space_parts_get_cell_index_mapper(void *map_data, int nr_parts, pos_z); #endif + /* Compute minimal mass */ + min_mass = min(min_mass, hydro_get_mass(p)); + /* Update the position */ p->x[0] = pos_x; p->x[1] = pos_y; @@ -1129,6 +1135,9 @@ void space_parts_get_cell_index_mapper(void *map_data, int nr_parts, for (int k = 0; k < s->nr_cells; k++) if (cell_counts[k]) atomic_add(&data->cell_counts[k], cell_counts[k]); free(cell_counts); + + /* Write back thee minimal part mass */ + atomic_min_f(&s->min_part_mass, min_mass); } /** @@ -1161,6 +1170,9 @@ void space_gparts_get_cell_index_mapper(void *map_data, int nr_gparts, if (cell_counts == NULL) error("Failed to allocate temporary cell count buffer."); + /* Init the local minimal part mass */ + float min_mass = FLT_MAX; + for (int k = 0; k < nr_gparts; k++) { /* Get the particle */ @@ -1192,6 +1204,9 @@ void space_gparts_get_cell_index_mapper(void *map_data, int nr_gparts, pos_z); #endif + /* Compute minimal mass */ + if (gp->type == swift_type_dark_matter) min_mass = min(min_mass, gp->mass); + /* Update the position */ gp->x[0] = pos_x; gp->x[1] = pos_y; @@ -1202,6 +1217,9 @@ void space_gparts_get_cell_index_mapper(void *map_data, int nr_gparts, for (int k = 0; k < s->nr_cells; k++) if (cell_counts[k]) atomic_add(&data->cell_counts[k], cell_counts[k]); free(cell_counts); + + /* Write back thee minimal part mass */ + atomic_min_f(&s->min_gpart_mass, min_mass); } /** @@ -1234,6 +1252,9 @@ void space_sparts_get_cell_index_mapper(void *map_data, int nr_sparts, if (cell_counts == NULL) error("Failed to allocate temporary cell count buffer."); + /* Init the local minimal part mass */ + float min_mass = FLT_MAX; + for (int k = 0; k < nr_sparts; k++) { /* Get the particle */ @@ -1265,6 +1286,9 @@ void space_sparts_get_cell_index_mapper(void *map_data, int nr_sparts, pos_z); #endif + /* Compute minimal mass */ + min_mass = min(min_mass, sp->mass); + /* Update the position */ sp->x[0] = pos_x; sp->x[1] = pos_y; @@ -1275,11 +1299,16 @@ void space_sparts_get_cell_index_mapper(void *map_data, int nr_sparts, for (int k = 0; k < s->nr_cells; k++) if (cell_counts[k]) atomic_add(&data->cell_counts[k], cell_counts[k]); free(cell_counts); + + /* Write back thee minimal part mass */ + atomic_min_f(&s->min_spart_mass, min_mass); } /** * @brief Computes the cell index of all the particles. * + * Also computes the minimal mass of all #part. + * * @param s The #space. * @param ind The array of indices to fill. * @param cell_counts The cell counters to update. @@ -1291,6 +1320,9 @@ void space_parts_get_cell_index(struct space *s, int *ind, int *cell_counts, const ticks tic = getticks(); + /* Re-set the minimal mass */ + s->min_part_mass = FLT_MAX; + /* Pack the extra information */ struct index_data data; data.s = s; @@ -1309,6 +1341,8 @@ void space_parts_get_cell_index(struct space *s, int *ind, int *cell_counts, /** * @brief Computes the cell index of all the g-particles. * + * Also computes the minimal mass of all dark-matter #gpart. + * * @param s The #space. * @param gind The array of indices to fill. * @param cell_counts The cell counters to update. @@ -1320,6 +1354,9 @@ void space_gparts_get_cell_index(struct space *s, int *gind, int *cell_counts, const ticks tic = getticks(); + /* Re-set the minimal mass */ + s->min_gpart_mass = FLT_MAX; + /* Pack the extra information */ struct index_data data; data.s = s; @@ -1338,6 +1375,8 @@ void space_gparts_get_cell_index(struct space *s, int *gind, int *cell_counts, /** * @brief Computes the cell index of all the s-particles. * + * Also computes the minimal mass of all #spart. + * * @param s The #space. * @param sind The array of indices to fill. * @param cell_counts The cell counters to update. @@ -1349,6 +1388,9 @@ void space_sparts_get_cell_index(struct space *s, int *sind, int *cell_counts, const ticks tic = getticks(); + /* Re-set the minimal mass */ + s->min_spart_mass = FLT_MAX; + /* Pack the extra information */ struct index_data data; data.s = s; diff --git a/src/space.h b/src/space.h index e8ad1929fa3902221ac24de790547e7e9aa5163d..2f17854e38cf26c24fb463bddb25d995ec144b69 100644 --- a/src/space.h +++ b/src/space.h @@ -146,6 +146,15 @@ struct space { /*! The top-level FFT task */ struct task *grav_top_level; + /*! Minimal mass of all the #part */ + float min_part_mass; + + /*! Minimal mass of all the dark-matter #gpart */ + float min_gpart_mass; + + /*! Minimal mass of all the #spart */ + float min_spart_mass; + /*! General-purpose lock for this space. */ swift_lock_type lock;