Search
Close
Search
Advanced Search
Search Menu
AI Discovery Assistant

Abstract

The dynamics of black hole (BH) seeds in high-redshift galaxies is key to understand their ability to grow via accretion and to pair in close binaries during galactic mergers. To properly follow the dynamics of BHs we develop a physically motivated model to capture unresolved dynamical friction from stars, dark matter, and gas. We first validate the model and then we use it to investigate the dynamics of seed BHs born at z ∼ 9 in dwarf proto-galaxies. We perform a suite of zoom cosmological simulations with spatial resolution as high as 10 pc and with a stellar and dark matter mass resolution of |$2\times 10^3 \, \, $| and |$2\times 10^5 \, \, \mathrm{ M}_{\odot }$|⁠, respectively. We first explore the dynamics of a seed BH in the galaxy where it is born and show that it is highly erratic if the seed mass is less than |$10^5\, \, \mathrm{ M}_{\odot }$|⁠. The dynamics is dominated by the stellar component, whose distribution is irregular and patchy, thus inducing stochasticity in the orbits: the BH may be anywhere in the proto-galaxy. When this dwarf merges into a larger galaxy, it is paramount to simulate the process with very high spatial and mass resolution in order to correctly account for the stripping of the stellar envelope of the satellite BH. The outcome of the encounter could be either a tight binary or, at least temporary, a wandering BH, leading to multiple BHs in a galaxy, each inherited from a different merger.

galaxies: evolution, galaxies: kinematics and dynamics

1 INTRODUCTION

The high-redshift Universe is the birthplace of the seeds of the supermassive black holes (BHs) observed in today’s galaxy centre (Kormendy & Ho 2013a). A variety of different physical mechanisms for seed formation have been proposed (Woods et al. 2018, and references therein), but observational constraints are hampered, since the seeds are predicted to have relatively low masses (⁠|$10^2\!-\!10^5 \, \, \mathrm{ M}_{\odot }$|⁠) and form at high redshift (z > 6), making their electromagnetic emission faint (Reines & Comastri 2016).

The seeds build-up their mass via accretion of gas and stars, or via mergers with other BHs (e.g. Volonteri, Haardt & Madau 2003). When BHs merge, they emit gravitational waves, and detection of such waves provides a complementary way of probing BH seeds (Sesana, Volonteri & Haardt 2007a; Barausse 2012; Dayal et al. 2018; Hartwig, Agarwal & Regan 2018; Ricarte & Natarajan 2018). For BHs with masses in the range |$10^4\!-\!10^7\, \, \mathrm{ M}_{\odot }$| the gravitational waves have frequency around mHz, and they are therefore primary targets for LISA, which can detect BHs with such masses out to z > 20 (Amaro-Seoane et al. 2017).

However, before coalescing by emission of gravitational waves, which can merge BHs of |$10^4\!-\!10^7\, \, \mathrm{ M}_{\odot }$| in less than a Hubble time once their separation is ∼10−4–10−2 pc, BHs have a long journey (Begelman, Blandford & Rees 1980). They are initially separated by tens of kpc and sit in the centre of separate galaxies, which eventually merge. Then, the long process of dynamical friction (Chandrasekhar 1943) begins, driving BHs towards the centre of the galaxy remnant, until they form a binary when their separation is pc-scale (e.g. Mayer et al. 2007; Pfister et al. 2017, and references therein). Once the binary has formed, scattering with stars (e.g. Quinlan 1996; Sesana, Haardt & Madau 2007b; Khan et al. 2012; Vasiliev, Antonini & Merritt 2015), interactions with massive or circumbinary discs (e.g. Dotti et al. 2007; Haiman, Kocsis & Menou 2009; Goicovic et al. 2016) or even three-body scattering with another incoming BH (see Bonetti et al. 2018, and references therein) are invoked to bridge the final gap to where emission of gravitational waves becomes efficient.

Cosmological simulations are excellent tools to study the properties of BH evolution over cosmic time, since they can track the joint evolution of BHs and of the galaxies they are embedded in (Tremmel et al. 2018a). Large-volume simulations provide good statistics, having a large number of galaxies and BHs in their boxes, but lack of mass and spatial resolution means that not even the formation of BH binaries can be resolved. Zoom simulations can have much higher resolution, but they allow for the study of a limited number of galaxies and BHs. In this paper, we present a model to better track the dynamics of BHs, validate it and show its limitations. We then use our model in high-resolution zoomed cosmological simulations to study the yet unexplored dynamics of BHs of mass |$10^4\!-\!10^5\, \mathrm{ M}_{\odot } ,$| in a cosmological context, primary targets for the LISA observatory (Amaro-Seoane et al. 2017) .

2 DYNAMICAL FRICTIgON IN NUMERICAL SIMULATIONS

Due to their high mass, BHs attract surrounding material, gas, stars, and dark matter, which create an overdensity lagging their passage. This overdensity drags and decelerates the moving BH: this phenomenon is referred to as dynamical friction (Chandrasekhar 1943; Chapon, Mayer & Teyssier 2013). To resolve the resulting force in a numerical simulation, Pfister et al. (2017) have shown that the spatial resolution, or the softening, should be smaller than the influence radius
(1)
where M• is the mass of the BH and σ is the velocity dispersion of material (gas, stars, or dark matter) around the BH. This is because the typical size of the drag, partly causing dynamical friction, has a typical size of the same order as rinf (Colpi, Mayer & Governato 1999). In cosmological simulations, the typical resolution is |${\sim } 100 \, \mathrm{ pc }\!-\!1 \, \mathrm{ kpc }$|⁠, much larger than the pc-scale needed to resolve rinf for a |$10^7 \, \mathrm{ M}_{\odot }$| BH in a Milky Way like galaxy. Therefore, we must remove by hand the momentum that a BH would lose through dynamical friction if we were able to resolve the phenomenon. In this section, we first describe how we implement unresolved dynamical friction in the adaptive mesh refinement code ramses (Teyssier 2002) for collisionless particles (stars and dark matter); the code already includes a correction for dynamical friction from gas (Dubois et al. 2012).

We follow an approach similar to Tremmel et al. (2015), although we include not only the contribution to dynamical friction from slow moving particles but also from fast-moving particles, which can play an important role when the density profile becomes shallow (Antonini & Merritt 2012; Dosopoulou & Antonini 2017).

We measure all the quantities needed to estimate dynamical friction in a sphere |$\mathcal {S}$| centred on the BH with a radius 4Δx, where Δx corresponds to the minimum grid size. We chose |$\mathcal {S}$| to be consistent with the already existing implementation for gas accretion, feedback, and dynamical friction (Dubois et al. 2012).

We report here equation (30) from Chandrasekhar (1943). This gives an analytical estimate of the amount of momentum that must be removed to BHs due to dynamical friction:
(2)
where we denote as M• the mass of the BH, as |$\vec{v_\bullet }$| (with magnitude v•) the relative velocity of the BH with respect to the velocity of the background, and |$\tilde{\vec{v}}$| defined below in equation (4); lnΛ = ln(4Δx/rdef) is the Coulomb logarithm (this expression is justified below); and f is the distribution function:
(3)
Here |$\vec{v_i}$| (with magnitude vi) is the relative velocity of particle i with respect to the velocity of the background, mi is the mass of particle i and δ is the Dirac function.
The velocity of the background, |$\tilde{\vec{v}}$|⁠, is simply the mass-weighted velocity of all particles (except the BH particle) enclosed in |$\mathcal {S}$|⁠:
(4)
where M is the total mass enclosed in |$\mathcal {S}$|⁠. We stress here that the background velocity is computed for stars and dark matter separately, the reason is that dynamical friction assumes particles with similar masses, which is a reasonable assumption if we consider an assembly of stars, and an assembly of dark matter particles, but not if we consider stars and dark matter particles together. Therefore we compute the contribution from dark matter, |$\vec{a}_{\rm DF,DM}$|⁠, and stars, |$\vec{a}_{\text{DF},\star }$| separately.
We justify here the expression above for the Coulomb logarithm lnΛ = ln(4Δx/rdef). In the classical derivation of dynamical friction (Chandrasekhar 1943) the Coulomb logarithm represents the ratio between the ‘minimum’ and ‘maximum’ impact parameters that affect the velocity change. The minimum impact parameter represents that a deflection of 90° is required, which in the Keplerian case is the deflection radius:
(5)
(6)
while the maximum impact parameter is the distance at which the stellar density becomes sufficiently ‘smaller’ than around the BH to become insignificant in modifying its velocity. In our case, gravity is computed self-consistently by the code outside |$\mathcal {S}$|⁠; therefore, the integration must be stopped at 4Δx if we do not want to double count dynamical friction. This naturally leads to lnΛ = ln(4Δx/rdef). Furthermore, as explained in Beckmann et al. (2017), using subgrid models when resolution is sufficient to account for dynamical friction can lead to incorrect results. For this reason, when 4Δx ≤ rdef, we set |$\vec{a}_\text{DF}$| to 0.

3 ADDITIONAL PHYSICS: GALAXIES AND BLACK HOLES

ramses follows the evolution of the gas using the second-order MUSCL-Hancock scheme for the Euler equations. The approximate Harten–Lax–Van Leer Contact (Toro 1997) Riemann solver with a MinMod total variation diminishing scheme to reconstruct the interpolated variables from their cell-centred values is used to compute the unsplit Godunov fluxes at cell interfaces. An equation of state of perfect gas composed of monoatomic particles with adiabatic index γ = 5/3 is assumed to close the full set of fluid equations. Collisionless particles (dark matter, stellar, and BH particles) are evolved using a particle-mesh solver with a cloud-in-cell interpolation. The size of the cloud-in-cell interpolation is that of the local cell for BHs and stars; however, dark matter particles can only project their mass on the grid down to a minimum cell size of ΔxDM > Δx (as these particles are usually larger in mass than stars or gas, we smooth their mass distribution to reduce their contribution to shot noise). When cloud-in-cell interpolation is used, therefore, even if the mass of the dark matter particle is larger than the BH mass, since the dark matter distribution is smoothed, scattering off dark matter particles becomes unimportant.

Gas is allowed to cool by hydrogen and helium with a contribution from metals using cooling curves from Sutherland & Dopita (1993) for temperatures above |$10^4 \, \mathrm{ K }$|⁠. For gas below |$10^4 \, \mathrm{ K }$| and down to our minimum temperature of |$10 \, \mathrm{ K }$|⁠, we use the fitting functions of Rosen & Bregman (1995).

Star formation is stochastically sampled from a random Poisson distribution (Rasera & Teyssier 2006): at each time-step Δt, in each cell of size Δx containing a gas mass Mgas, the mass of newly formed stars, M⋆, new, follows a Schmidt law:
(7)
where |$t_\mathrm{ff}=\sqrt{3\pi /32\mathrm{G}\rho _\mathrm{gas}}$| is the free-fall time, ρgas  =  Mgas/Δx3 is the gas density in the cell, and ϵ depends on the local turbulence of the gas, as detailed in Trebitsch et al. (2018).

For the feedback of supernovae, we use the mechanical feedback described in Kimm & Cen (2014), in which star particles older than 5 Myr release |$\eta _\mathrm{SN}\times 10^{50}\, \mathrm{ erg }\,\, \mathrm{ M}_{\odot }^{ -1}$|⁠, where ηSN = 0.2. The amount of energy and momentum deposited depends on local properties of the gas (density and metallicity) so that it captures either the Sedov or the supernovae-plough expansion phase of the explosion.

We use the model of BHs described in Dubois et al. (2012), where accretion is computed using the Bondi–Hoyle–Littleton formalism capped at the Eddington luminosity. AGN feedback consists of a dual-mode approach, where thermal energy, corresponding to 15 per cent of the bolometric luminosity (with radiative efficiency of ϵr = 0.1), is injected at high accretion rates (luminosity above 0.01 the Eddington luminosity); otherwise feedback is modelled with a bipolar jet with a velocity of |$10^4\, \rm km\, s^{-1}$| and an efficiency of 100 per cent. We slightly modify the implementation of BH dynamics: in the original ramses version, the mass of the BH is deposited on to the so-called ‘cloud’ particles, which uniformly pave a sphere of 4Δx radius on a grid of Δx/2 inter-cloud distance. This has the effect of smoothing the density, and therefore, when two BHs pass close by, their potential is shallower than it should be, and this delays the formation of the binary. We simply deposit all the mass of BHs on to their central cloud particle and then perform the cloud-in-cell to obtain more accurate dynamics, while using the rest of cloud particles to compute the Bondi–Hoyle–Littleton accretion rate on to the BH.

We include the dynamical friction implementation as described in Section 2 when necessary, and dynamical friction from gas was already included (Dubois, Volonteri & Silk 2014a) using equation (12) from Ostriker (1999). Two BHs are allowed to merge when they are separated by less than 4Δx and the kinetic energy of the binary is lower than the gravitational energy.

Due to our inability to resolve the cold and dense regions of the interstellar medium, gas dynamical friction, and gas accretion, which depend linearly on the gas density, ρgas can be underestimated in simulations. To correct for this lack of resolution, Booth & Schaye (2009); Dubois et al. (2014a) compute the expressions obtained analytically (Chandrasekhar gas dynamical friction and Bondi accretion) and boost them by:
(8)
where ρth is a free parameter, similar to that used for star formation, which is linked to the Jeans length and depends on resolution. This parameter is typically calibrated via the phase diagram of gas, and it is |${\sim } 1\, \, \mathrm{ amu }\, \mathrm{ cm }^{-3}$| for Δx = 100 pc and |${\sim } 50\, \, \mathrm{ amu }\, \mathrm{ cm }^{-3}$| for Δx = 10 pc. ξ can differ for accretion and dynamical friction, for the rest of the paper, we will use ξ = α when we relate to boosting accretion, and ξ = β when we refer to gas dynamical friction. Booth & Schaye (2009) performed a parameter study and found that α = 2 is the optimal value to recover the BH mass–galaxy mass relation (Kormendy & Ho 2013b). Although the expression of these boosts is not physically motivated, they have been proven to give excellent match with observations, as shown by large cosmological simulations (Dubois et al. 2014b). For this work we will either use ξ = 2 either ξ = 0 (no boost), and in Section 4.2 we also vary ρth within an order of magnitude of the typical value expected for the chosen resolution. We briefly explore the effects of the boost in Section 4.2.

4 VALIDATION OF THE DYNAMICAL FRICTION IMPLEMENTATION

4.1 Isolated dark matter halo

In order to compare the dynamical friction time-scale with analytical estimates (Lacey & Cole 1993; Colpi et al. 1999; Taffoni et al. 2003), we test our implementation following the dynamics of a BH moving in a dark matter halo.

The dark matter halo, initialized with dice (Perret 2016), follows a Navarro, Frank, and White (NFW; Navarro, Frenk & White 1997) profile with a total virial mass |$M_\text{vir}=2\times 10^{11}\, \mathrm{ M}_{\odot }$|⁠, a concentration parameter of 4 and a virial radius |$R_\text{vir}=45\, \mathrm{ kpc }$|⁠, typical of redshift 3. We set the total spin parameter to 0.04 consistent with the average spin parameter of cosmological dark matter haloes (Bullock et al. 2001), which only mildly evolves between z = 3 and today (Mu noz-Cuartas et al. 2011; Ahn et al. 2014). The BH mass is set to |$10^8 \, \mathrm{ M}_{\odot }$|⁠, it is initially |$5\, \mathrm{ kpc }$| away from the centre with a tangential velocity of |$57 \, \mathrm{ km }\, \mathrm{ s }^{-1}$|⁠, corresponding to 50 per cent of the circular velocity. In the simulation, the influence radius varies between 10 and 100 pc: it is at best resolved by two-cell elements, therefore dynamics is generally not treated properly and dynamical friction must be added ad hoc with our subgrid model when necessary.

As done in Tremmel et al. (2015), we can estimate the goodness of the method by comparing it to the analytical estimate of the ‘sinking time’, τDF, defined as the time it will take for a satellite to sink to a target, using equation (12) from Taffoni et al. (2003):
(9)
(10)
where Mvir is the virial mass of the target, vc is the circular velocity at the virial radius, G the gravitational constant, rc is the radius at which a test particle moving in the potential of the target has the same energy as the satellite, Ms is the mass of the satellite, J is the specific angular momentum of the satellite in the frame of the target, Jc is the specific angular momentum of circular orbit at rc and α depends on Ms, Mvir, Rvir, and rc and is given by equation (15) from Taffoni et al. (2003). In this case, the target is the halo and the satellite is the BH. Using this approach, we find that the BH should sink in the potential well of the halo in |$600 \, \mathrm{ Myr }$|⁠.

We perform two simulations which only differ by the presence (PD), or not (NoDrag), of dynamical friction on to the BH using our subgrid model. In both cases, the size of the box is |$100\, \mathrm{ kpc }$|⁠, slightly larger than 2Rvir and we allow refinement from levels 7–11, leading to a maximum physical resolution of |$\Delta x = \Delta x_{\rm DM} = 50 \, \mathrm{ pc }$|⁠, similar to what simulations reach in cosmological zooms (Dubois et al. 2014a). The refinement is done using a quasi-Lagrangian criterion: a cell is refined if its mass exceeds 8 × mDM, where mDM is the mass of dark matter particles, and we refine at maximum level up to 4Δx around the BH. We set the mass of dark matter particles to |$10^5 \, \mathrm{ M}_{\odot }$|⁠, in good agreement with the value suggested by Power et al. (2003): |$m_\text{DM}=M_\text{vir}\left(R_\text{vir}/\Delta x \right)^{-2}\sim 2\times 10^5 \, \mathrm{ M}_{\odot }$|⁠.

We show in Fig. 1, for both simulations, the distance between the BH and the centre of the halo; we also include for comparison the analytical estimate given by equation (10). The result is quite clear and in agreement with Tremmel et al. (2015): adding unresolved dynamical friction contributes to recover sinking times estimated analytically. In the following section, we set in a more realistic problem where a BH sinks in a galaxy including not only dark matter but also gas, stars, and many associated processes (cooling, star formation, supernovae feedback).

Distance of the BH to the centre of the halo as a function of time. The vertical dashed line is the analytical estimate of the total sinking time from Taffoni et al. (2003). With our prescription, the sinking time is much shorter than without and in very good agreement with the analytical estimate. The resolution in these simulations is 50 pc and the BH to dark matter mass particle ratio set equal to 1000. See section 4.1 for details.
Figure 1.

Distance of the BH to the centre of the halo as a function of time. The vertical dashed line is the analytical estimate of the total sinking time from Taffoni et al. (2003). With our prescription, the sinking time is much shorter than without and in very good agreement with the analytical estimate. The resolution in these simulations is 50 pc and the BH to dark matter mass particle ratio set equal to 1000. See section 4.1 for details.

Open in new tabDownload slide

4.2 Isolated galaxy

We run a suite of simulations (see Table 1) of a BH sinking in the potential well of an idealized isolated galaxy. Our suite contains low-resolution (⁠|$\Delta x = \Delta x _{\rm DM} = 50\, \mathrm{ pc }$|⁠) simulations, similar to what high-resolution zoomed cosmological simulations can reach today. Thus it is a good test to see how our implementation will act in this context. Contrary to the dark matter halo case, we do not have analytical estimates to provide a benchmark. To overcome this issue, we run a high-resolution test (⁠|$\Delta x = \Delta x _{\rm DM} = 1\, \mathrm{ pc }$|⁠) to perform the comparison. The set-up is chosen such that, with 50 pc resolution, during the sinking, the deflection radius, as defined in equation (5), is not always resolved (see Fig. 2). In this case, dynamics is not properly treated and dynamical friction must be added ad hoc with our subgrid model. Conversely, with 1 pc resolution the deflection radius is always resolved during the sinking and dynamical friction is well captured by the gravity solver of ramses, thus providing the correct dynamics.

Deflection radius (solid line) and resolution of the different simulations (dashed lines) of idealized isolated galaxies. In the low-resolution case the deflection radius is not always resolved, leading to incorrect dynamics of the BH and the need to add unresolved dynamical friction. In the high resolution run the deflection radius is always resolved and dynamical friction is self-consistently captured by the gravity solver. All quantities shown as a function of time.
Figure 2.

Deflection radius (solid line) and resolution of the different simulations (dashed lines) of idealized isolated galaxies. In the low-resolution case the deflection radius is not always resolved, leading to incorrect dynamics of the BH and the need to add unresolved dynamical friction. In the high resolution run the deflection radius is always resolved and dynamical friction is self-consistently captured by the gravity solver. All quantities shown as a function of time.

Open in new tabDownload slide
Table 1.
Open in new tab

Different simulations performed for the isolated galaxy test, with their name and the use, or not, of our new model. [1 and 2] Typical density and exponent that we used to boost gas friction, following equation (8). [3] Resolution of the simulation.

NameDynamical friction|$\rho _\text{th} ^{[1]}$|β[2]Δx[3]
amu cc−1pc
NoDragNoneXX50
GnB_nPDGas (no boost)X050
GB0.1_nPDGas (boost)0.1250
GB1_nPDGas (boost)1250
GB15_nPDGas (boost)15250
GnB_PDGas (no boost)+stars+dark matterX050
HRNoneXX0.76
NameDynamical friction|$\rho _\text{th} ^{[1]}$|β[2]Δx[3]
amu cc−1pc
NoDragNoneXX50
GnB_nPDGas (no boost)X050
GB0.1_nPDGas (boost)0.1250
GB1_nPDGas (boost)1250
GB15_nPDGas (boost)15250
GnB_PDGas (no boost)+stars+dark matterX050
HRNoneXX0.76
Table 1.
Open in new tab

Different simulations performed for the isolated galaxy test, with their name and the use, or not, of our new model. [1 and 2] Typical density and exponent that we used to boost gas friction, following equation (8). [3] Resolution of the simulation.

NameDynamical friction|$\rho _\text{th} ^{[1]}$|β[2]Δx[3]
amu cc−1pc
NoDragNoneXX50
GnB_nPDGas (no boost)X050
GB0.1_nPDGas (boost)0.1250
GB1_nPDGas (boost)1250
GB15_nPDGas (boost)15250
GnB_PDGas (no boost)+stars+dark matterX050
HRNoneXX0.76
NameDynamical friction|$\rho _\text{th} ^{[1]}$|β[2]Δx[3]
amu cc−1pc
NoDragNoneXX50
GnB_nPDGas (no boost)X050
GB0.1_nPDGas (boost)0.1250
GB1_nPDGas (boost)1250
GB15_nPDGas (boost)15250
GnB_PDGas (no boost)+stars+dark matterX050
HRNoneXX0.76

We initialize with DICE an ideal galaxy at redshift three with a total virial mass of |$2\times 10^{11}\, \mathrm{ M}_{\odot }$| and a spin parameter of 0.04. The galaxy is composed of four components.

  • A dark matter halo with a mass of |$1.95\times 10^{11} \, \mathrm{ M}_{\odot }$|⁠, slightly lighter than in Section 4.1. It has a virial radius of |$45\, \mathrm{ kpc }$| and the density follows an NFW profile with a concentration parameter of 4.

  • A gas disc with a total mass of |$2.4\times 10^9 \, \mathrm{ M}_{\odot }$|⁠. The density follows an exponential disc + sech-z profile with a scale radius of |$1.28\, \mathrm{ kpc }$| and an aspect ratio of 1:10. We impose an initial constant absolute metallicity and temperature of 10−3 and |$10^5\, \mathrm{ K }$|⁠, respectively.

  • A stellar disc with a total mass of |$1.6\times 10^9 \, \mathrm{ M}_{\odot }$|⁠. The density follows an exponential disc + sech-z profile with a scale radius of |$1.28\, \mathrm{ kpc }$| and an aspect ratio of 1:10. We impose an initial constant absolute metallicity of 10−3. Additionally, to avoid unphysical initial starbursts regularly found in ideal simulations (Capelo et al. 2015), we give an age distribution to stellar particles to mimic a 5 |$\, \mathrm{ M}_{\odot } \, \mathrm{ yr }^{-1}$| star formation rate.

  • A stellar bulge with a total mass of |$8\times 10^8 \, \mathrm{ M}_{\odot }$|⁠. The density follows a Hernquist profile (Hernquist 1990) with a scale radius of |$0.128\, \mathrm{ kpc }$|⁠. We impose a constant absolute metallicity of 2 × 10−4 (5 times smaller than in the disc to mimic the older age of stars in the bulge). Similarly, we give an age to stellar particles to mimic a 0.5 |$\, \mathrm{ M}_{\odot } \, \mathrm{ yr }^{-1}$| star formation rate.

In the low-resolution simulations (50 pc), the mass of dark matter particles is set to |$10^6~\, \mathrm{ M}_{\odot }$| and that of star particles to |$2\times 10^4 \, \mathrm{ M}_{\odot }$|⁠. In the high-resolution simulations (1 pc) the mass of dark matter particles is set to |$5\times 10^4 \, \mathrm{ M}_{\odot }$| and that of star particles to |$2\times 10^3 \, \mathrm{ M}_{\odot }$|⁠. In both cases the size of the box is 100 kpc and we allow for refinement from levels 7–11 in the low-resolution simulations and from 7 to 17 in the high resolution one, refining the mesh when |$M_\mathrm{DM}^\mathrm{cell}+ 10 M_\mathrm{ b}^\mathrm{cell} \ge 8 m_\mathrm{DM}$|⁠, where |$M_\mathrm{DM}^{\rm cell}$| and |$M_\mathrm{ b}^\mathrm{cell}$| are, respectively, the mass of dark matter and baryons in the cell. Maximum refinement is enforced within 4Δx around the BH.

After initializing this galaxy, we switch on cooling, star formation, supernovae feedback (see Section 3), and let the galaxy relax for 100 Myr. At that point, a BH with mass |$10^7 \, \mathrm{ M}_{\odot }$| is placed in the z = 0 plane, 1 kpc from the centre and with a tangential velocity of 21 km s−1, corresponding to 30 per cent of the circular velocity. Accretion and feedback from the BH are not included in order to keep the BH mass constant and isolate the effects of dynamical friction. We include dynamical friction with different implementations: from collisionless particles and gas without boost (i.e.  no free parameters), or only from gas, with or without a boost factor. The simulation properties and set-up are summarized in Table 1.

We show in Fig. 3, for all our simulations, the distance between the BH and the centre of the galaxy as a function of time. We first stress the difference between low-resolution simulations with gravity only, i.e. without including the dynamical friction model (NoDrag, blue line), and the simulations at high resolution where the deflection radius is resolved (HR, black line). In agreement with the results of Pfister et al. (2017), resolving at least the deflection radius is mandatory to properly capture the dynamics of the BH in the dynamical friction phase.

Distance of the BH to the centre of the galaxy as a function of time and resolutions of the different simulations (dashed lines). The good behaviour of our model is confirmed by the agreement between the solid orange curve (low resolution, use of our model) and the solid black one (high resolution). See section 4.2 for details.
Figure 3.

Distance of the BH to the centre of the galaxy as a function of time and resolutions of the different simulations (dashed lines). The good behaviour of our model is confirmed by the agreement between the solid orange curve (low resolution, use of our model) and the solid black one (high resolution). See section 4.2 for details.

Open in new tabDownload slide

We now compare simulations where we vary ρth (GB0.1_nPD, GB1_nPD, GB15_nPD) but we do not include dynamical friction from stars and dark matter. As expected, the lower ρth, the larger the boost, the faster the BH sinks. The choice of ρth must be performed accurately: if ρth is too low, BHs can get caught in a passing clump and either follow the clump outside the galaxy centre, or remain artificially in a dense environment where accretion is triggered, resulting in an overestimate of the mass of BHs. If ρth is too high, instead, the correction to dynamical friction is insufficient and the orbital decay is delayed. In this particular case, ρth between 0.1 and 1|$\, \mathrm{ amu }\, \mathrm{ cm }^{-3}$| is the best value to recover the high-resolution results, but the exact value may depend on additional factors such as the gas fraction (50 per cent in our case).

We finish with the simulation where the influence radius is not always resolved, but in which we include subgrid dynamical friction from stars, dark matter, and gas (without any boost) following our implementation (GnB_PD, orange line). This implementation does not contain any free parameters and avoids the arbitrary choice of ρth. This simulation is in excellent agreement with the high-resolution simulation (HR, thin black line), confirming the good behaviour of our model in a realistic, although idealized, galaxy.

4.3 Limits of the model: low-mass black holes

In this section, we explore the limits of our implementation when a BH has a mass so low that two-body interactions with star and dark matter particles significantly perturb its dynamics.

We run simulations similar to those described in Section 4.1 but decreasing the mass of the BH down to the mass of dark matter particles (⁠|$m_{\rm DM}=10^5\, \mathrm{ M}_{\odot }$|⁠). To contain computational costs, we also change the orbital parameters of the BH such as the analytical estimates from Taffoni et al. (2003); however, τDF, remains a few Gyrs. We list the parameters of the simulations in Table 2.

Table 2.
Open in new tab

Different simulations we perform to test the limits of our model in terms of particle mass ratio. We indicate the different mass ratio between the BH and dark matter particles, the initial distance of the BH from the centre of the halo, the initial velocity of the BH, and the analytical estimate for the time the BH should take to reach the centre of the halo from Taffoni et al. (2003). In all cases, we run a simulation with (PD) and without (NoDrag) our model.

M•/mDMd0v0τDF
|$\, \mathrm{ kpc }$||$\, \mathrm{ km }\, \mathrm{ s }^{-1}$||$\, \mathrm{ Gyr }$|
1028.22.38
1165.26
M•/mDMd0v0τDF
|$\, \mathrm{ kpc }$||$\, \mathrm{ km }\, \mathrm{ s }^{-1}$||$\, \mathrm{ Gyr }$|
1028.22.38
1165.26
Table 2.
Open in new tab

Different simulations we perform to test the limits of our model in terms of particle mass ratio. We indicate the different mass ratio between the BH and dark matter particles, the initial distance of the BH from the centre of the halo, the initial velocity of the BH, and the analytical estimate for the time the BH should take to reach the centre of the halo from Taffoni et al. (2003). In all cases, we run a simulation with (PD) and without (NoDrag) our model.

M•/mDMd0v0τDF
|$\, \mathrm{ kpc }$||$\, \mathrm{ km }\, \mathrm{ s }^{-1}$||$\, \mathrm{ Gyr }$|
1028.22.38
1165.26
M•/mDMd0v0τDF
|$\, \mathrm{ kpc }$||$\, \mathrm{ km }\, \mathrm{ s }^{-1}$||$\, \mathrm{ Gyr }$|
1028.22.38
1165.26

We show in Fig. 4 the distance of the BH to the centre of the halo as a function of time. It is clear that our model works very well when BHs have a mass larger than 10 times the mass of particles causing dynamical friction. If the mass of the BH is similar to that of particles causing dynamical friction; however, it is scattered through two-body interactions and the model becomes less reliable, as also noted by Tremmel et al. (2015).

Different simulations performed to test the effect of reducing the mass ratio between the BH and dark matter particles. We indicate the use (PD) or not (NoDrag) of our prescription for dynamical friction. If the BH mass is similar to the dark matter particle mass, the efficacy of the model becomes limited. See section 4.3 for details.
Figure 4.

Different simulations performed to test the effect of reducing the mass ratio between the BH and dark matter particles. We indicate the use (PD) or not (NoDrag) of our prescription for dynamical friction. If the BH mass is similar to the dark matter particle mass, the efficacy of the model becomes limited. See section 4.3 for details.

Open in new tabDownload slide

In Sections 4.2 and 5.1, the mass of dark matter particles is larger than that of BHs. However, we use cloud-in-cell interpolation to smooth the dark matter distribution, and we ensure that the mass of star particles, which are the main source of dynamical friction, is lower than the mass of BHs.

5 COSMOLOGICAL SIMULATIONS

5.1 Set-up

We now move to the full cosmological context, endeavouring to study the dynamical behaviour of seed BHs in high-redshift galaxies. We run a suite of cosmological simulations, with the code ramses. We zoom-in on one halo using different prescriptions for the dynamics of BHs. As we are interested in understanding the evolution of BHs in typical galaxies, we chose a halo with a minor/major merger rate comparable to the mean evolution obtained by Fakhouri, Ma & Boylan-Kolchin (2010) in this mass range. The physics is similar to that of the simulations described in Section 3 but for the refinement strategy: we refine if |$M_\mathrm{DM}^\mathrm{cell}+ (\Omega _m / \Omega _\mathrm{ b} -1) M_\mathrm{ b}^\mathrm{cell} \ge 8 m_\mathrm{DM}$|⁠, where MDM and |$M_\mathrm{ b}^\mathrm{cell}$| are, respectively, the mass of dark matter and baryons in the cell, and Ωm and Ωb are the total matter and baryon density. The minimum cell size, Δx is kept roughly constant in proper physical size with redshift: an additional level of refinement is added every time the expansion factor, aexp, decreases by a factor of two, such that the maximum level, lmax, is reached at aexp = 0.8. For simplicity, we further assume that |$\Delta x=L_\textrm{box}/2^{l_\mathrm{max}}$|⁠, where Lbox is the size of the box at redshift 0. Concerning the subgrid physics of BHs (see Section 3) we use β = 2 to boost accretion, gas friction is not boosted (α = 0), and the value of ρth depends on resolution. The specifications of each simulations are described in Table 3.

Table 3.
Open in new tab

Properties of the suite of cosmological simulations performed. [1] Typical density and exponent that we used to boost accretion, following equation (8). [2–5] Spatial/mass resolution of the simulation. [6] Seed mass of BHs.

NameParticle|$\rho _\text{th} ^{[1]}$|Δx[2]|$\Delta x^{[3]}_{\rm DM}$||$m_\star ^{[4]}$||$m_\text{DM} ^{[5]}$||$M_\bullet ^{[6]}$|
dynamical frictionamu cc−1pcpc|$\, \mathrm{ M}_{\odot }$||$\, \mathrm{ M}_{\odot }$||$\, \mathrm{ M}_{\odot }$|
LR_PD_BH1e4✓57223002 × 1042 × 106104
MR_PD_BH1e4✓103623002 × 1042 × 106104
HR_PD_BH1e4✓5095722 × 1032 × 105104
HR_nPD_BH1e4✗5095722 × 1032 × 105104
HR_PD_BH1e5✓5095722 × 1032 × 105105
HR_nPD_BH1e5✗5095722 × 1032 × 105105
NameParticle|$\rho _\text{th} ^{[1]}$|Δx[2]|$\Delta x^{[3]}_{\rm DM}$||$m_\star ^{[4]}$||$m_\text{DM} ^{[5]}$||$M_\bullet ^{[6]}$|
dynamical frictionamu cc−1pcpc|$\, \mathrm{ M}_{\odot }$||$\, \mathrm{ M}_{\odot }$||$\, \mathrm{ M}_{\odot }$|
LR_PD_BH1e4✓57223002 × 1042 × 106104
MR_PD_BH1e4✓103623002 × 1042 × 106104
HR_PD_BH1e4✓5095722 × 1032 × 105104
HR_nPD_BH1e4✗5095722 × 1032 × 105104
HR_PD_BH1e5✓5095722 × 1032 × 105105
HR_nPD_BH1e5✗5095722 × 1032 × 105105
Table 3.
Open in new tab

Properties of the suite of cosmological simulations performed. [1] Typical density and exponent that we used to boost accretion, following equation (8). [2–5] Spatial/mass resolution of the simulation. [6] Seed mass of BHs.

NameParticle|$\rho _\text{th} ^{[1]}$|Δx[2]|$\Delta x^{[3]}_{\rm DM}$||$m_\star ^{[4]}$||$m_\text{DM} ^{[5]}$||$M_\bullet ^{[6]}$|
dynamical frictionamu cc−1pcpc|$\, \mathrm{ M}_{\odot }$||$\, \mathrm{ M}_{\odot }$||$\, \mathrm{ M}_{\odot }$|
LR_PD_BH1e4✓57223002 × 1042 × 106104
MR_PD_BH1e4✓103623002 × 1042 × 106104
HR_PD_BH1e4✓5095722 × 1032 × 105104
HR_nPD_BH1e4✗5095722 × 1032 × 105104
HR_PD_BH1e5✓5095722 × 1032 × 105105
HR_nPD_BH1e5✗5095722 × 1032 × 105105
NameParticle|$\rho _\text{th} ^{[1]}$|Δx[2]|$\Delta x^{[3]}_{\rm DM}$||$m_\star ^{[4]}$||$m_\text{DM} ^{[5]}$||$M_\bullet ^{[6]}$|
dynamical frictionamu cc−1pcpc|$\, \mathrm{ M}_{\odot }$||$\, \mathrm{ M}_{\odot }$||$\, \mathrm{ M}_{\odot }$|
LR_PD_BH1e4✓57223002 × 1042 × 106104
MR_PD_BH1e4✓103623002 × 1042 × 106104
HR_PD_BH1e4✓5095722 × 1032 × 105104
HR_nPD_BH1e4✗5095722 × 1032 × 105104
HR_PD_BH1e5✓5095722 × 1032 × 105105
HR_nPD_BH1e5✗5095722 × 1032 × 105105

5.1.1 Initial conditions

We assume a ΛCDM cosmology with total matter density Ωm = 0.3089, baryon density Ωb = 0.0486, dark energy density |$\Omega _\Lambda = 0.6911$|⁠, amplitude of the matter power spectrum σ8 = 0.8159, ns = 0.9667 spectral index, and Hubble constant |$H_0 = 67.74 \, \mathrm{ km }\, \mathrm{ s }^{-1} \, \mathrm{ Mpc }^{-1}$| consistent with the Planck data (Planck Collaboration XIII, 2016). The initial conditions are produced with music (Hahn & Abel 2013). The box size of the simulations is |$L_\mathrm{box} = 73.8 \, \mathrm{ Mpc }$|⁠, with a coarse grid of 2563 dark matter particles corresponding to a dark matter mass resolution of |$m_\mathrm{DM,coarse} = 3 \times 10^9 \, \mathrm{ M}_{\odot }$|⁠. A high-resolution region is defined around a halo of |$M_\mathrm{vir} = 10^{12} \, \mathrm{ M}_{\odot }$| at z = 2 that contains only high-resolution dark matter particles (see Table 2 for the mass of high-resolution dark matter particles in each simulation) within 2 rvir (⁠|$r_\mathrm{vir} = 100 \, \mathrm{ kpc }$|⁠). The halo is a progenitor of a group of galaxies whose mass is |$M_\mathrm{vir} = 7\times 10^{12} \, \mathrm{ M}_{\odot }$| at z = 0.

5.1.2 Finding haloes and galaxies

We construct catalogues of haloes and galaxies using the adaptahop halo finder (Aubert, Pichon & Colombi 2004), which uses an SPH-like kernel to compute densities at the location of each particle and partitions the ensemble of particles into subhaloes based on saddle points in the density field. Haloes contain at least 200 dark matter particles. Galaxies are identified in the same way, and contain at least 200 stellar particles. We then construct a merger tree for haloes and galaxies with treemaker (Tweed et al. 2009).

5.1.3 Estimate of the sinking time

We consider once again the ‘sinking time’, τDF, defined as the time it takes for a satellite to sink to a target, using equation (10). To compute τDF for a BH in its own galaxy (Section 5.2), we consider that the satellite is the BH, for which we have the dynamical properties, and we consider that the target is the galaxy, for which we compute the different properties with the halo finder.

To compute τDF for a BH during a galaxy merger (Section 5.3), we have to take into account that Ms evolves. Initially, the BH is surrounded by its own galaxy, which is itself surrounded by a halo, and it is the system BH + galaxy + halo that undergoes dynamical friction. Therefore, we must match BHs to galaxies and galaxies to haloes to have the corrected satellite mass, i.e. Ms is similar to the mass of the halo.

In a second phase, the dark matter halo and outer stellar layers of the secondary galaxy disperse into that of the primary, and the BH remains surrounded only by a fraction of the initial stellar mass, and we identify the evolving Ms via the halo finder. Finally, the BH remains naked, and Ms is the BH mass. To give an order of magnitude for this final phase, in the early universe, where galaxies have velocity dispersion as small as tens of |$\, \mathrm{ km }\, \mathrm{ s }^{-1}$|⁠, unless the BH is very massive (⁠|${\gtrsim } 10^5 \, \mathrm{ M}_{\odot }$|⁠), or surrounded by a bound dense stellar cluster, acting as if |$\, M_\mathrm{ s}\,$| is larger, the sinking time is longer than Gyr if the distance to the centre if larger than |${\sim }100 \, \mathrm{ pc }$|⁠, which is likely to be the case if the BH is scattered due to anisotropies of the galaxy, either when it is in isolation or during mergers.

5.2 Dynamics of a seed black hole in its own galaxy

We focus on a satellite galaxy which merges with the main galaxy when the age of the Universe is about 1 Gyr. In Fig. 5, we show snapshots at the beginning of the interaction between the main galaxy, on which the figure is centred, and the satellite, to the top left of the main galaxy. This satellite hosts a BH and we study its dynamics while the galaxy is in relative isolation. This case is interesting because it explores the prospects for a seed BH to remain surrounded by dense cold gas available for growing the BH and make it observable as a faint AGN.

Top left: stellar density (black: $10^{-4} \, \mathrm{ M}_{\odot } \, \mathrm{ pc }^{-3}$, white: $1 \, \mathrm{ M}_{\odot } \, \mathrm{ pc }^{-3}$), centred on the main galaxy with an 80 comoving kpc box size, of LR_PD_BH1e4. Top right: MR_PD_BH1e4. Bottom left: HR_PD_BH1e4. Bottom right: HR_PD_BH1e5. The panels show the exact same galaxy at the same time to highlight the effects of resolution. In Section 5.2, we discuss the dynamics of a BH in the satellite galaxy on the top-left corner of each panel (the BH is highlighted in red and its ID is 2 in the four panels). In Section 5.3, we discuss instead the interaction between this BH and the main BH in the central galaxy (also highlighted in red, and its ID is 1 in the four panels).
Figure 5.

Top left: stellar density (black: |$10^{-4} \, \mathrm{ M}_{\odot } \, \mathrm{ pc }^{-3}$|⁠, white: |$1 \, \mathrm{ M}_{\odot } \, \mathrm{ pc }^{-3}$|⁠), centred on the main galaxy with an 80 comoving kpc box size, of LR_PD_BH1e4. Top right: MR_PD_BH1e4. Bottom left: HR_PD_BH1e4. Bottom right: HR_PD_BH1e5. The panels show the exact same galaxy at the same time to highlight the effects of resolution. In Section 5.2, we discuss the dynamics of a BH in the satellite galaxy on the top-left corner of each panel (the BH is highlighted in red and its ID is 2 in the four panels). In Section 5.3, we discuss instead the interaction between this BH and the main BH in the central galaxy (also highlighted in red, and its ID is 1 in the four panels).

Open in new tabDownload slide

We start by studying how the different sources of friction (dark matter, stars, and gas) contribute to the dynamical evolution. Fig. 6 presents the density in gas and stars around the BH (we do not include dark matter since its contribution is negligible). Gas is more chaotic than stars, but stars themselves do not provide a constant acceleration because they are also irregularly distributed. Beyond the sheer inhomogeneity, gas can shock, cool, inflow, and outflow making its dynamical friction contribution unpredictable a priori. The presence of satellites also perturbs the BH orbit when it is far from the centre, see e.g. Fig. 5: in a typical high-redshift environment a BH feels acceleration coming from different directions.

Top panel: mean stellar density within 4Δx around the BH in the satellite galaxy, as a function of time for a subset of the simulations listed in Table 3, as noted in the inset. Bottom panel: ratio of the stellar density and gas density within 4Δx around the satellite BH, as a function of time.
Figure 6.

Top panel: mean stellar density within 4Δx around the BH in the satellite galaxy, as a function of time for a subset of the simulations listed in Table 3, as noted in the inset. Bottom panel: ratio of the stellar density and gas density within 4Δx around the satellite BH, as a function of time.

Open in new tabDownload slide

Moving to how this affects the BH’s orbits, we show τDF as a function of time in the bottom panel of Fig. 7, computed for the different simulations, using the method described in Section 5.1.3. We also show in the top panel of Fig. 7 the distance of the BH to the centre of its host galaxy.

Top panel: distance of the BH from the centre of its host galaxy, before interaction with another larger galaxy, as a function of time for all simulations listed in Table 3. Bottom panel: sinking time τDF for the secondary BH with respect to its host galaxy, computed using equation (10) (replacing Ms by the mass of the BH in this equation), as a function of time.
Figure 7.

Top panel: distance of the BH from the centre of its host galaxy, before interaction with another larger galaxy, as a function of time for all simulations listed in Table 3. Bottom panel: sinking time τDF for the secondary BH with respect to its host galaxy, computed using equation (10) (replacing Ms by the mass of the BH in this equation), as a function of time.

Open in new tabDownload slide

First, we find that, as long as the seeding mass of the BH is |$10^4 \, \mathrm{ M}_{\odot }$|⁠, all the simulations, independently of the resolution and the different models used for the BH dynamics, show a similar trend: the sinking time is, at least, 1–10 Gyr. Since in all cases vc slowly increases from 7 to 30 |$\, \mathrm{ km }\, \mathrm{ s }^{-1}$| and M• remains close to the BH seed mass, the reason of this large τDF is the dependency of the sinking time with the distance of the centre of the galaxy, which is shown in the top panel of Fig. 7. Even in HR_PD_BH1e4 and HR_nPD_BH1e4, where the BH is five times heavier than the star particles (those mostly contributing to the dynamical friction here) and forms at |${\sim } 70 \, \mathrm{ pc }$| from the centre, it is rapidly ejected and remains hundreds of pc away from the centre. Clumps and anisotropies are observed both in the stellar and gas central distributions. Due to such irregularities in the underlying galaxy, the BH undergoes a physically motivated random walk out of the centre of the potential well, as it also happens in lower redshift dwarfs (Bellovary et al. 2019). When the BH is more massive, |$10^5 \, \mathrm{ M}_{\odot }$|⁠, it remains in the centre of its host, with a sinking time less than 100 Myr. |$10^5 \, \mathrm{ M}_{\odot }$| seems therefore to be the minimum requirement to imagine that a BH is well stabilized in the centre of its host. BHs with masses lower than |$10^5 \, \mathrm{ M}_{\odot }$| are scattered within the galaxy due to irregularities of the gas/stellar potential and oscillate around the centre of their host galaxies, remaining far from the dense gas regions, therefore we expect them to have low accretion rates (Smith et al. 2018) and be difficult to observe.

5.3 Formation of a black hole binary in a high-redshift galaxy merger

We now focus on the same satellite galaxy, and follow the dynamical evolution of its BH during and after its host infalls into the halo of the larger galaxy. It is typically after this kind of event, when the galaxy remnant has settled and the massive BHs have sunk to the centre of the potential well, which massive BH binaries form.

We show in Fig. 8 τDF as a function of time, for all the simulations (bottom panel). We see that, initially, when the BH is still embedded in the satellite galaxy (solid line), its dynamics is the same for all simulations: the large-scale dynamics is independent of the subgrid model we use. However, what happens following the disruption of the satellite galaxy (dotted line) differs significantly from one simulation to the other: in some cases, the satellite BH sinks towards the centre and ‘merges’ (we recall that BHs are allowed to merge when they are separated by less than 4Δx and the kinetic energy of the binary is lower than the gravitational energy, but the real merger happens below our resolution) with the central BH of the main galaxy (the subsequent evolution is shown as a dashed line), in other cases, the BH stalls hundreds of pc away from the centre. We also show in the top panel of Fig. 8 the distance of the satellite BH to the central galaxy it is sinking in.

Top panel: distance of the BH originally in the satellite galaxy from the centre of the main galaxy as a function of time for a subset of the simulations listed in Table 3, as noted in the inset. If the BH in the satellite galaxy merges with the BH of the central galaxy, we show its subsequent evolution with dashed lines. Bottom panel: sinking time τDF, computed using equation (10), as a function of time. We show the different phases: when the BH is still surrounded by material (solid line), when the BH is naked (dotted line; note the rapid increase in the sinking time because of the drop in Ms in equation 10), and when the BH has merged with the BH of the central galaxy (dashed line).
Figure 8.

Top panel: distance of the BH originally in the satellite galaxy from the centre of the main galaxy as a function of time for a subset of the simulations listed in Table 3, as noted in the inset. If the BH in the satellite galaxy merges with the BH of the central galaxy, we show its subsequent evolution with dashed lines. Bottom panel: sinking time τDF, computed using equation (10), as a function of time. We show the different phases: when the BH is still surrounded by material (solid line), when the BH is naked (dotted line; note the rapid increase in the sinking time because of the drop in Ms in equation 10), and when the BH has merged with the BH of the central galaxy (dashed line).

Open in new tabDownload slide

We first compare the simulations HR_PD_BH1e5 - HR_nPD_BH1e5, and HR_PD_BH1e4 - HR_nPD_BH1e4, which differ only by the use or not of our subgrid model for dynamical friction from stars and dark matter. Fig. 7 shows that the model does not help in keeping BHs in the centre, as discussed in Section 5.2: the galaxy is so chaotic that BHs wander no matter the implementation. When the galaxy is more settled; however, as it is the case when the satellite BH falls into the main galaxy, we see the effects of our model (see Fig. 8). When our prescription is used, the BH remains closer to the centre; none the less the BHs do not merge as would happen if the BHs were artificially repositioned at the centre of mass of the halo, as is sometimes done in cosmological simulations (e.g. Vogelsberger et al. 2013; Schaye et al. 2015).

We now focus on simulations with |$10^4 \, \mathrm{ M}_{\odot }$| seeds. After its stellar and gaseous envelope has been dispersed (dotted line), the BH should take 1–100 Gyr to sink towards the centre of the galaxy, and indeed, it stalls at ∼hundreds of pc. This is in agreement with our understanding of dynamical friction: it is a very long process if the mass of the BH is low. The presence of a nuclear star cluster could speed-up the process (Biernacki, Teyssier & Bleuler 2017), increasing the mass experiencing dynamical friction, but due to our limited resolution, such compact structures of typical size of a few pc to ∼ten of pc are not captured here (Georgiev et al. 2016), and the envelope of the BH is rapidly stripped (dotted line). In the medium-resolution case (MR_PD_BH1e4) the BH in the larger galaxy has also been scattered, similarly to what happened for the case studied in Section 5.2. Accidentally, the two BHs merge while they are both off centre and the remnant of this merger remains hundreds of pc away from the centre. If we admit that this merger is physical, it is interesting to note that mergers of light seeds BHs are possible, though the dynamics is highly erratic. Multiple BHs in galaxies, each inherited from a different merger, are generically expected (e.g. Governato, Colpi & Maraschi 1994; Schneider et al. 2002; Volonteri & Perna 2005; Bonetti et al. 2018; Tremmel et al. 2018b).

Finally, we compare HR_PD_BH1e5 and HR_PD_BH1e4 which differ only by the seed mass of the BH. In HR_PD_BH1e5, the BH being more massive, it remains surrounded by a dense stellar concentration which does not disrupt (no dotted line), increasing even more the effective |$\, M_\mathrm{ s}\,$| and resulting in a smooth decay to the centre of the main galaxy and a BH merger.

These experiments make us believe that |$\lt 10^4\, \mathrm{ M}_{\odot }$| seed BHs are less likely to contribute to the merging population observable by LISA than larger mass seed BHs. This does not exclude that these low-mass BHs may eventually sink in the centre of galaxies and contribute to the massive BH population, but the presence of a dense stellar cluster or of bound gas on scales not resolved in this study, which would make the effective |$\, M_\mathrm{ s}\,$| larger, appears to be crucial (e.g. Callegari et al. 2009). Off-centre mergers, happening by chance, as in MR_nGB_PD_BH1e4 could also contribute.

6 CONCLUSIONS

We present a model to correct the dynamics of BHs in the ramses code, which is currently the only code to include a physically motivated model for dynamical friction on to BHs from gas, stars, and dark matter. We use this model in a suite of cosmological simulations to understand the dynamics of seed BHs (⁠|$10^4\!-\!10^5 \, \mathrm{ M}_{\odot }$|⁠) in high-redshift galaxies and during galaxy mergers. We summarize our findings below:

  • dynamical friction from stars has generically a more stabilizing effect than dynamical friction from gas, which can shock and is subject to inflows and outflows. In high-redshift galaxies, however, the stellar distribution is irregular and does not necessarily provide a smooth distribution within which BHs can decay undisturbed. The presence of satellite galaxies can also perturb the orbit of a BH.

  • From the results of our best resolution cosmological simulation, BHs with masses of the order of |$10^4 \, \mathrm{ M}_{\odot }$| are subject to the fluctuations of the underlying stellar gravitational potential, which leads to a random walk-type of trajectory. This appears to be unique of a high-z environment in which substructures undergo rapid evolution. If BHs were to be seeded in nuclear star clusters, or had masses of |$10^5 \, \mathrm{ M}_{\odot }$| or higher, they would be well stabilized galaxy centres.

  • Similarly, following galaxy mergers, if the mass of BHs in satellite galaxies is |$\sim 10^4 \, \mathrm{ M}_{\odot }$|⁠, it is unlikely that they participate in the merging population, although off-centre mergers can occur fortuitously. If seed BHs have larger masses, |${\sim } 10^5 \, \mathrm{ M}_{\odot }$|⁠, or they are embedded in dense bound stellar or gaseous envelopes, they can smoothly reach the centre of the larger galaxy and merge with the companion BH.

ACKNOWLEDGEMENTS

MV and HP acknowledge support from the European Research Council (Project no. 614199, ‘BLACK’). HP acknowledges support from the COST Association (CA16104 Gravitational waves, black holes, and fundamental physics), and thanks the University of Milano Bicocca for hosting him. This work was granted access to the HPC resources of under the allocations A0020406955 and A0040406955 made by GENCI. This work has made use of the Horizon Cluster hosted by the Institut d’Astrophysique de Paris; we thank Stephane Rouberol for running smoothly this cluster for us. Finally, we acknowledge Darren Croton for his useful comments.

REFERENCES

Ahn
J.
,
Kim
J.
,
Shin
J.
,
Kim
S.
,
Choi
Y.-Y.
,
2014
,
J. Korean Astron. Soc.
,
47
:
10.5303/JKAS.2014.47.2.77

Crossref

 

maro-Seoane
P.
et al. .,
2017
, preprint (arXiv:1702.00786)

Antonini
F.
,
Merritt
D.
,
2012
,
ApJ
,
745
,
83
10.1088/0004-637X/745/1/83

Crossref
Search ADS
 

Aubert
D.
,
Pichon
C.
,
Colombi
S.
,
2004
,
MNRAS
,
352
,
376
10.1111/j.1365-2966.2004.07883.x

Crossref
Search ADS
 

Barausse
E.
,
2012
,
MNRAS
,
423
,
2533
10.1111/j.1365-2966.2012.21057.x

Crossref
Search ADS
 

Beckmann
R. S.
et al. .,
2017
,
MNRAS
,
472
,
949
10.1093/mnras/stx1831

Crossref
Search ADS
 

Begelman
M. C.
,
Blandford
R. D.
,
Rees
M. J.
,
1980
,
Nature
,
287
,
307
10.1038/287307a0

Crossref
Search ADS
 

Bellovary
J. M.
,
Cleary
C. E.
,
Munshi
F.
,
Tremmel
M.
,
Christensen
C. R.
,
Brooks
A.
,
Quinn
T. R.
,
2019
,
MNRAS
,
482
,
2913
10.1093/mnras/sty2842

Crossref

 

Biernacki
P.
,
Teyssier
R.
,
Bleuler
A.
,
2017
,
MNRAS
,
469
,
295
10.1093/mnras/stx845

Crossref
Search ADS
 

Bonetti
M.
,
Haardt
F.
,
Sesana
A.
,
Barausse
E.
,
2018
,
MNRAS
,
477
,
3910
10.1093/mnras/sty896

Crossref
Search ADS
 

Booth
C. M.
,
Schaye
J.
,
2009
,
MNRAS
,
398
,
53
10.1111/j.1365-2966.2009.15043.x

Crossref
Search ADS
 

Bullock
J. S.
,
Dekel
A.
,
Kolatt
T. S.
,
Kravtsov
A. V.
,
Klypin
A. A.
,
Porciani
C.
,
Primack
J. R.
,
2001
,
ApJ
,
555
,
240
10.1086/321477

Crossref
Search ADS
 

Callegari
S.
,
Mayer
L.
,
Kazantzidis
S.
,
Colpi
M.
,
Governato
F.
,
Quinn
T.
,
Wadsley
J.
,
2009
,
ApJ
,
696
,
L89
10.1088/0004-637X/696/1/L89

Crossref
Search ADS
 

Capelo
P. R.
,
Volonteri
M.
,
Dotti
M.
,
Bellovary
J. M.
,
Mayer
L.
,
Governato
F.
,
2015
,
MNRAS
,
447
,
2123

Crossref
Search ADS
 

Chandrasekhar
S.
,
1943
,
ApJ
,
97
,
255
10.1086/144517

Crossref
Search ADS
 

Chapon
D.
,
Mayer
L.
,
Teyssier
R.
,
2013
,
MNRAS
,
429
,
3114
10.1093/mnras/sts568

Crossref
Search ADS
 

Colpi
M.
,
Mayer
L.
,
Governato
F.
,
1999
,
ApJ
,
525
,
720

Crossref
Search ADS
 

Dayal
P.
,
Rossi
E. M.
,
Shiralilou
B.
,
Piana
O.
,
Choudhury
T. R.
,
Volonteri
M.
,
2018
, preprint (arXiv e-prints)

Dosopoulou
F.
,
Antonini
F.
,
2017
,
ApJ
,
840
,
31
10.3847/1538-4357/aa6b58

Crossref
Search ADS
 

Dotti
M.
,
Colpi
M.
,
Haardt
F.
,
Mayer
L.
,
2007
,
MNRAS
,
379
,
956
10.1111/j.1365-2966.2007.12010.x

Crossref
Search ADS
 

Dubois
Y.
et al. .,
2014b
,
MNRAS
,
444
,
1453
10.1093/mnras/stu1227

Crossref
Search ADS
 

Dubois
Y.
,
Devriendt
J.
,
Slyz
A.
,
Teyssier
R.
,
2012
,
MNRAS
,
420
,
2662
10.1111/j.1365-2966.2011.20236.x

Crossref
Search ADS
 

Dubois
Y.
,
Volonteri
M.
,
Silk
J.
,
2014a
,
MNRAS
,
440
,
1590
10.1093/mnras/stu373

Crossref
Search ADS
 

Fakhouri
O.
,
Ma
C.-P.
,
Boylan-Kolchin
M.
,
2010
,
MNRAS
,
406
,
2267
10.1111/j.1365-2966.2010.16859.x

Crossref
Search ADS
 

Georgiev
I. Y.
,
Böker
T.
,
Leigh
N.
,
Lützgendorf
N.
,
Neumayer
N.
,
2016
,
MNRAS
,
457
,
2122
10.1093/mnras/stw093

Crossref
Search ADS
 

Goicovic
F. G.
,
Cuadra
J.
,
Sesana
A.
,
Stasyszyn
F.
,
Amaro-Seoane
P.
,
Tanaka
T. L.
,
2016
,
MNRAS
,
455
,
1989

Crossref
Search ADS
 

Governato
F.
,
Colpi
M.
,
Maraschi
L.
,
1994
,
MNRAS
,
271
,
317

Crossref
Search ADS
 

Hahn
O.
,
Abel
T.
,
2013
,
Astrophysics Source Code Library
,
record ascl:1311.011

Google Scholar

OpenURL Placeholder Text

 

Haiman
Z.
,
Kocsis
B.
,
Menou
K.
,
2009
,
ApJ
,
700
,
1952
10.1088/0004-637X/700/2/1952

Crossref
Search ADS
 

Hartwig
T.
,
Agarwal
B.
,
Regan
J. A.
,
2018
,
MNRAS
,
479
,
L23
10.1093/mnrasl/sly091

Crossref
Search ADS
 

Hernquist
L.
,
1990
,
ApJ
,
356
,
359
10.1086/168845

Crossref
Search ADS
 

Khan
F. M.
,
Berentzen
I.
,
Berczik
P.
,
Just
A.
,
Mayer
L.
,
Nitadori
K.
,
Callegari
S.
,
2012
,
ApJ
,
756
,
30
10.1088/0004-637X/756/1/30

Crossref
Search ADS
 

Kimm
T.
,
Cen
R.
,
2014
,
ApJ
,
788
,
121
10.1088/0004-637X/788/2/121

Crossref
Search ADS
 

Kormendy
J.
,
Ho
L. C.
,
2013a
,
ARA&A
,
51
,
511
10.1146/annurev-astro-082708-101811

Crossref
Search ADS
 

Kormendy
J.
,
Ho
L. C.
,
2013b
,
ARAA
,
51
,
511
10.1146/annurev-astro-082708-101811

Crossref
Search ADS
 

Lacey
C.
,
Cole
S.
,
1993
,
MNRAS
,
262
,
627
10.1093/mnras/262.3.627

Crossref
Search ADS
 

Mayer
L.
,
Kazantzidis
S.
,
Madau
P.
,
Colpi
M.
,
Quinn
T.
,
Wadsley
J.
,
2007
,
Science
,
316
,
1874
10.1126/science.1141858

Crossref
Search ADS
PubMed
 

Mu noz-Cuartas
J. C.
,
Macciò
A. V.
,
Gottlöber
S.
,
Dutton
A. A.
,
2011
,
MNRAS
,
411
,
584
10.1111/j.1365-2966.2010.17704.x

Crossref
Search ADS
 

Navarro
J. F.
,
Frenk
C. S.
,
White
S. D. M.
,
1997
,
ApJ
,
490
,
493
10.1086/304888

Crossref
Search ADS
 

Ostriker
E. C.
,
1999
,
ApJ
,
513
,
252
10.1086/306858

Crossref
Search ADS
 

Perret
V.
,
2016
,
Astrophysics Source Code Library
,
record ascl:1607.002

Google Scholar

OpenURL Placeholder Text

 

Pfister
H.
,
Lupi
A.
,
Capelo
P. R.
,
Volonteri
M.
,
Bellovary
J. M.
,
Dotti
M.
,
2017
,
MNRAS
,
471
,
3646
10.1093/mnras/stx1853

Crossref
Search ADS
 

Planck Collaboration XIII
,
2016
,
A&A
,
594
,
A13
10.1051/0004-6361/201525830

Crossref
Search ADS
 

Power
C.
,
Navarro
J. F.
,
Jenkins
A.
,
Frenk
C. S.
,
White
S. D. M.
,
Springel
V.
,
Stadel
J.
,
Quinn
T.
,
2003
,
MNRAS
,
338
,
14
10.1046/j.1365-8711.2003.05925.x

Crossref
Search ADS
 

Quinlan
G. D.
,
1996
,
New Astron.
,
1
,
35
10.1016/S1384-1076(96)00003-6

Crossref
Search ADS
 

Rasera
Y.
,
Teyssier
R.
,
2006
,
A&A
,
445
,
1
10.1051/0004-6361:20053116

Crossref
Search ADS
 

Reines
A. E.
,
Comastri
A.
,
2016
,
Publ. Astron. Soc. Aust.
,
33
,
e054
10.1017/pasa.2016.46

Crossref
Search ADS
 

Ricarte
A.
,
Natarajan
P.
,
2018
,
MNRAS
,
481
,
3278
10.1093/mnras/sty2448

Crossref
Search ADS
 

Rosen
A.
,
Bregman
J. N.
,
1995
,
ApJ
,
440
,
634
10.1086/175303

Crossref
Search ADS
 

Schaye
J.
et al. .,
2015
,
MNRAS
,
446
,
521
10.1093/mnras/stu2058

Crossref
Search ADS
 

Schneider
R.
,
Ferrara
A.
,
Natarajan
P.
,
Omukai
K.
,
2002
,
ApJ
,
571
,
30
10.1086/339917

Crossref
Search ADS
 

Sesana
A.
,
Volonteri
M.
,
Haardt
F.
,
2007a
,
MNRAS
,
377
,
1711
10.1111/j.1365-2966.2007.11734.x

Crossref
Search ADS
 

Sesana
A.
,
Haardt
F.
,
Madau
P.
,
2007b
,
ApJ
,
660
,
546
10.1086/513016

Crossref
Search ADS
 

Smith
B. D.
,
Regan
J. A.
,
Downes
T. P.
,
Norman
M. L.
,
O'Shea
B. W.
,
Wise
J. H.
,
2018
,
MNRAS
,
480
,
3762

Crossref
Search ADS
 

Sutherland
R. S.
,
Dopita
M. A.
,
1993
,
ApJS
,
88
,
253
10.1086/191823

Crossref
Search ADS
 

Taffoni
G.
,
Mayer
L.
,
Colpi
M.
,
Governato
F.
,
2003
,
MNRAS
,
341
,
434
10.1046/j.1365-8711.2003.06395.x

Crossref
Search ADS
 

Teyssier
R.
,
2002
,
A&A
,
385
,
337
10.1051/0004-6361:20011817

Crossref
Search ADS
 

Toro
E. F.
,
1997
,
Riemann Solvers and Numerical Methods for Fluid Dynamics
.
Springer Verlag
,
Berlin, Heidelberg

Google Scholar

OpenURL Placeholder Text

 

Trebitsch
M.
,
Volonteri
M.
,
Dubois
Y.
,
Madau
P.
,
2018
,
MNRAS
,
478
,
5607
10.1093/mnras/sty1406

Crossref
Search ADS
 

Tremmel
M.
,
Governato
F.
,
Volonteri
M.
,
Quinn
T. R.
,
2015
,
MNRAS
,
451
,
1868
10.1093/mnras/stv1060

Crossref
Search ADS
 

Tremmel
M.
,
Governato
F.
,
Volonteri
M.
,
Quinn
T. R.
,
Pontzen
A.
,
2018a
,
MNRAS
,
475
,
4967
10.1093/mnras/sty139

Crossref
Search ADS
 

Tremmel
M.
,
Governato
F.
,
Volonteri
M.
,
Pontzen
A.
,
Quinn
T. R.
,
2018b
,
ApJ
,
857
,
L22
10.3847/2041-8213/aabc0a

Crossref
Search ADS
 

Tweed
D.
,
Devriendt
J.
,
Blaizot
J.
,
Colombi
S.
,
Slyz
A.
,
2009
,
A&A
,
506
,
647
10.1051/0004-6361/200911787

Crossref
Search ADS
 

Vasiliev
E.
,
Antonini
F.
,
Merritt
D.
,
2015
,
ApJ
,
810
,
49
10.1088/0004-637X/810/1/49

Crossref
Search ADS
 

Vogelsberger
M.
,
Genel
S.
,
Sijacki
D.
,
Torrey
P.
,
Springel
V.
,
Hernquist
L.
,
2013
,
MNRAS
,
436
,
3031
10.1093/mnras/stt1789

Crossref
Search ADS
 

Volonteri
M.
,
Perna
R.
,
2005
,
MNRAS
,
358
,
913
10.1111/j.1365-2966.2005.08832.x

Crossref
Search ADS
 

Volonteri
M.
,
Haardt
F.
,
Madau
P.
,
2003
,
ApJ
,
582
,
559
10.1086/344675

Crossref
Search ADS
 

Woods
T. E.
et al. .,
2018
, preprint (arXiv e-prints)

© 2019 The Author(s) Published by Oxford University Press on behalf of the Royal Astronomical Society
This article is published and distributed under the terms of the Oxford University Press, Standard Journals Publication Model (https://dbpia.nl.go.kr/journals/pages/open_access/funder_policies/chorus/standard_publication_model)
Issue Section:
Article
Download all slides