diff --git a/examples/EAGLE_ICs/EAGLE_12/eagle_12.yml b/examples/EAGLE_ICs/EAGLE_12/eagle_12.yml
index eccf8529098eea1e157619d07d3c588bb30ef289..ee28df4164cea519d2136080714d8121ab8ed1fc 100644
--- a/examples/EAGLE_ICs/EAGLE_12/eagle_12.yml
+++ b/examples/EAGLE_ICs/EAGLE_12/eagle_12.yml
@@ -165,6 +165,7 @@ EAGLEFeedback:
EAGLEAGN:
subgrid_seed_mass_Msun: 1.5e5 # Black hole subgrid mass at creation time in solar masses.
max_eddington_fraction: 1. # Maximal allowed accretion rate in units of the Eddington rate.
+ viscuous_alpha: 1e6 # Normalisation constant of the Bondi viscuous time-scale accretion reduction term
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.
diff --git a/examples/EAGLE_ICs/EAGLE_25/eagle_25.yml b/examples/EAGLE_ICs/EAGLE_25/eagle_25.yml
index a960af6dd9ff949053c5a053ae29a7bd97e19d9c..230e5320b91968fd46ddb021cc68e1fddb3d9758 100644
--- a/examples/EAGLE_ICs/EAGLE_25/eagle_25.yml
+++ b/examples/EAGLE_ICs/EAGLE_25/eagle_25.yml
@@ -166,6 +166,7 @@ EAGLEFeedback:
EAGLEAGN:
subgrid_seed_mass_Msun: 1.5e5 # Black hole subgrid mass at creation time in solar masses.
max_eddington_fraction: 1. # Maximal allowed accretion rate in units of the Eddington rate.
+ viscuous_alpha: 1e6 # Normalisation constant of the Bondi viscuous time-scale accretion reduction term
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.
diff --git a/examples/EAGLE_ICs/EAGLE_50/eagle_50.yml b/examples/EAGLE_ICs/EAGLE_50/eagle_50.yml
index b31a77d87669d5786ab1fdb904e1f6683b9dece9..9e95d625c6fccf3e84dbcef08033d24bf106e20f 100644
--- a/examples/EAGLE_ICs/EAGLE_50/eagle_50.yml
+++ b/examples/EAGLE_ICs/EAGLE_50/eagle_50.yml
@@ -166,6 +166,7 @@ EAGLEFeedback:
EAGLEAGN:
subgrid_seed_mass_Msun: 1.5e5 # Black hole subgrid mass at creation time in solar masses.
max_eddington_fraction: 1. # Maximal allowed accretion rate in units of the Eddington rate.
+ viscuous_alpha: 1e6 # Normalisation constant of the Bondi viscuous time-scale accretion reduction term
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.
diff --git a/examples/EAGLE_low_z/EAGLE_12/eagle_12.yml b/examples/EAGLE_low_z/EAGLE_12/eagle_12.yml
index b7009a7d3b06087d875c10aa1861a8456dd3f4e0..3ac39899c432f291f62107abbe4e1cca31c53c21 100644
--- a/examples/EAGLE_low_z/EAGLE_12/eagle_12.yml
+++ b/examples/EAGLE_low_z/EAGLE_12/eagle_12.yml
@@ -157,6 +157,7 @@ EAGLEFeedback:
EAGLEAGN:
subgrid_seed_mass_Msun: 1.5e5 # Black hole subgrid mass at creation time in solar masses.
max_eddington_fraction: 1. # Maximal allowed accretion rate in units of the Eddington rate.
+ viscuous_alpha: 1e6 # Normalisation constant of the Bondi viscuous time-scale accretion reduction term
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.
diff --git a/examples/EAGLE_low_z/EAGLE_25/eagle_25.yml b/examples/EAGLE_low_z/EAGLE_25/eagle_25.yml
index 5a753cdf0da41dd4ec102eff351a83658c4761ee..9fb2243887625998b7e9a322aac0932720dfba0a 100644
--- a/examples/EAGLE_low_z/EAGLE_25/eagle_25.yml
+++ b/examples/EAGLE_low_z/EAGLE_25/eagle_25.yml
@@ -164,6 +164,7 @@ EAGLEFeedback:
EAGLEAGN:
subgrid_seed_mass_Msun: 1.5e5 # Black hole subgrid mass at creation time in solar masses.
max_eddington_fraction: 1. # Maximal allowed accretion rate in units of the Eddington rate.
+ viscuous_alpha: 1e6 # Normalisation constant of the Bondi viscuous time-scale accretion reduction term
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.
diff --git a/examples/EAGLE_low_z/EAGLE_50/eagle_50.yml b/examples/EAGLE_low_z/EAGLE_50/eagle_50.yml
index e36eb4ad8d814b7d37006578a2c2796f3ec934b9..d36d68bdf68810beb9edcdf97c95cedf02f818ac 100644
--- a/examples/EAGLE_low_z/EAGLE_50/eagle_50.yml
+++ b/examples/EAGLE_low_z/EAGLE_50/eagle_50.yml
@@ -156,6 +156,7 @@ EAGLEFeedback:
EAGLEAGN:
subgrid_seed_mass_Msun: 1.5e5 # Black hole subgrid mass at creation time in solar masses.
max_eddington_fraction: 1. # Maximal allowed accretion rate in units of the Eddington rate.
+ viscuous_alpha: 1e6 # Normalisation constant of the Bondi viscuous time-scale accretion reduction term
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.
diff --git a/examples/EAGLE_low_z/EAGLE_6/eagle_6.yml b/examples/EAGLE_low_z/EAGLE_6/eagle_6.yml
index cfe76f51efb0920ffdb092d6af7db0104db0ee85..b4cd37e079b5131dec72c7b240fc7a98eded6236 100644
--- a/examples/EAGLE_low_z/EAGLE_6/eagle_6.yml
+++ b/examples/EAGLE_low_z/EAGLE_6/eagle_6.yml
@@ -166,6 +166,7 @@ EAGLEFeedback:
EAGLEAGN:
subgrid_seed_mass_Msun: 1.5e5 # Black hole subgrid mass at creation time in solar masses.
max_eddington_fraction: 1. # Maximal allowed accretion rate in units of the Eddington rate.
+ viscuous_alpha: 1e6 # Normalisation constant of the Bondi viscuous time-scale accretion reduction term
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.
diff --git a/examples/parameter_example.yml b/examples/parameter_example.yml
index eec64eb7e7c03640e344ec37419b06e3b3c1b663..4c20710b52fda4bd1b71adb2c89a57eedc000e6d 100644
--- a/examples/parameter_example.yml
+++ b/examples/parameter_example.yml
@@ -396,6 +396,7 @@ EAGLEFeedback:
EAGLEAGN:
subgrid_seed_mass_Msun: 1.5e5 # Black hole subgrid mass at creation time in solar masses.
max_eddington_fraction: 1. # Maximal allowed accretion rate in units of the Eddington rate.
+ viscuous_alpha: 1e6 # Normalisation constant of the Bondi viscuous time-scale accretion reduction term
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.
diff --git a/src/black_holes/EAGLE/black_holes.h b/src/black_holes/EAGLE/black_holes.h
index ccb97e015b50046e10b4966548bf59d55d73c526..c359467f739368ec305b9ccfd8ec4e66c3467940 100644
--- a/src/black_holes/EAGLE/black_holes.h
+++ b/src/black_holes/EAGLE/black_holes.h
@@ -82,6 +82,9 @@ __attribute__((always_inline)) INLINE static void black_holes_init_bpart(
bp->velocity_gas[0] = 0.f;
bp->velocity_gas[1] = 0.f;
bp->velocity_gas[2] = 0.f;
+ bp->circular_velocity_gas[0] = 0.f;
+ bp->circular_velocity_gas[1] = 0.f;
+ bp->circular_velocity_gas[2] = 0.f;
bp->ngb_mass = 0.f;
bp->num_ngbs = 0;
}
@@ -136,6 +139,9 @@ __attribute__((always_inline)) INLINE static void black_holes_end_density(
bp->velocity_gas[0] *= h_inv_dim;
bp->velocity_gas[1] *= h_inv_dim;
bp->velocity_gas[2] *= h_inv_dim;
+ bp->circular_velocity_gas[0] *= h_inv_dim;
+ bp->circular_velocity_gas[1] *= h_inv_dim;
+ bp->circular_velocity_gas[2] *= h_inv_dim;
const float rho_inv = 1.f / bp->rho_gas;
@@ -144,6 +150,9 @@ __attribute__((always_inline)) INLINE static void black_holes_end_density(
bp->velocity_gas[0] *= rho_inv;
bp->velocity_gas[1] *= rho_inv;
bp->velocity_gas[2] *= rho_inv;
+ bp->circular_velocity_gas[0] *= rho_inv * h_inv;
+ bp->circular_velocity_gas[1] *= rho_inv * h_inv;
+ bp->circular_velocity_gas[2] *= rho_inv * h_inv;
}
/**
@@ -290,12 +299,13 @@ __attribute__((always_inline)) INLINE static void black_holes_prepare_feedback(
const double num_ngbs_to_heat = props->num_ngbs_to_heat;
const double delta_T = props->AGN_delta_T_desired;
const double delta_u = delta_T * props->temp_to_u_factor;
+ const double alpha_visc = props->alpha_visc;
/* (Subgrid) mass of the BH (internal units) */
const double BH_mass = bp->subgrid_mass;
/* Convert the quantities we gathered to physical frame (all internal units)
- */
+ * Note: for the velocities this means peculiar velocities */
const double gas_rho_phys = bp->rho_gas * cosmo->a3_inv;
const double gas_c_phys = bp->sound_speed_gas * cosmo->a_factor_sound_speed;
const double gas_v_peculiar[3] = {bp->velocity_gas[0] * cosmo->a_inv,
@@ -306,6 +316,11 @@ __attribute__((always_inline)) INLINE static void black_holes_prepare_feedback(
bp->v[1] * cosmo->a_inv,
bp->v[2] * cosmo->a_inv};
+ const double gas_v_circular[3] = {
+ bp->circular_velocity_gas[0] * cosmo->a_inv,
+ bp->circular_velocity_gas[1] * cosmo->a_inv,
+ bp->circular_velocity_gas[2] * cosmo->a_inv};
+
/* Difference in peculiar velocity between the gas and the BH
* Note that there is no need for a Hubble flow term here. We are
* computing the gas velocity at the position of the black hole. */
@@ -316,13 +331,31 @@ __attribute__((always_inline)) INLINE static void black_holes_prepare_feedback(
v_diff_peculiar[1] * v_diff_peculiar[1] +
v_diff_peculiar[2] * v_diff_peculiar[2];
+ /* Norm of the cirecular velocity of the gas around the BH */
+ const double tangential_velocity2 = gas_v_circular[0] * gas_v_circular[0] +
+ gas_v_circular[0] * gas_v_circular[0] +
+ gas_v_circular[0] * gas_v_circular[0];
+ const double tangential_velocity = sqrt(tangential_velocity2);
+
/* We can now compute the Bondi accretion rate (internal units) */
const double gas_c_phys2 = gas_c_phys * gas_c_phys;
const double denominator2 = v_diff_norm2 + gas_c_phys2;
const double denominator_inv = 1. / sqrt(denominator2);
- const double Bondi_rate = 4. * M_PI * G * G * BH_mass * BH_mass *
- gas_rho_phys * denominator_inv * denominator_inv *
- denominator_inv;
+ double Bondi_rate = 4. * M_PI * G * G * BH_mass * BH_mass * gas_rho_phys *
+ denominator_inv * denominator_inv * denominator_inv;
+
+ /* Compute the reduction factor from Rosas-Guevara et al. (2015) */
+ const double Bondi_radius = G * BH_mass / gas_c_phys2;
+ const double Bondi_time = Bondi_radius / gas_c_phys;
+ const double r_times_v_tang = Bondi_radius * tangential_velocity;
+ const double r_times_v_tang_3 =
+ r_times_v_tang * r_times_v_tang * r_times_v_tang;
+ const double viscous_time = 2. * M_PI * r_times_v_tang_3 /
+ (1e-6 * alpha_visc * G * G * BH_mass * BH_mass);
+ const double f_visc = max(Bondi_time / viscous_time, 1.);
+
+ /* Limit the Bondi rate by the Bondi viscuous time ratio */
+ Bondi_rate *= f_visc;
/* Compute the Eddington rate (internal units) */
const double Eddington_rate =
diff --git a/src/black_holes/EAGLE/black_holes_iact.h b/src/black_holes/EAGLE/black_holes_iact.h
index 27723f00c92e0eafd011cf29bdd017215c0a776c..955803786486e4f0edddd5df2a194685d2762ba9 100644
--- a/src/black_holes/EAGLE/black_holes_iact.h
+++ b/src/black_holes/EAGLE/black_holes_iact.h
@@ -84,6 +84,12 @@ __attribute__((always_inline)) INLINE static void runner_iact_nonsym_bh_density(
bi->velocity_gas[1] += mj * vj[1] * wi;
bi->velocity_gas[2] += mj * vj[2] * wi;
+ /* Contribution to the circular valocity */
+ const float dv[3] = {bi->v[0] - vj[0], bi->v[1] - vj[1], bi->v[2] - vj[2]};
+ bi->circular_velocity_gas[0] += mj * wi * (dx[1] * dv[2] - dx[2] * dv[1]);
+ bi->circular_velocity_gas[1] += mj * wi * (dx[0] * dv[2] - dx[2] * dv[0]);
+ bi->circular_velocity_gas[2] += mj * wi * (dx[0] * dv[1] - dx[1] * dv[0]);
+
#ifdef DEBUG_INTERACTIONS_BH
/* Update ngb counters */
if (si->num_ngb_density < MAX_NUM_OF_NEIGHBOURS_BH)
diff --git a/src/black_holes/EAGLE/black_holes_part.h b/src/black_holes/EAGLE/black_holes_part.h
index 45be4461ceea1beefaa84ce4705280b4a9b45491..598edffecae228a47f57224be9727ee50ccc4074 100644
--- a/src/black_holes/EAGLE/black_holes_part.h
+++ b/src/black_holes/EAGLE/black_holes_part.h
@@ -94,6 +94,9 @@ struct bpart {
/*! Smoothed velocity (peculiar) of the gas surrounding the black hole */
float velocity_gas[3];
+ /*! Curl of the velocity field around the black hole */
+ float circular_velocity_gas[3];
+
/*! Total mass of the gas neighbours. */
float ngb_mass;
diff --git a/src/black_holes/EAGLE/black_holes_properties.h b/src/black_holes/EAGLE/black_holes_properties.h
index 6cb0b4ea7056c19b7c131a9d226675b5c888b5ce..f2e17aef9f71abd6ef1f36946442d568e6e94bbe 100644
--- a/src/black_holes/EAGLE/black_holes_properties.h
+++ b/src/black_holes/EAGLE/black_holes_properties.h
@@ -63,6 +63,9 @@ struct black_holes_props {
/*! Feedback coupling efficiency of the black holes. */
float epsilon_f;
+ /*! Normalisation of the viscuous angular momentum accretion reduction */
+ float alpha_visc;
+
/* ---- Properties of the feedback model ------- */
/*! Temperature increase induced by AGN feedback (Kelvin) */
@@ -144,6 +147,7 @@ INLINE static void black_holes_props_init(struct black_holes_props *bp,
parser_get_param_float(params, "EAGLEAGN:radiative_efficiency");
bp->epsilon_f =
parser_get_param_float(params, "EAGLEAGN:coupling_efficiency");
+ bp->alpha_visc = parser_get_param_float(params, "EAGLEAGN:viscuous_alpha");
/* Feedback parameters ---------------------------------- */