diff --git a/examples/EAGLE_low_z/EAGLE_6/eagle_6.yml b/examples/EAGLE_low_z/EAGLE_6/eagle_6.yml index fb885fa44f119bb3e296947997f71581d5ca0f48..313a5a384324eecd56455ea22bbc96d147982d6b 100644 --- a/examples/EAGLE_low_z/EAGLE_6/eagle_6.yml +++ b/examples/EAGLE_low_z/EAGLE_6/eagle_6.yml @@ -87,3 +87,30 @@ EAGLECooling: He_reion_z_centre: 3.5 He_reion_z_sigma: 0.5 He_reion_eV_p_H: 2.0 +# EAGLE star formation parameters +EAGLEStarFormation: + EOS_density_norm_H_p_cm3: 0.1 # Physical density used for the normalisation of the EOS assumed for the star-forming gas in Hydrogen atoms per cm^3. + EOS_temperature_norm_K: 8000 # Temperature om the polytropic EOS assumed for star-forming gas at the density normalisation in Kelvin. + EOS_gamma_effective: 1.3333333 # Slope the of the polytropic EOS assumed for the star-forming gas. + gas_fraction: 0.3 # The gas fraction used internally by the model. + KS_normalisation: 1.515e-4 # The normalization of the Kennicutt-Schmidt law in Msun / kpc^2 / yr. + KS_exponent: 1.4 # The exponent of the Kennicutt-Schmidt law. + KS_min_over_density: 57.7 # The over-density above which star-formation is allowed. + KS_high_density_threshold_H_p_cm3: 1e3 # Hydrogen number density above which the Kennicut-Schmidt law changes slope in Hydrogen atoms per cm^3. + KS_high_density_exponent: 2.0 # Slope of the Kennicut-Schmidt law above the high-density threshold. + KS_temperature_margin_dex: 0.5 # Logarithm base 10 of the maximal temperature difference above the EOS allowed to form stars. + threshold_norm_H_p_cm3: 0.1 # Normalisation of the metal-dependant density threshold for star formation in Hydrogen atoms per cm^3. + threshold_Z0: 0.002 # Reference metallicity (metal mass fraction) for the metal-dependant threshold for star formation. + threshold_slope: -0.64 # Slope of the metal-dependant star formation threshold + threshold_max_density_H_p_cm3: 10.0 # Maximal density of the metal-dependant density threshold for star formation in Hydrogen atoms per cm^3. + +# Parameters for the EAGLE "equation of state" +EAGLEEntropyFloor: + Jeans_density_threshold_H_p_cm3: 0.1 # Physical density above which the EAGLE Jeans limiter entropy floor kicks in expressed in Hydrogen atoms per cm^3. + Jeans_over_density_threshold: 10. # Overdensity above which the EAGLE Jeans limiter entropy floor can kick in. + Jeans_temperature_norm_K: 8000 # Temperature of the EAGLE Jeans limiter entropy floor at the density threshold expressed in Kelvin. + Jeans_gamma_effective: 1.3333333 # Slope the of the EAGLE Jeans limiter entropy floor + Cool_density_threshold_H_p_cm3: 1e-5 # Physical density above which the EAGLE Cool limiter entropy floor kicks in expressed in Hydrogen atoms per cm^3. + Cool_over_density_threshold: 10. # Overdensity above which the EAGLE Cool limiter entropy floor can kick in. + Cool_temperature_norm_K: 8000 # Temperature of the EAGLE Cool limiter entropy floor at the density threshold expressed in Kelvin. + Cool_gamma_effective: 1. # Slope the of the EAGLE Cool limiter entropy floor diff --git a/src/active.h b/src/active.h index 89955ae430aa0a685cdc513b5c11ac44f6b093e8..4bd3ab3bd6f18e152b9dafc76215ec3e008a415a 100644 --- a/src/active.h +++ b/src/active.h @@ -412,7 +412,6 @@ __attribute__((always_inline)) INLINE static int cell_is_starting_stars( return (c->stars.ti_beg_max == e->ti_current); } - /** * @brief Is this particle starting its time-step now ? * diff --git a/src/cell.c b/src/cell.c index 8f63a3c51492cd3b0f65a5d55f79fe54aa4e293b..de9d16f8639ded787b5a7c0b22042d796d6d057e 100644 --- a/src/cell.c +++ b/src/cell.c @@ -3404,7 +3404,7 @@ int cell_unskip_gravity_tasks(struct cell *c, struct scheduler *s) { * @return 1 If the space needs rebuilding. 0 otherwise. */ int cell_unskip_stars_tasks(struct cell *c, struct scheduler *s) { - + struct engine *e = s->space->e; const int nodeID = e->nodeID; int rebuild = 0; diff --git a/src/cell.h b/src/cell.h index 529993a7433385656a0964867ddae1bf1f0a6e7c..40d525c7fd49a5ffd858fc5b28e8ca69b17d6b3c 100644 --- a/src/cell.h +++ b/src/cell.h @@ -209,7 +209,7 @@ struct pcell_step { /*! Maximal integer end-of-timestep in this cell (stars) */ integertime_t ti_end_max; - + /*! Maximal distance any #part has travelled since last rebuild */ float dx_max_part; @@ -554,7 +554,7 @@ struct cell { /*! Maximum beginning of (integer) time step in this cell for star tasks. */ integertime_t ti_beg_max; - + /*! Number of #spart updated in this cell. */ int updated; diff --git a/src/engine_maketasks.c b/src/engine_maketasks.c index d88bd304d0f99c13a040ce926215244555161b03..af7ecf1e9c930b1bc5e6f6a4e4bb4cb2608346fb 100644 --- a/src/engine_maketasks.c +++ b/src/engine_maketasks.c @@ -225,21 +225,21 @@ void engine_addtasks_send_stars(struct engine *e, struct cell *ci, /* Check if any of the density tasks are for the target node. */ for (l = ci->stars.density; l != NULL; l = l->next) if (l->t->ci->nodeID == nodeID || - (l->t->cj != NULL && l->t->cj->nodeID == nodeID)) + (l->t->cj != NULL && l->t->cj->nodeID == nodeID)) break; /* If so, attach send tasks. */ if (l != NULL) { - + if (t_feedback == NULL) { - + /* Make sure this cell is tagged. */ cell_ensure_tagged(ci); /* Create the tasks and their dependencies? */ t_feedback = scheduler_addtask(s, task_type_send, task_subtype_spart, - ci->mpi.tag, 0, ci, cj); - + ci->mpi.tag, 0, ci, cj); + /* The send_stars task should unlock the super_cell's kick task. */ scheduler_addunlock(s, t_feedback, ci->hydro.super->stars.stars_out); @@ -407,7 +407,7 @@ void engine_addtasks_recv_hydro(struct engine *e, struct cell *c, for (struct link *l = c->stars.density; l != NULL; l = l->next) { scheduler_addunlock(s, t_rho, l->t); } - + /* Recurse? */ if (c->split) for (int k = 0; k < 8; k++) @@ -443,7 +443,7 @@ void engine_addtasks_recv_stars(struct engine *e, struct cell *c, /* Create the tasks. */ t_feedback = scheduler_addtask(s, task_type_recv, task_subtype_spart, - c->mpi.tag, 0, c, NULL); + c->mpi.tag, 0, c, NULL); } c->mpi.stars.recv = t_feedback; @@ -451,7 +451,7 @@ void engine_addtasks_recv_stars(struct engine *e, struct cell *c, for (struct link *l = c->stars.density; l != NULL; l = l->next) { scheduler_addunlock(s, l->t, t_feedback); } - + for (struct link *l = c->stars.feedback; l != NULL; l = l->next) { scheduler_addunlock(s, t_feedback, l->t); } @@ -866,13 +866,13 @@ void engine_make_hierarchical_tasks_hydro(struct engine *e, struct cell *c) { /* Stars */ if (with_stars) { - c->stars.drift = scheduler_addtask(s, task_type_drift_spart, + c->stars.drift = scheduler_addtask(s, task_type_drift_spart, task_subtype_none, 0, 0, c, NULL); - if(!with_feedback) { - scheduler_addunlock(s, c->stars.drift, c->super->kick2); - } + if (!with_feedback) { + scheduler_addunlock(s, c->stars.drift, c->super->kick2); + } } - + /* Subgrid tasks: cooling */ if (with_cooling) { @@ -1780,20 +1780,19 @@ void engine_make_extra_hydroloop_tasks_mapper(void *map_data, int num_elements, #endif if (with_feedback) { - scheduler_addunlock(sched, ci->hydro.super->stars.sorts, - t_star_density); - scheduler_addunlock(sched, ci->hydro.super->hydro.sorts, - t_star_density); - - if (ci->hydro.super != cj->hydro.super) { + scheduler_addunlock(sched, ci->hydro.super->stars.sorts, + t_star_density); + scheduler_addunlock(sched, ci->hydro.super->hydro.sorts, + t_star_density); + + if (ci->hydro.super != cj->hydro.super) { scheduler_addunlock(sched, cj->hydro.super->stars.sorts, t_star_density); scheduler_addunlock(sched, cj->hydro.super->hydro.sorts, t_star_density); - - } + } } - + if (ci->nodeID == nodeID) { scheduler_addunlock(sched, t_force, ci->hydro.super->hydro.end_force); @@ -2027,18 +2026,18 @@ void engine_make_extra_hydroloop_tasks_mapper(void *map_data, int num_elements, #endif if (with_feedback) { - scheduler_addunlock(sched, ci->hydro.super->stars.sorts, - t_star_density); - scheduler_addunlock(sched, ci->hydro.super->hydro.sorts, - t_star_density); - if (ci->hydro.super != cj->hydro.super) { - scheduler_addunlock(sched, cj->hydro.super->stars.sorts, - t_star_density); - scheduler_addunlock(sched, cj->hydro.super->hydro.sorts, - t_star_density); - } - } - + scheduler_addunlock(sched, ci->hydro.super->stars.sorts, + t_star_density); + scheduler_addunlock(sched, ci->hydro.super->hydro.sorts, + t_star_density); + if (ci->hydro.super != cj->hydro.super) { + scheduler_addunlock(sched, cj->hydro.super->stars.sorts, + t_star_density); + scheduler_addunlock(sched, cj->hydro.super->hydro.sorts, + t_star_density); + } + } + if (ci->nodeID == nodeID) { scheduler_addunlock(sched, t_force, ci->hydro.super->hydro.end_force); diff --git a/src/engine_marktasks.c b/src/engine_marktasks.c index a2f8b9a506f99fb41c75bc18f80b7ef2c4a146b7..e3b53861dd3c39137c380a45eba3bb4f55924562 100644 --- a/src/engine_marktasks.c +++ b/src/engine_marktasks.c @@ -73,16 +73,16 @@ void engine_activate_stars_mpi(struct engine *e, struct scheduler *s, scheduler_activate(s, ci->mpi.hydro.recv_xv); /* If the local cell is active, more stuff will be needed. */ - //MATTHIEU - //scheduler_activate_send(s, cj->mpi.stars.send, ci_nodeID); + // MATTHIEU + // scheduler_activate_send(s, cj->mpi.stars.send, ci_nodeID); /* If the local cell is active, send its ti_end values. */ scheduler_activate_send(s, cj->mpi.send_ti, ci_nodeID); } if (ci_active_stars) { - //MATTHIEU - //scheduler_activate(s, ci->mpi.stars.recv); + // MATTHIEU + // scheduler_activate(s, ci->mpi.stars.recv); /* If the foreign cell is active, we want its ti_end values. */ scheduler_activate(s, ci->mpi.recv_ti); @@ -104,16 +104,16 @@ void engine_activate_stars_mpi(struct engine *e, struct scheduler *s, scheduler_activate(s, cj->mpi.hydro.recv_xv); /* If the local cell is active, more stuff will be needed. */ - //MATTHIEU - //scheduler_activate_send(s, ci->mpi.stars.send, cj_nodeID); + // MATTHIEU + // scheduler_activate_send(s, ci->mpi.stars.send, cj_nodeID); /* If the local cell is active, send its ti_end values. */ scheduler_activate_send(s, ci->mpi.send_ti, cj_nodeID); } if (cj_active_stars) { - //MATTHIEU - //scheduler_activate(s, cj->mpi.stars.recv); + // MATTHIEU + // scheduler_activate(s, cj->mpi.stars.recv); /* If the foreign cell is active, we want its ti_end values. */ scheduler_activate(s, cj->mpi.recv_ti); diff --git a/src/error.h b/src/error.h index e99d1c56fb99b51bfce98ca978b823dbe29794d3..181a32c78b4ef09708eeade002157157895ddde5 100644 --- a/src/error.h +++ b/src/error.h @@ -54,7 +54,7 @@ extern int engine_rank; fprintf(stderr, "[%04i] %s %s:%s():%i: " s "\n", engine_rank, \ clocks_get_timesincestart(), __FILE__, __FUNCTION__, __LINE__, \ ##__VA_ARGS__); \ - swift_abort(1); \ + swift_abort(1); \ }) #else #define error(s, ...) \ diff --git a/src/runner.c b/src/runner.c index 27d67ddd876b58c3b2ee5b3e4e3e90e7d8f97ac3..3e58f1c7e3b960d4689815eda612c5c31b4f6267 100644 --- a/src/runner.c +++ b/src/runner.c @@ -1816,7 +1816,9 @@ void runner_do_kick1(struct runner *r, struct cell *c, int timer) { TIMER_TIC; /* Anything to do here? */ - if (!cell_is_starting_hydro(c, e) && !cell_is_starting_gravity(c, e) && !cell_is_starting_stars(c, e)) return; + if (!cell_is_starting_hydro(c, e) && !cell_is_starting_gravity(c, e) && + !cell_is_starting_stars(c, e)) + return; /* Recurse? */ if (c->split) { @@ -1998,7 +2000,9 @@ void runner_do_kick2(struct runner *r, struct cell *c, int timer) { TIMER_TIC; /* Anything to do here? */ - if (!cell_is_active_hydro(c, e) && !cell_is_active_gravity(c, e) && !cell_is_active_stars(c, e)) return; + if (!cell_is_active_hydro(c, e) && !cell_is_active_gravity(c, e) && + !cell_is_active_stars(c, e)) + return; /* Recurse? */ if (c->split) { @@ -2386,13 +2390,13 @@ void runner_do_timestep(struct runner *r, struct cell *c, int timer) { ti_stars_end_min = min(ti_current + ti_new_step, ti_stars_end_min); ti_stars_end_max = max(ti_current + ti_new_step, ti_stars_end_max); - ti_gravity_end_min = min(ti_current + ti_new_step, ti_gravity_end_min); - ti_gravity_end_max = max(ti_current + ti_new_step, ti_gravity_end_max); + ti_gravity_end_min = min(ti_current + ti_new_step, ti_gravity_end_min); + ti_gravity_end_max = max(ti_current + ti_new_step, ti_gravity_end_max); /* What is the next starting point for this cell ? */ ti_stars_beg_max = max(ti_current, ti_stars_beg_max); - ti_gravity_beg_max = max(ti_current, ti_gravity_beg_max); - + ti_gravity_beg_max = max(ti_current, ti_gravity_beg_max); + /* star particle is inactive but not inhibited */ } else { @@ -2403,16 +2407,16 @@ void runner_do_timestep(struct runner *r, struct cell *c, int timer) { get_integer_time_end(ti_current, sp->time_bin); const integertime_t ti_beg = - get_integer_time_begin(ti_current + 1, sp->time_bin); - + get_integer_time_begin(ti_current + 1, sp->time_bin); + ti_stars_end_min = min(ti_end, ti_stars_end_min); ti_stars_end_max = max(ti_end, ti_stars_end_max); - ti_gravity_end_min = min(ti_end, ti_gravity_end_min); - ti_gravity_end_max = max(ti_end, ti_gravity_end_max); - + ti_gravity_end_min = min(ti_end, ti_gravity_end_min); + ti_gravity_end_max = max(ti_end, ti_gravity_end_max); + /* What is the next starting point for this cell ? */ ti_stars_beg_max = max(ti_beg, ti_stars_beg_max); - ti_gravity_beg_max = max(ti_beg, ti_gravity_beg_max); + ti_gravity_beg_max = max(ti_beg, ti_gravity_beg_max); } } } else { @@ -2421,7 +2425,7 @@ void runner_do_timestep(struct runner *r, struct cell *c, int timer) { for (int k = 0; k < 8; k++) { if (c->progeny[k] != NULL) { struct cell *restrict cp = c->progeny[k]; - + /* Recurse */ runner_do_timestep(r, cp, 0); diff --git a/src/space.c b/src/space.c index 08480866395436f272c9956b8a36421f7a97fc01..a9bf4ceb832a92d085853895335e5b5428330961 100644 --- a/src/space.c +++ b/src/space.c @@ -331,9 +331,11 @@ void space_regrid(struct space *s, int verbose) { const struct cell *c = &s->cells_top[s->local_cells_with_particles_top[k]]; if (c->hydro.h_max > h_max) { + message("hydro h_max=%e", c->hydro.h_max); h_max = c->hydro.h_max; } if (c->stars.h_max > h_max) { + message("stars h_max=%e", c->stars.h_max); h_max = c->stars.h_max; } } @@ -2773,8 +2775,8 @@ void space_split_recursive(struct space *s, struct cell *c, ti_gravity_end_max = max(ti_gravity_end_max, cp->grav.ti_end_max); ti_gravity_beg_max = max(ti_gravity_beg_max, cp->grav.ti_beg_max); ti_stars_end_min = min(ti_stars_end_min, cp->stars.ti_end_min); - ti_stars_end_max = min(ti_stars_end_max, cp->stars.ti_end_max); - ti_stars_beg_max = min(ti_stars_beg_max, cp->stars.ti_beg_max); + ti_stars_end_max = min(ti_stars_end_max, cp->stars.ti_end_max); + ti_stars_beg_max = min(ti_stars_beg_max, cp->stars.ti_beg_max); /* Increase the depth */ if (cp->maxdepth > maxdepth) maxdepth = cp->maxdepth;