From 4d92cc87eb321ea236b948dfc4a69c0d7ee25517 Mon Sep 17 00:00:00 2001
From: Matthieu Schaller <schaller@strw.leidenuniv.nl>
Date: Wed, 15 May 2019 22:22:18 +0100
Subject: [PATCH] Made a separate function for the seeding of black holes.
---
examples/EAGLE_low_z/EAGLE_6/eagle_6.yml | 23 +++--
examples/main.c | 2 +-
examples/main_fof.c | 2 +-
src/engine.c | 5 +-
src/fof.c | 106 ++++++++++++++++-------
src/fof.h | 15 ++--
6 files changed, 106 insertions(+), 47 deletions(-)
diff --git a/examples/EAGLE_low_z/EAGLE_6/eagle_6.yml b/examples/EAGLE_low_z/EAGLE_6/eagle_6.yml
index 1d2cf42f6a..f726ebaeda 100644
--- a/examples/EAGLE_low_z/EAGLE_6/eagle_6.yml
+++ b/examples/EAGLE_low_z/EAGLE_6/eagle_6.yml
@@ -64,13 +64,14 @@ SPH:
# Parameters for the Friends-Of-Friends algorithm
FOF:
- basename: fof_output # Filename for the FOF outputs.
- min_group_size: 20 # The minimum no. of particles required for a group. Defaults to 20 if unspecified.
- linking_length_ratio: 0.2 # Linking length in units of the main inter-particle separation.
- absolute_linking_length: -1. # (Optional) Absolute linking length (in internal units)
- run_freq: 100 # (Optional) The no. of steps between each FOF search. Defaults to 2000.
- group_id_default: 2147483647 # (Optional) Sets the group ID of particles in groups below the minimum size. Defaults to 2^31 - 1 if unspecified. Has to be positive.
- group_id_offset: 1 # (Optional) Sets the offset of group ID labeling. Defaults to 1 if unspecified.
+ basename: fof_output # Filename for the FOF outputs.
+ min_group_size: 256 # The minimum no. of particles required for a group.
+ linking_length_ratio: 0.2 # Linking length in units of the main inter-particle separation.
+ black_hole_seed_halo_mass: 1.5e10 # Minimal halo mass in which to seed a black hole (in solar masses).
+ absolute_linking_length: -1. # (Optional) Absolute linking length (in internal units)
+ run_freq: 100 # (Optional) The no. of steps between each FOF search.
+ group_id_default: 2147483647 # (Optional) Sets the group ID of particles in groups below the minimum size. Defaults to 2^31 - 1 if unspecified. Has to be positive.
+ group_id_offset: 1 # (Optional) Sets the offset of group ID labeling. Defaults to 1 if unspecified.
# Parameters related to the initial conditions
InitialConditions:
@@ -129,3 +130,11 @@ EAGLEEntropyFloor:
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
+
+# EAGLE AGN model
+EAGLEAGN:
+ max_eddington_fraction: 1. # Maximal allowed accretion rate in units of the Eddington rate.
+ radiative_efficiency: 0.1 # Fraction of the accreted mass that gets radiated.
+ coupling_efficiency: 0.15 # Fraction of the radiated energy that couples to the gas in feedback events.
+ AGN_delta_T_K: 3.16228e8 # Change in temperature to apply to the gas particle in an AGN feedback event in Kelvin.
+ AGN_num_ngb_to_heat: 1. # Target number of gas neighbours to heat in an AGN feedback event.
diff --git a/examples/main.c b/examples/main.c
index b334b128af..caf70fbe79 100644
--- a/examples/main.c
+++ b/examples/main.c
@@ -831,7 +831,7 @@ int main(int argc, char *argv[]) {
/* Initialise the FOF properties */
bzero(&fof_properties, sizeof(struct fof_props));
- if (with_fof) fof_init(&fof_properties, params);
+ if (with_fof) fof_init(&fof_properties, params, &prog_const, &us);
/* Be verbose about what happens next */
if (myrank == 0) message("Reading ICs from file '%s'", ICfileName);
diff --git a/examples/main_fof.c b/examples/main_fof.c
index 2e8aa0e930..eeeca3ba96 100644
--- a/examples/main_fof.c
+++ b/examples/main_fof.c
@@ -916,7 +916,7 @@ int main(int argc, char *argv[]) {
/* Initialise the FOF properties */
bzero(&fof_properties, sizeof(struct fof_props));
- if (with_fof) fof_init(&fof_properties, params);
+ if (with_fof) fof_init(&fof_properties, params, &prog_const, &us);
/* Construct the engine policy */
int engine_policies = ENGINE_POLICY | engine_policy_steal;
diff --git a/src/engine.c b/src/engine.c
index 590daef4e1..b9bb4aadc6 100644
--- a/src/engine.c
+++ b/src/engine.c
@@ -4874,6 +4874,7 @@ void engine_unpin(void) {
* @param Ngas total number of gas particles in the simulation.
* @param Ngparts total number of gravity particles in the simulation.
* @param Nstars total number of star particles in the simulation.
+ * @param Nblackholes total number of black holes in the simulation.
* @param policy The queuing policy to use.
* @param verbose Is this #engine talkative ?
* @param reparttype What type of repartition algorithm are we using ?
@@ -4891,6 +4892,7 @@ void engine_unpin(void) {
* @param cooling_func The properties of the cooling function.
* @param starform The #star_formation model of this run.
* @param chemistry The chemistry information.
+ * @param fof_properties The #fof_props.
*/
void engine_init(struct engine *e, struct space *s, struct swift_params *params,
long long Ngas, long long Ngparts, long long Nstars,
@@ -6304,7 +6306,8 @@ void engine_fof(struct engine *e) {
/* Perform FOF search over foreign particles and
* find groups which require black hole seeding. */
- fof_search_tree(e->fof_properties, e->s);
+ fof_search_tree(e->fof_properties, e->s, /*dump_results=*/1,
+ /*seed_black_holes=*/1);
/* Reset flag. */
e->run_fof = 0;
diff --git a/src/fof.c b/src/fof.c
index 470cfcdf6a..f1563a97ea 100644
--- a/src/fof.c
+++ b/src/fof.c
@@ -50,6 +50,15 @@
#define UNION_BY_SIZE_OVER_MPI (1)
#define FOF_COMPRESS_PATHS_MIN_LENGTH (2)
+/**
+ * @brief Properties of a group used for black hole seeding
+ */
+enum fof_halo_seeding_props {
+ fof_halo_has_no_gas = -1LL,
+ fof_halo_has_black_hole = -2LL,
+ fof_halo_has_too_low_mass = -3LL
+};
+
#ifdef WITH_MPI
/* MPI types used for communications */
@@ -68,8 +77,12 @@ size_t node_offset;
*
* @param props the #fof_props structure to fill.
* @param params the parameter file parser.
+ * @param phys_const The physical constants in internal units.
+ * @param us The internal unit system.
*/
-void fof_init(struct fof_props *props, struct swift_params *params) {
+void fof_init(struct fof_props *props, struct swift_params *params,
+ const struct phys_const *phys_const,
+ const struct unit_system *us) {
/* Base name for the FOF output file */
parser_get_param_string(params, "FOF:basename", props->base_name);
@@ -112,6 +125,13 @@ void fof_init(struct fof_props *props, struct swift_params *params) {
if (props->l_x_ratio <= 0. && props->l_x_ratio != -1.)
error("The FOF linking length can't be negative!");
+ /* Read the minimal halo mass for black hole seeding */
+ props->seed_halo_mass =
+ parser_get_param_double(params, "FOF:black_hole_seed_halo_mass");
+
+ /* Convert to internal units */
+ props->seed_halo_mass *= phys_const->const_solar_mass;
+
#if defined(WITH_MPI) && defined(UNION_BY_SIZE_OVER_MPI)
if (engine_rank == 0)
message(
@@ -624,8 +644,8 @@ __attribute__((always_inline)) INLINE static void fof_compute_send_recv_offsets(
* @brief Perform a FOF search using union-find on a given leaf-cell
*
* @param props The properties fof the FOF scheme.
- * @param s The #space where the particles are.
* @param l_x2 The square of the FOF linking length.
+ * @param space_gparts The start of the #gpart array in the #space structure.
* @param c The #cell in which to perform FOF.
*/
void fof_search_self_cell(const struct fof_props *props, const double l_x2,
@@ -702,9 +722,10 @@ void fof_search_self_cell(const struct fof_props *props, const double l_x2,
* @brief Perform a FOF search using union-find between two cells
*
* @param props The properties fof the FOF scheme.
+ * @param dim The dimension of the simulation volume.
* @param l_x2 The square of the FOF linking length.
* @param periodic Are we using periodic BCs?
- * @param gparts The start of the #gpart array in the #space structure.
+ * @param space_gparts The start of the #gpart array in the #space structure.
* @param ci The first #cell in which to perform FOF.
* @param cj The second #cell in which to perform FOF.
*/
@@ -954,7 +975,7 @@ void fof_search_pair_cells_foreign(
* @param dim The dimension of the space.
* @param search_r2 the square of the FOF linking length.
* @param periodic Are we using periodic BCs?
- * @param gparts The start of the #gpart array in the #space structure.
+ * @param space_gparts The start of the #gpart array in the #space structure.
* @param ci The first #cell in which to perform FOF.
* @param cj The second #cell in which to perform FOF.
*/
@@ -1069,9 +1090,10 @@ void rec_fof_search_pair_foreign(
* @brief Recursively perform a union-find FOF on a cell.
*
* @param props The properties fof the FOF scheme.
- * @param s The #space where the particles are.
* @param dim The dimension of the space.
+ * @param space_gparts The start of the #gpart array in the #space structure.
* @param search_r2 the square of the FOF linking length.
+ * @param periodic Are we using periodic BCs?
* @param c The #cell in which to perform FOF.
*/
void rec_fof_search_self(const struct fof_props *props, const double dim[3],
@@ -1237,7 +1259,7 @@ void fof_unpack_group_mass_mapper(hashmap_key_t key, hashmap_value_t *value,
}
#endif
-/*
+/**
* @brief Calculates the total mass of each group above min_group_size and finds
* the densest particle for black hole seeding.
*/
@@ -1585,28 +1607,36 @@ void fof_calc_group_mass(struct fof_props *props, const struct space *s,
/* JSW TODO: Parallelise with threadpool*/
for (size_t i = 0; i < nr_gparts; i++) {
- /* Only check groups above the minimum size and mass threshold. */
- if (gparts[i].group_id != group_id_default &&
- group_mass[gparts[i].group_id - group_id_offset] > seed_halo_mass) {
+ const size_t index = gparts[i].group_id - group_id_offset;
+
+ /* Only check groups above the minimum mass threshold. */
+ if (gparts[i].group_id != group_id_default) {
+
+ if (group_mass[index] > seed_halo_mass) {
- const size_t index = gparts[i].group_id - group_id_offset;
+ /* Find the densest gas particle.
+ * Account for groups that already have a black hole and groups that
+ * contain no gas. */
+ if (gparts[i].type == swift_type_gas &&
+ max_part_density_index[index] != fof_halo_has_black_hole) {
- /* Find the densest gas particle.
- * Account for groups that already have a black hole and groups that
- * contain no gas. */
- if (gparts[i].type == swift_type_gas &&
- max_part_density_index[index] != fof_halo_has_black_hole) {
+ const size_t gas_index = -gparts[i].id_or_neg_offset;
+ const float rho_com = hydro_get_comoving_density(&parts[gas_index]);
- /* Update index if a denser gas particle is found. */
- if (parts[-gparts[i].id_or_neg_offset].rho > max_part_density[index]) {
- max_part_density_index[index] = -gparts[i].id_or_neg_offset;
- max_part_density[index] = parts[-gparts[i].id_or_neg_offset].rho;
+ /* Update index if a denser gas particle is found. */
+ if (rho_com > max_part_density[index]) {
+ max_part_density_index[index] = gas_index;
+ max_part_density[index] = rho_com;
+ }
}
- }
- /* If there is already a black hole in the group we don't need to create a
- new one. */
- else if (gparts[i].type == swift_type_black_hole) {
- max_part_density_index[index] = fof_halo_has_black_hole;
+ /* If there is already a black hole in the group we don't need to create
+ a new one. */
+ else if (gparts[i].type == swift_type_black_hole) {
+ max_part_density_index[index] = fof_halo_has_black_hole;
+ }
+
+ } else {
+ max_part_density_index[index] = fof_halo_has_too_low_mass;
}
}
}
@@ -1727,6 +1757,9 @@ void fof_find_foreign_links_mapper(void *map_data, int num_elements,
#endif
}
+void fof_seed_black_holes(const struct fof_props *props, struct space *s,
+ int num_groups, struct group_length *group_sizes) {}
+
/* Dump FOF group data. */
void fof_dump_group_data(const struct fof_props *props,
const char *out_file_name, struct space *s,
@@ -1802,6 +1835,7 @@ void fof_dump_group_data(const struct fof_props *props,
* @brief Search foreign cells for links and communicate any found to the
* appropriate node.
*
+ * @param props the properties of the FOF scheme.
* @param s Pointer to a #space.
*/
void fof_search_foreign_cells(struct fof_props *props, const struct space *s) {
@@ -2213,9 +2247,17 @@ void fof_search_foreign_cells(struct fof_props *props, const struct space *s) {
#endif /* WITH_MPI */
}
-/* Perform a FOF search on gravity particles using the cells and applying the
- * Union-Find algorithm.*/
-void fof_search_tree(struct fof_props *props, struct space *s) {
+/**
+ * @brief Perform a FOF search on gravity particles using the cells and applying
+ * the Union-Find algorithm.
+ *
+ * @param props The properties of the FOF scheme.
+ * @param s The #space containing the particles.
+ * @param dump_results Do we want to write the group catalogue to a file?
+ * @param seed_black_holes Do we want to seed black holes in haloes?
+ */
+void fof_search_tree(struct fof_props *props, struct space *s,
+ const int dump_results, const int seed_black_holes) {
const size_t nr_gparts = s->nr_gparts;
const size_t min_group_size = props->min_group_size;
@@ -2531,8 +2573,14 @@ void fof_search_tree(struct fof_props *props, struct space *s) {
clocks_from_ticks(getticks() - tic_seeding), clocks_getunit());
/* Dump group data. */
- fof_dump_group_data(props, output_file_name, s, num_groups_local,
- high_group_sizes);
+ if (dump_results) {
+ fof_dump_group_data(props, output_file_name, s, num_groups_local,
+ high_group_sizes);
+ }
+
+ if (seed_black_holes) {
+ fof_seed_black_holes(props, s, num_groups_local, high_group_sizes);
+ }
/* Free the left-overs */
swift_free("fof_high_group_sizes", high_group_sizes);
diff --git a/src/fof.h b/src/fof.h
index 7ff19219e5..c30782afa2 100644
--- a/src/fof.h
+++ b/src/fof.h
@@ -26,16 +26,12 @@
#include "align.h"
#include "parser.h"
-enum fof_halo_seeding_props {
- fof_halo_has_no_gas = -1LL,
- fof_halo_has_black_hole = -2LL,
- fof_halo_has_too_low_mass = -3LL
-};
-
/* Avoid cyclic inclusions */
struct gpart;
struct space;
struct engine;
+struct unit_system;
+struct phys_const;
/* MPI message required for FOF. */
struct fof_mpi {
@@ -163,11 +159,14 @@ struct cell_pair_indices {
} SWIFT_STRUCT_ALIGN;
/* Function prototypes. */
-void fof_init(struct fof_props *props, struct swift_params *params);
+void fof_init(struct fof_props *props, struct swift_params *params,
+ const struct phys_const *phys_const,
+ const struct unit_system *us);
void fof_create_mpi_types(void);
void fof_allocate(const struct space *s, const long long total_nr_DM_particles,
struct fof_props *props);
-void fof_search_tree(struct fof_props *props, struct space *s);
+void fof_search_tree(struct fof_props *props, struct space *s,
+ const int dump_results, const int seed_black_holes);
void rec_fof_search_self(const struct fof_props *props, const double dim[3],
const double search_r2, const int periodic,
const struct gpart *const space_gparts,
--
GitLab