diff --git a/src/black_holes/EAGLE/#black_holes_part.h# b/src/black_holes/EAGLE/#black_holes_part.h#
deleted file mode 100644
index 85b2599a0e545a9807353fcb7d637a9f36b6a0cf..0000000000000000000000000000000000000000
--- a/src/black_holes/EAGLE/#black_holes_part.h#
+++ /dev/null
@@ -1,92 +0,0 @@
-/*******************************************************************************
- * This file is part of SWIFT.
- * Copyright (c) 2019 Matthieu Schaller (schaller@strw.leidenuniv.nl)
- *
- * This program is free software: you can redistribute it and/or modify
- * it under the terms of the GNU Lesser General Public License as published
- * by the Free Software Foundation, either version 3 of the License, or
- * (at your option) any later version.
- *
- * This program is distributed in the hope that it will be useful,
- * but WITHOUT ANY WARRANTY; without even the implied warranty of
- * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
- * GNU General Public License for more details.
- *
- * You should have received a copy of the GNU Lesser General Public License
- * along with this program. If not, see <http://www.gnu.org/licenses/>.
- *
- ******************************************************************************/
-#ifndef SWIFT_EAGLE_BLACK_HOLE_PART_H
-#define SWIFT_EAGLE_BLACK_HOLE_PART_H
-
-/* Some standard headers. */
-#include <stdlib.h>
-
-/**
- * @brief Particle fields for the black hole particles.
- *
- * All quantities related to gravity are stored in the associate #gpart.
- */
-struct bpart {
-
- /*! Particle ID. */
- long long id;
-
- /*! Pointer to corresponding gravity part. */
- struct gpart* gpart;
-
- /*! Particle position. */
- double x[3];
-
- /* Offset between current position and position at last tree rebuild. */
- float x_diff[3];
-
- /*! Particle velocity. */
- float v[3];
-
- /*! Black hole mass */
- float mass;
-
- /* Particle cutoff radius. */
- float h;
-
- /*! Particle time bin */
- timebin_t time_bin;
-
- struct {
-
- /* Number of neighbours. */
- float wcount;
-
- /* Number of neighbours spatial derivative. */
- float wcount_dh;
-
- } density;
-
-#ifdef SWIFT_DEBUG_CHECKS
-
- /* Time of the last drift */
- integertime_t ti_drift;
-
- /* Time of the last kick */
- integertime_t ti_kick;
-
-#endif
-
-#ifdef DEBUG_INTERACTIONS_BLACK_HOLES
- /*! Number of interactions in the density SELF and PAIR */
- int num_ngb_density;
-
- /*! List of interacting particles in the density SELF and PAIR */
- long long ids_ngbs_density[MAX_NUM_OF_NEIGHBOURS_BLACK_HOLES];
-
- /*! Number of interactions in the force SELF and PAIR */
- int num_ngb_force;
-
- /*! List of interacting particles in the force SELF and PAIR */
- long long ids_ngbs_force[MAX_NUM_OF_NEIGHBOURS_BLACK_HOLES];
-#endif
-
-} SWIFT_STRUCT_ALIGN;
-
-#endif /* SWIFT_EAGLE_BLACK_HOLE_PART_H */
diff --git a/src/black_holes/EAGLE/.#black_holes_part.h b/src/black_holes/EAGLE/.#black_holes_part.h
deleted file mode 120000
index 76c963c625e58adde9803f0f530cc7ece95fff2a..0000000000000000000000000000000000000000
--- a/src/black_holes/EAGLE/.#black_holes_part.h
+++ /dev/null
@@ -1 +0,0 @@
-matthieu@laptop.2863:1556535537
\ No newline at end of file
diff --git a/src/black_holes/EAGLE/black_holes.h b/src/black_holes/EAGLE/black_holes.h
index 67333955e8216401573b0841414e617e5adc1079..e3e2d02cd051d3ecc75c1ddb63b7e5f68f51898a 100644
--- a/src/black_holes/EAGLE/black_holes.h
+++ b/src/black_holes/EAGLE/black_holes.h
@@ -127,6 +127,14 @@ __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;
+
+ const float rho_inv = 1.f / bp->rho_gas;
+
+ /* Finish the calculation by undoing the mass smoothing */
+ bp->sound_speed_gas *= rho_inv;
+ bp->velocity_gas[0] *= rho_inv;
+ bp->velocity_gas[1] *= rho_inv;
+ bp->velocity_gas[2] *= rho_inv;
}
/**
@@ -157,8 +165,16 @@ __attribute__((always_inline)) INLINE static void black_holes_prepare_feedback(
struct bpart* restrict bp, const struct phys_const* constants,
const struct cosmology* cosmo, const double dt) {
+ /* Gather some physical constants (all in internal units) */
const double G = constants->const_newton_G;
const double c = constants->const_speed_light_c;
+ const double proton_mass = constants->const_proton_mass;
+ const double sigma_Thomson = constants->const_thomson_cross_section;
+
+ /* Gather the parameters of the model */
+ const double f_Edd = 1.;
+ const double epsilon_r = 0.1;
+ const double epsilon_f = 0.15;
/* (Subgrid) mass of the BH (internal units) */
const double BH_mass = bp->subgrid_mass;
@@ -175,7 +191,7 @@ __attribute__((always_inline)) INLINE static void black_holes_prepare_feedback(
bp->v[1] * cosmo->a_inv,
bp->v[2] * cosmo->a_inv};
- /* Difference in velocity between the gas and the BH */
+ /* Difference in peculiar velocity between the gas and the BH */
const double v_diff_peculiar[3] = {gas_v_peculiar[0] - bh_v_peculiar[0],
gas_v_peculiar[1] - bh_v_peculiar[1],
gas_v_peculiar[2] - bh_v_peculiar[2]};
@@ -191,13 +207,18 @@ __attribute__((always_inline)) INLINE static void black_holes_prepare_feedback(
gas_rho_phys * denominator_inv * denominator_inv *
denominator_inv;
+ /* Compute the Eddington rate */
+ const double Eddington_rate =
+ 4. * M_PI * G * BH_mass * proton_mass / (epsilon_r * c * sigma_Thomson);
+
+ /* Limit the accretion rate to the Eddington fraction */
+ const double accr_rate = max(Bondi_rate, f_Edd * Eddington_rate);
+
/* Factor in the radiative efficiency */
- const double epsilon_r = 0.1;
- const double mass_rate = (1. - epsilon_r) * Bondi_rate;
- const double luminosity = epsilon_r * Bondi_rate * c * c;
+ const double mass_rate = (1. - epsilon_r) * accr_rate;
+ const double luminosity = epsilon_r * accr_rate * c * c;
/* Integrate forward in time */
- const double epsilon_f = 0.15;
bp->subgrid_mass += mass_rate * dt;
bp->energy_reservoir += luminosity * epsilon_f * dt;
}
diff --git a/src/black_holes/EAGLE/black_holes_iact.h b/src/black_holes/EAGLE/black_holes_iact.h
index c93fdc6de21141eebbe3f298f645d433d7cd2dbf..507a4a4f32df58fa011682dacfb28a6c2e2bfd15 100644
--- a/src/black_holes/EAGLE/black_holes_iact.h
+++ b/src/black_holes/EAGLE/black_holes_iact.h
@@ -66,15 +66,15 @@ __attribute__((always_inline)) INLINE static void runner_iact_nonsym_bh_density(
const float cj = hydro_get_comoving_soundspeed(pj);
/* Contribution to the smoothed sound speed */
- bi->sound_speed_gas += cj * wi;
+ bi->sound_speed_gas += mj * cj * wi;
/* Neighbour peculiar drifted velocity */
- const float v[3] = {pj->v[0], pj->v[1], pj->v[2]};
+ const float vj[3] = {pj->v[0], pj->v[1], pj->v[2]};
/* Contribution to the smoothed velocity */
- bi->velocity_gas[0] += v[0] * wi;
- bi->velocity_gas[1] += v[1] * wi;
- bi->velocity_gas[2] += v[2] * wi;
+ bi->velocity_gas[0] += mj * vj[0] * wi;
+ bi->velocity_gas[1] += mj * vj[1] * wi;
+ bi->velocity_gas[2] += mj * vj[2] * wi;
#ifdef DEBUG_INTERACTIONS_BH
/* Update ngb counters */