Add Journal to My Library
Monthly Notices of the Royal Astronomical Society
, Volume 476 (2) – May 1, 2018

20 pages

/lp/ou_press/active-galactic-nucleus-outflows-in-galaxy-discs-1blHFvp0zo

- Publisher
- Oxford University Press
- Copyright
- © 2018 The Author(s) Published by Oxford University Press on behalf of the Royal Astronomical Society
- ISSN
- 0035-8711
- eISSN
- 1365-2966
- D.O.I.
- 10.1093/mnras/sty229
- Publisher site
- See Article on Publisher Site

Abstract Galactic outflows, driven by active galactic nuclei (AGNs), play a crucial role in galaxy formation and in the self-regulated growth of supermassive black holes (BHs). AGN feedback couples to and affects gas, rather than stars, and in many, if not most, gas-rich galaxies cold gas is rotationally supported and settles in a disc. We present a 2D analytical model for AGN-driven outflows in a gaseous disc and demonstrate the main improvements, compared to existing 1D solutions. We find significant differences for the outflow dynamics and wind efficiency. The outflow is energy-driven due to inefficient cooling up to a certain AGN luminosity (∼1043 erg s−1 in our fiducial model), above which the outflow remains momentum-driven in the disc up to galactic scales. We reproduce results of 3D simulations that gas is preferentially ejected perpendicular to the disc and find that the fraction of ejected interstellar medium is lower than in 1D models. The recovery time of gas in the disc, defined as the free-fall time from the radius to which the AGN pushes the ISM at most, is remarkably short, of the order 1 Myr. This indicates that AGN-driven winds cannot suppress BH growth for long. Without the inclusion of supernova feedback, we find a scaling of the BH mass with the halo velocity dispersion of MBH ∝ σ4.8. black hole physics, methods: analytical, galaxies: active, quasars: general, quasars: supermassive black holes 1 INTRODUCTION The majority of galaxies in the local and distant universe host supermassive black holes (SMBHs) in their centres. The accretion discs around SMBHs are nature's most efficient engines to convert gravitational energy of the infalling gas into radiation. These extreme sources of energy are expected to have a strong feedback effect on the galactic gas content and star formation. AGN feedback may also self-regulate the growth of the central SMBH and its inclusion in models of galaxy formation improves the match between the simulated and the observed galaxy luminosity function for massive galaxies (Bower et al. 2006; Croton et al. 2006; Somerville et al. 2008; Hirschmann et al. 2014; Croton et al. 2016). There are several observations that indicate a correlation between the mass of the central SMBH and large-scale properties of the host galaxy, such as the stellar velocity dispersion σ, luminosity, or the bulge mass (Magorrian et al. 1998; Ferrarese & Merritt 2000; Gebhardt et al. 2000; Tremaine et al. 2002; Marconi & Hunt 2003; Häring & Rix 2004; Gültekin et al. 2009). These relations imply a possible co-evolution of the black hole (BH) and its host galaxy, and AGN feedback has been suggested to be responsible for regulating accretion on the BH and star formation in the galaxy, possibly guiding this correlation (see Heckman & Kauffmann 2011, for a review). The co-evolution, however, is still debated and there are many open questions regarding the nature of AGN feedback (recent reviews on the observational evidence and the theoretical modelling of AGN feedback can be found in Fabian 2012; King & Pounds 2015). Is this relation between the BH mass and galactic properties still valid in galaxies at higher redshift or with lower masses (Greene, Ho & Barth 2008; Baldassare et al. 2015; Reines & Volonteri 2015; Habouzit, Volonteri & Dubois 2017)? How and in which form is the energy from the AGN on sub-pc scales communicated to the interstellar medium (ISM) on kpc scales? Under which conditions can gas escape the gravitational potential of a galaxy? How could the first BHs grow to masses above 109 M⊙ within the first billion years of the universe, although strong radiative feedback is expected to suppress further gas inflow on to the BH? Our current understanding of AGN-driven outflows is based on observations at various wavelengths. The dynamics of the ISM, which is a first indicator of AGN feedback, can be traced by dust, CO, C+, or C ii emission with radio telescopes, such as ALMA. The observations of gas kinematics of high-z galaxies is one of ALMA's main science drivers and it can reach sub-kpc resolution at z = 6 (Wang et al. 2013; Trakhtenbrot et al. 2017; Decarli et al. 2018). Radio observations reveal large molecular outflows with velocities of 100 km s−1 up to a few times 1000 km s−1 (Aalto et al. 2012; Cicone et al. 2014). Recently, Maiolino et al. (2017) report star formation inside a galactic outflow, another strong indication for the presence of an outflow of cold gas. In Mrk231, the closest and best-studied quasar, we observe velocities of ∼1000 km s−1 and mass outflow rates of the order ∼1000 M⊙ yr−1 (Feruglio et al. 2010; Rupke & Veilleux 2011; Sturm et al. 2011). These outflows have to be powered by the AGN, because SN-driven outflows can only account for outflow velocities up to ∼600 km s−1 (Martin 2005; Sharma & Nath 2013). Besides these molecular outflows, the Fe K lines in the X-ray reveal highly ionized AGN outflows with mildly relativistic velocities of ∼0.1c (Chartas et al. 2002; Pounds et al. 2003; Reeves, O'Brien & Ward 2003; Cappi 2006; Gofford et al. 2013). Some ultrafast outflows even have velocities up to ∼0.3c, but the majority has velocities around ∼0.1c (King & Pounds 2015). The consensual picture of AGN-driven outflows consists of an inner, line-driven disc wind with mildly relativistic velocities (also called ‘nuclear wind’). This highly ionized wind is shocked once it encounters the denser ISM and drives a forward shock into the ISM, which sweeps up most of the gas. As we discuss below, the details of the driving mechanism depend on the efficiency of cooling in the shocked wind: if the internal energy is radiated away, only the momentum input from the AGN drives the outflow (‘momentum-driven’). If cooling is inefficient, the hot shocked wind adiabatically expands and drives the shock into the ISM (‘energy-driven’). These different driving mechanisms have fundamental influence on the outflow dynamics and the thermal state of the outflow. Most of the thermal energy is carried by the shocked wind, and the shocked ISM accounts for the bulk part of the kinetic energy. Several groups have developed 1D analytical models of AGN-driven outflows in these two different regimes (Silk & Rees 1998; Fabian 1999; King 2003; Murray, Quataert & Thompson 2005; Silk & Nusser 2010; Faucher-Giguère & Quataert 2012; Ishibashi & Fabian 2012; Zubovas & King 2012). While these models can reproduce several observations, the assumption of a spherically symmetric gas distribution is too simplistic to capture the complexity of a realistic galaxy. Moreover, these models often assume an energy- or momentum-driven outflow, without self-consistently solving for the corresponding transition. For a momentum-driven outflow, you find a characteristic BH mass of King & Pounds (2015) \begin{eqnarray} M_\sigma = \frac{f_\mathrm{gas}\kappa }{\pi G} \sigma ^4, \end{eqnarray} (1)with the gas opacity κ, above which the AGN can eject the gas out of the gravitational potential of its host galaxy and hence prevents further BH growth. The ∝ σ4 scaling is independent of the underlying density distribution (McQuillin & McLaughlin 2012). If, on the other hand, we assume an energy-driven outflow, this characteristic mass and its scaling with the velocity dispersion is given by \begin{eqnarray} M_\sigma = \frac{11 f_\mathrm{gas} \kappa }{\eta \pi G^2 c} \sigma ^5, \end{eqnarray} (2)where η is proportional to the mechanical luminosity (Haehnelt, Natarajan & Rees 1998; Silk & Rees 1998; Fabian 2012; McQuillin & McLaughlin 2013). Observationally, different values for the exponent α of the M–σ relation, MBH ∝ σα, are proposed, ranging from α = 3.7 (Gebhardt et al. 2000) to α = 5.6 (McConnell & Ma 2013) with other values within this range (see e.g. Ferrarese & Merritt 2000; Tremaine et al. 2002; Gültekin et al. 2009; Kormendy & Ho 2013). Numerical simulations of AGN-driven outflows in 3D demonstrate that the outflow takes the path of least resistance and propagates preferentially towards the steepest density gradient, i.e. along the poles of the host galaxy (Gaspari et al. 2011a; Gabor & Bournaud 2013, 2014; Gaspari, Brighenti & Ruszkowski 2013; Costa et al. 2014; Costa, Sijacki & Haehnelt 2015; Roos et al. 2015; Bieri et al. 2017; Barai et al. 2018). This results in highly anisotropic outflows as found by Costa et al. (2014) and observationally confirmed by Cicone et al. (2015). However, 3D simulations remain expensive, have explored limited parameter space, and have not yet focused on low-mass or high-redshift galaxies. Our aim is to develop a flexible, but general model to study AGN feedback in galaxies where rotational support in the gas component is non-negligible. We derive a 2D analytical model for AGN-driven outflows: two dimensions are sufficient to capture the geometry of a disc-like galaxy, which allows simultaneously an AGN-feeding inflow and an ISM-ejecting outflow, but it is still simple enough to be treated analytically and flexible enough to explore a large parameter space with this model. We include a self-consistent transition from momentum- to energy-driving and derive observables, such as velocities, mass outflow rates, and momentum loading. 2 METHODOLOGY In this section, we describe the basic equations that govern our 2D analytical approach, justify necessary approximations, and introduce our physically motivated galaxy model. 2.1 Galaxy model We present and motivate the individual components and characteristic radii of our analytical approach. 2.1.1 Gas density profile of the disc Our main focus is to model AGN-driven outflows in low-mass galaxies, which represent the first galaxies at high redshift or dwarf galaxies at lower redshift (Penny et al. 2018). We assume that the gas distribution in these galaxies can be approximated by a 2D axisymmetric density profile. The surface brightness of many spiral galaxies can be described by an exponential profile (Sérsic 1963; Freeman 1970; Courteau, de Jong & Broeils 1996). Smit et al. (2018) recently report the observation of two z ≈ 6.8 galaxies with a clear velocity gradient in C ii, which might suggest the presence of a turbulent, rotation-dominated disc. Based on imaging spectroscopy of the H α emission line, Genzel et al. (2017) confirm that the observed velocity profiles of high-z disc galaxies favour a thick exponential disc profile and Pawlik, Milosavljević & Bromm (2011) also find a disc structure for the first galaxy in their 3D hydrodynamical simulations. We therefore assume an exponential profile for the surface density \begin{eqnarray} \Sigma (R) = \Sigma _0 \exp (-R/R_0), \end{eqnarray} (3)where R0 is the scaling radius, given by (Mo, Mao & White 1998) \begin{eqnarray} R_0 = \frac{\lambda }{\sqrt{2}} \left( \frac{j_d}{m_d} \right) R_\mathrm{vir}, \end{eqnarray} (4)with jd being the fraction of the specific angular momentum of the halo that goes into the disc, md the fraction of the halo mass that settles into the disc, and λ the spin parameter of the halo. For typical haloes, we assume jd ≈ md ≈ λ ≈ 0.05 and hence \begin{eqnarray} R_0 \approx 0.035 R_\mathrm{vir}. \end{eqnarray} (5)The scale height of an isothermal, self-gravitating disc is given by \begin{eqnarray} H(R) = \frac{c_s ^2}{\pi G \Sigma (R)} = \frac{c_s ^2}{\pi G \Sigma _0} \exp (R/R_0), \end{eqnarray} (6)where cs is the sound speed of the gas. The derivation of this relation assumes a locally constant surface density and the real scale height might be smaller in regions where the disc surface density changes significantly with radius. For the total density profile of the disc, we require \begin{eqnarray} \Sigma (R) = \int _{- \infty } ^{\infty } \rho (R,h)\, \mathrm{d}h. \end{eqnarray} (7)By assuming \begin{eqnarray} \rho (R,h) \propto \cosh ^{-2} (h/H(R)) \end{eqnarray} (8)we find \begin{eqnarray} \rho (R,h) = \frac{\Sigma _0}{2 H(R)} \exp (-R/R_0) \cosh ^{-2} (h/H(R)). \end{eqnarray} (9)The gas mass inside a radius R for the exponential disc is given by \begin{eqnarray} M(<R) = 2 \pi \Sigma _0 R_0 \left[ R_0 - \exp (-R/R_0) (R+R_0) \right]. \end{eqnarray} (10) For reasonable values of the spin parameter λ, the second term is negligible. Assuming that mdMh is the mass of the disc within Rvir, we can calculate the normalization of the surface density as \begin{eqnarray} \Sigma _0 = \frac{m_d M_h}{\pi (\lambda R_\mathrm{vir})^2}, \end{eqnarray} (11)where Mh is the mass of the halo. 2.1.2 Dark matter Following Costa et al. (2014), we assume a spherically symmetric Hernquist profile (Hernquist 1990) of the form \begin{eqnarray} \rho _{\text{DM}} = (1-f_\mathrm{gas})\frac{M_h}{2 \pi }\frac{a}{r(r+a)^3} \end{eqnarray} (12)and \begin{eqnarray} M_{\text{DM}} (<r) = (1-f_\mathrm{gas}) M_h \frac{r^2}{(r+a)^2}, \end{eqnarray} (13)with fgas = Ωb/Ωm = 0.17 and the scaling radius a = 0.1Rvir. With small r we denote the radius in spherical geometry in contract to capital R that we use for the radius in the disc plane in cylindrical coordinates. 2.1.3 1D reference model To highlight the differences of a 2D model, we also compare our new model to a spherical 1D gas distribution for which we assume that the radial gas profiles follows a Hernquist profile with \begin{eqnarray} \rho _\mathrm{gas} = \frac{f_\mathrm{gas}M_h}{2 \pi }\frac{a}{r(r+a)^3}, \end{eqnarray} (14)where a = 0.1Rvir is the scaling radius. 2.1.4 Fiducial galaxy model For our fiducial model, we assume typical parameters of a z = 6 galaxy with a halo mass of Mh = 108 M⊙. The BH has a mass of MBH = 105 M⊙ and shines at an Eddington ratio of fEdd = 0.3. We assume that the AGN shines for a time of ton = 10 Myr with a constant luminosity that is given by the BH mass and the Eddington ratio. The resulting surface density normalization and disc scale height at the scale radius are Σ0 = 0.28 g cm−2 and H(R0) = 2.6 pc, respectively. In Section 3.2, we vary all these model parameters independently to analyse their influence on the outflow dynamics. We do not include a stellar population in the galaxy because its gravitational potential has no direct influence on the outflow dynamics. Galaxies at z ≥ 6 have a stellar mass fraction of at most $$10\,\,\rm{per\,\,cent}$$, relative to their baryon content (O'Shea et al. 2015; Trebitsch et al. 2017). Low-mass galaxies with a halo mass of Mh ≈ 108 M⊙, such as our fiducial model, have even lower stellar mass fractions of $$0.1\rm{-}1\,\,{per\,\,cent}$$ (Behroozi, Wechsler & Conroy 2013; Trebitsch et al. 2017). In one test scenario, we include a stellar population as a Hernquist profile with a scaling radius of a* = 5 pc (Park et al. 2016) with a stellar mass fraction of $$10\,\,\rm{per\,\,cent}$$. Even under this conservative assumption, the time evolution of the AGN-driven outflow is not distinguishable from another model, in which we reduce the gas mass by 10 per cent. This is because the dynamics of the outflow is dominated by the inertia of the swept-up gas and not by the gravity of the external potential. Phrased differently, the inclusion of a static stellar potential does not affect the outflow and only the implicated lower gas mass has a small effect. A smaller gas mass fraction, e.g. due to star formation at lower redshift, would decrease the surface density normalization. Also, gas in high-z galaxies is turbulent and has significant pressure support (Genzel et al. 2008; Förster Schreiber et al. 2009; Safranek-Shrader et al. 2012), which would additionally increase the disc scale height. These effects would lower the ISM density in the disc plane and permit the AGN-driven wind to propagate to larger radii. 2.1.5 Initial and final radii We start the integration of the outward-driven shock at an initial radius Rmin. The physical interpretation of this radius is the innermost distance where the disc wind starts to sweep up the interstellar medium (ISM). We use the self-gravity radius of the accretion disc as a proxy for Rmin and we express it in units of the Schwarzschild radius: \begin{eqnarray} R_\mathrm{SS} = \frac{2GM_\mathrm{BH}}{c^2} \approx 3\times 10^5 \frac{M_\mathrm{BH}}{\,\mathrm{M}_{\odot }}\, \mathrm{cm}. \end{eqnarray} (15)For standard values of the disc viscosity and Eddington ratio, the outer radius of the disc (i.e. self-gravity radius) and hence the inner radius for our integration is given by (Shakura & Sunyaev 1973; Goodman & Tan 2004; King, Pringle & Hofmann 2008; Dotti et al. 2010) \begin{eqnarray} R_\mathrm{min} \approx 10^{5} R_\mathrm{SS}, \end{eqnarray} (16)which yields consistent results for all the tested parameter combinations and for BH masses, spanning many orders of magnitude. We integrate the equations out to arbitrary radii, but the physically interesting regime is within the virial radius. Hence, we follow the outflow in each direction until it either becomes subsonic, or until it crosses the virial radius of the halo. Since the escape velocity at the virial radius is roughly equal to the sound speed, we assume that the wind escapes the gravitational potential of the halo if it reaches the virial radius and is still supersonic. 2.2 Quantifying the outflow We define a set of parameters that quantify the efficiency, dynamics, and nature of the outflow. We introduce a radius, Rmax, out to which the outflow has to push the gas at least to stop further gas accretion on to the central BH. Due to the higher gas density in the disc plane, it is a sufficient criterion for the outflow to reach Rmax in the disc plane. A sufficient criterion to prevent further gas accretion on to the central BH is Rmax = Rvir and a necessary criterion is Rmax > Rsoi, where Rsoi is the sphere of influence, within which the dynamics is dominated by the BH. Once pushed out to a certain radius, the question is if this gas can reach the BH in a sufficiently short amount of time to start refuelling the AGN. We hence relate Rmax to the recovery time tff, which we define as the minimum time the gas needs to fall back towards the central BH, after being pushed out to the radius Rmax. Expressing the free-fall time via the sound crossing time yields for Rmax in the disc plane \begin{eqnarray} t_\mathrm{ff} = 5\,\mathrm{Myr} \left( \frac{R_\mathrm{max}}{R_0} \right) \left( \frac{1+z}{7} \right)^{-3/2}. \end{eqnarray} (17)There is no unique criterion on the absolute value of Rmax or tff to shut off further gas accretion. We therefore present and discuss both relative to the scale radius of the disc and to the lifetime of the AGN. If tff is greater than the AGN lifetime, we can assume that the outflow is efficient enough to significantly suppress further mass accretion. The ejected mass Meject is the gas mass that is ejected from the gravitational potential of the halo, i.e. that passes the virial radius with a velocity larger than the escape velocity. The ejection angle θeject is the last angle (measured from the disc normal) for which gas can escape the halo. It defines the opening angle of a cone, which is cleared off gas by the AGN-driven wind. Gas outside this cone is pushed to Rmax(θ), but falls back towards the galactic centre over a time ≳ tff. The free-fall time in the disc plane is defined as the time which gas that has been pushed to a certain radius needs at least to fall back towards the central BH, once the shock becomes subsonic (see equation 17). This is an approximate quantification for how long an AGN can shut off its own gas supply. We also compare the mass accretion rate of the central BH: \begin{eqnarray} \dot{M}_\mathrm{acc} = \frac{L}{\eta c^2}, \end{eqnarray} (18)to the mass outflow rate: \begin{eqnarray} \dot{M}_\mathrm{out} = \sum _i M_\mathrm{shell,i} \frac{v_\mathrm{shell,i}}{R_\mathrm{shell,i}}, \end{eqnarray} (19)where the mass and velocity of the shell resolution elements are summed over all angles i. The mass outflow rate generally varies with time, but to make it more accessible for a direct comparison to the mass accretion rate, we define the mass outflow rate at the sphere of influence, which is one unique number. In theory, the shell crosses the sphere of influence at different times for different angles and we need to sum over shell elements at different ‘crossing-times’. However, in practice, the shell crosses the sphere of influence for all angles in the same time-step, because the density profile of the exponential disc is close to spherical in the centre and there is no preferred direction. Comparing these two mass rates quantifies the efficiency of the AGN to convert infalling, accreted mass into a mass outflow. The ratio Rperp/Rdisc quantifies the asymmetry of the outflow. At any time, Rdisc is the position of the shock front in the disc plane and Rperp perpendicular to it. If this ratio is equal to one, the wind propagates almost spherically symmetric. For values above unity, the outflow develops peanut-shaped lobes. This ratio illustrates the importance and improvement of a 2D treatment with respect to 1D wind solutions. We define the momentum boost, which is often referred to as ‘mechanical advantage’, as the ratio of the total momentum of the AGN-driven wind divided by the total momentum input by photons: \begin{eqnarray} \frac{\sum _i |{\boldsymbol p}_i|}{tL/c}. \end{eqnarray} (20)The vectorial sum would be equal to zero, but summing the absolute values of the momentum allows to quantify the efficiency of the AGN to accelerate the gas. In the momentum-driven regime we expect \begin{eqnarray} \sum _i |{\boldsymbol p}_i| \le tL/c, \end{eqnarray} (21)but the adiabatic expansion in the energy-driven case can yield mechanical advantages above unity. Observations find momentum boost in the range of ∼2–30 (Rupke & Veilleux 2011; Sturm et al. 2011) and 3D simulations yield momentum boosts of 1–30 (Cicone et al. 2014; Bieri et al. 2017). 2.3 Shock acceleration There are two systematically different types of outflows: if the shocked wind cools efficiently, the outflow is accelerated only by the momentum transfer and the outflow is called momentum-driven. A momentum-driven outflow develops a thin shock layer. For less efficient cooling, the shock-heated gas remains at temperatures of >107 K and the adiabatic expansion of the hot shocked wind accelerates the shock into the ISM. In this case, a thick shock layer develops and the outflow is labelled energy-driven, because the internal energy of the adiabatically expanding shock-heated wind drives the outflow. In a disc-like geometry, we expect the shock to propagate faster into the direction of lower density. This implies that a self-consistent 2D treatment of the wind propagation is important to capture all the relevant physics of the shock dynamics. We can assume that the shock front is locally plane-parallel and in the rest frame of the shock front, the Rankine–Hugoniot jump conditions for a plane parallel shock are as follows: \begin{eqnarray} \rho _1 u_1 &=& \rho _2 u_2, \end{eqnarray} (22) \begin{eqnarray} \rho _1 u_1 ^2 +P_1 &=& \rho _2 u_2 ^2 +P_2, \end{eqnarray} (23) \begin{eqnarray} \frac{1}{2} \rho _1 u_1 ^2 + \epsilon _1 +\frac{P_1}{\rho _1} = \frac{1}{2} \rho _2 u_2 ^2 + \epsilon _2 +\frac{P_2}{\rho _2}, \end{eqnarray} (24)where the index 1 refers to the pre-shock and the index 2 to the post-shock conditions of the density ρ, velocity u, pressure P, and the specific internal energy ε. For a strong shock with Mach number \begin{eqnarray} \mathcal {M}_1 = \left( \frac{\rho _1 u_1 ^2}{\gamma P_1} \right) ^{1/2} \gg 1, \end{eqnarray} (25)we find \begin{eqnarray} \frac{\rho _2}{\rho _1} \approx \frac{\gamma + 1}{\gamma - 1} = 4, \end{eqnarray} (26) \begin{eqnarray} P_2 \approx \frac{2}{\gamma +1} \rho _1 u_1 ^2 = \frac{4}{3} \rho _1 u_1 ^2, \end{eqnarray} (27) \begin{eqnarray} T_2 \approx \frac{2(\gamma -1)}{(\gamma +1)^2} \frac{m}{k_\mathrm{B}} u_1 ^2 = \frac{3}{16} \frac{m}{k_\mathrm{B}} u_1 ^2, \end{eqnarray} (28) \begin{eqnarray} \mathcal {M}_2 \approx \left( \frac{\gamma -1}{2 \gamma } \right) ^{1/2} \approx 0.45, \end{eqnarray} (29)where the last equalities are valid for an adiabatic coefficient of γ = 5/3 (Shull & Draine 1987). Hence, a shock converts supersonic gas into denser, slower moving, higher pressure, subsonic gas. The density of the shocked wind depends on the available cooling channels. For a radiative (momentum-driven) strong shock, the post-shock density is given by ρ2 = 4ρ1 (equation 26). In the energy-driven case, the density of the shocked wind depends on the thickness of the shock. We generally refer to the shock position as the contact discontinuity between the shocked wind and the shocked ISM, and the thickness of the shocked wind is given by \begin{eqnarray} \Delta r = r - r_\mathrm{sw}, \end{eqnarray} (30)where R and Rsw are the position of the contact discontinuity and of the wind shock, respectively (see illustrations of outflow dynamics in Faucher-Giguère & Quataert 2012; Zubovas & King 2012; King & Pounds 2015; Dashyan et al. 2018). The radius where the wind is shocked for an energy-driven outflow is given by Faucher-Giguère & Quataert (2012): \begin{eqnarray} r_\mathrm{sw} \approx r \sqrt{\frac{v}{v_\mathrm{in}}}, \end{eqnarray} (31)and the thickness of the shocked wind is Δr = r(1 − ($$v$$/$$v$$in)1/2). We do not use this thickness of the shock explicitly for the integration of the equations of motion, but we use it implicitly to calculate the rate of radiative cooling (Section 2.4.3). 2.3.1 Equation of motion: momentum-driven The shock is accelerated by the momentum input of the AGN, counteracted by gravity (King 2010; Ishibashi & Fabian 2014, 2015; King & Pounds 2015): \begin{eqnarray} \frac{\mathrm{d}}{\mathrm{d}t} \left[ M_\mathrm{shell}(r) \dot{r} \right] &=& \frac{f_\mathrm{edd} L_\mathrm{edd}}{c} \nonumber \\ &&-\,G \frac{M_\mathrm{shell}(r)(M_\mathrm{DM}(<r)+M_\mathrm{BH})}{r^2} - F_\mathrm{self\text{-}grav}, \nonumber\\ \end{eqnarray} (32)where fedd is the Eddington fraction, Mshell is the mass of the swept-up ISM, and \begin{eqnarray} L_\mathrm{edd} = \frac{4 \pi c G M_\mathrm{BH} m_p}{\sigma _T} \end{eqnarray} (33)is the Eddington luminosity. The term $$F_\mathrm{self\text{-}grav}$$ is not explicitly mentioned in most studies (King 2010; Ishibashi & Fabian 2014, 2015; King & Pounds 2015), but it accounts for the swept-up gas in the shell and the gravity contribution from the gas outside the shell. The gravitational force of the ambient gas is negligible, but self-gravity of the shell has to be included. The term L/c describes the total momentum per unit time that can be transferred to the gas via single scatter events. For Thomson scattering, we assume that the directions of the photons after one scattering are random and that the momentum transfer from secondary scattering events cancels out. However, there is a probability 0 < p ≤ 1 for the photons to escape freely without any interaction with the ISM: \begin{eqnarray} \ln (p)= - \frac{\sum _i \sigma _{t,i}}{\text{d} \Omega r^2} = - \frac{M_\mathrm{shell}}{m_p} \frac{\sigma _t}{d \Omega r^2} = - \kappa \Sigma _\mathrm{shell} = - \tau _\mathrm{shell}, \end{eqnarray} (34)with the opacity κ = σt/mp = 0.4 cm2 g−1. The fraction of the total force input L/c that actually couples to the gas is hence \begin{eqnarray} \frac{L}{c} \left( 1- e^{- \tau _\mathrm{shell}} \right). \end{eqnarray} (35)We do not explicitly account for this effect and show in Section 3.1 that it is negligibly small and irrelevant in most of the scenarios considered. However, see Ishibashi & Fabian (2015) for a more detailed discussion of this effect and for the influence of multiscattering of IR photons on dust grains. 2.3.2 Equation of motion: energy-driven In the energy conserving case, the shocked ISM is not accelerated by the direct momentum input but by the adiabatic expansion of the shocked wind. The momentum equation then reads \begin{eqnarray} \frac{\mathrm{d}}{\mathrm{d}t}[M_\mathrm{shell}(r)\dot{r}]+F_\mathrm{grav}=4 \pi r^2 P, \end{eqnarray} (36)where P is the pressure in the shocked wind, which can be calculated via energy conservation \begin{eqnarray} (\gamma -1) \frac{\mathrm{d}}{\mathrm{d}t} (PV) = \frac{\eta }{2} f_\mathrm{Edd} L_\mathrm{Edd} - P \frac{\mathrm{d}}{\mathrm{d}t}V - F_\mathrm{grav} \dot{r}, \end{eqnarray} (37)where V is the volume and η sets the mechanical luminosity, which can be expressed as η = $$v$$in/c (see Section 2.5). This energy conservation assumes that the thermal energy of the shocked wind is much higher than its kinetic energy, which is generally the case (Faucher-Giguère & Quataert 2012). 2.4 Cooling and heating In this section, we discuss the most important heating and cooling mechanisms of the gas. We demonstrate that inverse Compton cooling is the only relevant mechanism on the scales of interest (Ciotti & Ostriker 1997) and discuss why radiative cooling of the wind shock can be neglected. 2.4.1 Compton cooling The Compton cooling time of gas is given by \begin{eqnarray} t_\mathrm{compton} = \frac{3 m_e c}{8 \pi \sigma _T U_\mathrm{rad}} \frac{m_e c^2}{E}, \end{eqnarray} (38)with the radiation energy density \begin{eqnarray} U_\mathrm{rad} = \frac{f_\mathrm{Edd} L_\mathrm{Edd}}{4 \pi R^2 c} \end{eqnarray} (39)and the internal energy in the shocked gas \begin{eqnarray} E= \frac{9 m_p v_\mathrm{in}^2}{16}, \end{eqnarray} (40)where $$v$$in is the initial velocity of the disc wind before it shocks with the ISM (see Section 2.5), and mp and me are the masses of the proton and electron, respectively. The Compton cooling time is hence given by (see also Faucher-Giguère & Quataert 2012) \begin{eqnarray} t_\mathrm{compton} = \frac{2}{3 \pi f_\mathrm{edd}} \frac{c}{G M_\mathrm{BH}} \left( \frac{m_e}{m_p} \right) ^2 \left( \frac{v_\mathrm{in}}{c} \right) ^{-2} r^2. \end{eqnarray} (41) 2.4.2 Two temperature medium Inverse Compton cooling acts only on the electrons. Faucher-Giguère & Quataert (2012) analyse the cooling properties of a medium with different temperatures of the electrons and ions. They determine the characteristic time-scale to achieve thermal equilibrium (Te = Tp) to be \begin{eqnarray} t_{ei} = \frac{3 m_e m_p}{8 (2 \pi )^{1/2} n_p e^4 \ln \Lambda } \left( \frac{k_B T_e}{m_e} \frac{k_B T_p}{m_p} \right)^{3/2}, \end{eqnarray} (42)where e is the elementary charge and Λ ≈ 40. Assuming Tp = 10Te = 1010 K and a proton density of np = 1 cm−3 as typical post-shock conditions, this time-scale is of the order tei ≈ 1 Myr. Faucher-Giguère & Quataert (2012) show that this time-scale is longer than the Compton cooling time and might delay the cooling by up to two orders of magnitude. We do not include this effect of a two-temperature medium in our calculations and discuss its effect in Section 4.5. 2.4.3 Radiative cooling We only consider inverse Compton cooling in our calculation of the cooling time. To quantify the possible contribution of radiative cooling, we use the cooling function by Sutherland & Dopita (1993) in the functional form by Tozzi & Norman (2001). The post-shock temperature is given by equation (28) and we calculate the gas density based on equation (30) or equation (26), depending on whether the shock is energy- or momentum-driven, respectively. The radiative cooling time becomes \begin{eqnarray} t_\mathrm{rad} = \frac{9 v_\mathrm{in}^2 \mu m_\mathrm{p}^2}{32 \Lambda \rho _\mathrm{shock}}, \end{eqnarray} (43)where ρshock is the corresponding post-shock density and Λ the cooling function. The radial profile of the radiative cooling time and a comparison to other characteristic time-scales can be seen in Fig. 1. Within ∼0.01 pc, the velocity of the swept-up material is close to the initial wind velocity and therefore the shock would be infinitely thin (equation 31). This implies a very high density and an artificially short cooling time. In this case, we set the post-shock density to ρshock = 4ρ to avoid singularities. On all scales of interest, the radiative cooling time is always longer than the Compton cooling time and can be neglected. Figure 1. View largeDownload slide Radial profiles of characteristic time-scales for our fiducial model, where the flow time is defined as tflow = R/$$v$$. The different colours indicate different angles with respect to the disc normal, θ, and the Compton cooling time in this radial range is independent of the angle (grey). Figure 1. View largeDownload slide Radial profiles of characteristic time-scales for our fiducial model, where the flow time is defined as tflow = R/$$v$$. The different colours indicate different angles with respect to the disc normal, θ, and the Compton cooling time in this radial range is independent of the angle (grey). 2.5 Mechanical luminosity and initial wind velocity The post-shock temperature is set by the pre-shock wind velocity $$v$$in (equation 28). By assuming that the initial wind mass rate $$\dot{M}_\mathrm{in}$$ is equal to the accretion rate on to the BH $$\dot{M}_\mathrm{BH}$$, with a radiative efficiency of η ≈ 0.1 for the AGN, we obtain \begin{eqnarray} v_\mathrm{in} = \eta c \approx 30\,000\, \mathrm{km}\, \mathrm{s}^{-1}. \end{eqnarray} (44)The assumption $$\dot{M}_\mathrm{in}=\dot{M}_\mathrm{BH}$$ might not be accurate, but we use it as a working hypothesis, following e.g. Faucher-Giguère & Quataert (2012), Costa et al. (2014), and King & Pounds (2015), and we explore the dependence on $$v$$in in Section 3.2.3. The assumption of $$v$$in = 0.1c is supported by observations of the Doppler-shifted Fe K-band absorption lines, which are observed in many local AGNs (Pounds et al. 2003; Tombesi et al. 2010; Tombesi 2016). Ostriker et al. (2010) provide a more detailed and complete discussion of the momentum driving and mechanical luminosity. This assumption enters in the calculation of the cooling time and for the mechanical input of the AGN. If we assume that the initial mass outflow rate is only a fraction of the accretion rate: \begin{eqnarray} \dot{M}_\mathrm{in}=f_\mathrm{acc}\dot{M}_\mathrm{BH}, \end{eqnarray} (45)we would find $$v$$in = ηfaccc for the initial velocity and the energy equation would change to \begin{eqnarray} \frac{\eta }{2} L_\mathrm{AGN} \rightarrow \frac{\eta }{2} f_\mathrm{acc} L_\mathrm{AGN} = \frac{\eta }{2} f_\mathrm{acc} f_\mathrm{Edd} L_\mathrm{Edd}. \end{eqnarray} (46)Hence, in the parameter study of the outflow dynamics, facc and fEdd enter the equation of motion as a product, which allows us to keep one parameter fixed (facc = 1) and vary only the Eddington ratio. We assume that this initial line-driven disc wind is isotropic. This is supported by observations of ultrafast outflows in a large fraction of observed AGNs. Unless we see all AGNs from a very particular angle, these AGN-driven winds have large opening angles with an almost spherical driving mechanism (King & Pounds 2015). 2.6 Momentum- to energy-driven transition The main problem arises since the nature of the outflow depends on the dynamics of the shock and the shock dynamics depends on the nature of the outflow. We self-consistently check at any time if the outflow is energy- or momentum-driven (see Fig. 2). Figure 2. View largeDownload slide Characteristics times as a function of radius for different angles with respect to the disc normal for our fiducial model. The black dashed line shows the Compton cooling time (equation 41) and the solid lines show the flow time, for the self-consistent driving, i.e. energy-driven when tflow < tCompton. The dashed lines illustrate the dynamics of an outflow that is always momentum-driven. The transition from momentum- to energy-driven occurs at ∼2 × 10−3 pc and the momentum-driven solution is less efficient in driving an outflow. Figure 2. View largeDownload slide Characteristics times as a function of radius for different angles with respect to the disc normal for our fiducial model. The black dashed line shows the Compton cooling time (equation 41) and the solid lines show the flow time, for the self-consistent driving, i.e. energy-driven when tflow < tCompton. The dashed lines illustrate the dynamics of an outflow that is always momentum-driven. The transition from momentum- to energy-driven occurs at ∼2 × 10−3 pc and the momentum-driven solution is less efficient in driving an outflow. We expect the transition to happen at different radii for different angles with respect to the disc plane. This implies that there are phases where momentum- and energy-driving occur at the same time in a given galaxy. For the energy input by the AGN in driving the adiabatic expansion, we correct it by the solid angle dΩ of the energy-driven outflow, because for $$4\pi - \mathrm{d}\Omega$$ we expect the energy input to be radiated away in the momentum-driven shock front. However, the remaining question is, which volume represents the correct basis for the energy equation of the adiabatic expansion: the entire volume enclosed by the shock front or only the volume enclosed by the energy-driven shock front? A Gedanken experiment can help to settle this question: if we only had to consider the volume of the cone enclosed by the energy-driven shock (perpendicular to the disc plane, because we expect higher velocities in this direction), then what happens at the boundary of this cone? The adiabatically expanding volume does not know about whether it is enclosed by an energy- or momentum-driven shock front, but rather wants to expand in all directions. For the same reason, we have to adopt our treatment of the adiabatically expanding volume once the shock becomes subsonic: the gas still expands and can do so into the ISM with ∼cs, even if the shock front is dissolved by turbulence. For directions into which the shock becomes subsonic, we hence follow the propagation of the swept-up mass and the adiabatically expanding volume separately. We have tested the influence of these methods and find that the difference is small. 2.7 Validation of the model: comparison to 1D solution To test the robustness and convergence of the method, we first compare the results of our 2D analytical model to the results of Costa et al. (2014) and King (2005). We assume Mh = 1012 M⊙ for the halo mass, a Hernquist profile for the dark matter and the gas, where the gas mass is Mgas = fgasMh with fgas = 0.17 and the scaling radius a = 28 kpc. The virial radius is Rvir = 163 kpc, the AGN shines at the Eddington luminosity and for the BH we assume three different masses of MBH = 5 × 107, 108, 3 × 108 M⊙. We assume that there is no gas within the sphere of influence of the BH, which is ∼500 pc for these configurations. The results can be seen in Fig. 3. Figure 3. View largeDownload slide Comparison of our model with the results by King (2005), which Costa et al. (2014) can also reproduce with their 3D simulation, for different BH masses. Here, we assume the same spherical density distribution as Costa et al. (2014), but integrate the equations of motion with our 2D model. Our results agree well with the analytical solution by King (2005) and with the 3D simulations by Costa et al. (2014). Figure 3. View largeDownload slide Comparison of our model with the results by King (2005), which Costa et al. (2014) can also reproduce with their 3D simulation, for different BH masses. Here, we assume the same spherical density distribution as Costa et al. (2014), but integrate the equations of motion with our 2D model. Our results agree well with the analytical solution by King (2005) and with the 3D simulations by Costa et al. (2014). Our 2D analytical model is able to reproduce the outflow dynamics of the 1D model by King (2005) and of the 3D simulation by Costa et al. (2014) in a spherical gas profile. We then adopt the same functional forms for the gas distribution (Hernquist profile) and apply it to our fiducial model with Mh = 3 × 107 M⊙ and MBH = 105 M⊙. The direct comparison of the outflow dynamics in 1D and 2D is given in Fig. 4. For the same gas mass and AGN luminosity, the outflow in the 2D model can only eject mass within a cone with an opening angle of ∼45°, whereas the same outflow in the 1D model ejects all the gas out of the virial radius. The main reason for this lower efficiency in the 2D scenario is that the pressure in the 1D model cannot escape, but builds up over time. In the 2D scenario, however, the pressure decreases with time due to the fast expansion of the energy-driven wind perpendicular to the disc. Figure 4. View largeDownload slide Time evolution of the position of the shock front in polar coordinates. Top: 2D solution with the galactic disc plane at 0°. Bottom: 1D spherical solution. The concentric grey rings indicate the radius in log (r/pc) from 10 pc to 100 kpc and the magenta circle represents the virial radius at 980 pc. In the disc plane, the shock cannot propagate beyond a maximal radius of ∼20 pc and only gas at ≳47° escapes the virial radius. The 1D wind is slower at ≲30 kyr, but ejects more gas out of the virial radius at later times. Figure 4. View largeDownload slide Time evolution of the position of the shock front in polar coordinates. Top: 2D solution with the galactic disc plane at 0°. Bottom: 1D spherical solution. The concentric grey rings indicate the radius in log (r/pc) from 10 pc to 100 kpc and the magenta circle represents the virial radius at 980 pc. In the disc plane, the shock cannot propagate beyond a maximal radius of ∼20 pc and only gas at ≳47° escapes the virial radius. The 1D wind is slower at ≲30 kyr, but ejects more gas out of the virial radius at later times. 3 RESULTS 3.1 Standard case (fiducial parameters) We show the time evolution of the main dynamical quantities for our fiducial model in Fig. 5. Within ∼0.01 pc the swept-up mass and hence the gas inertia is negligibly small and the dynamics is dominated by the impinging disc wind, which has an initial velocity of $$v$$in = 0.1c in our fiducial model. Our model is insensitive to the choice of the initial radius of integration Rmin, as long as it is within this disc wind-dominated region. When the cooling time becomes longer than the flow time, the pressure rises by about 10 orders of magnitude and the outflow becomes energy-driven. The outflow slows down due to the inertia of the swept-up ISM, which was at rest first. Within ∼10 pc the density profile is almost spherically symmetric, but we can see a clear difference for the outflow beyond this radius for different angles with respect to the disc plane. While the shock front in the disc plane cannot escape the virial radius, the gas perpendicular to the disc plane gets ejected out of the galaxy. Figure 5. View largeDownload slide Time evolution of different quantities for our fiducial galaxy model. Top panel: position of the shock for different angles with respect to the disc. The virial and disc scaling radius are indicated by the black dashed lines. Middle panel: velocity for different angles. The wind becomes subsonic in the disc plane at around 0.5 Myr. Lower panel: the pressure of the shocked wind is isotropic and it rises rapidly once the gas can no longer cool efficiently. Figure 5. View largeDownload slide Time evolution of different quantities for our fiducial galaxy model. Top panel: position of the shock for different angles with respect to the disc. The virial and disc scaling radius are indicated by the black dashed lines. Middle panel: velocity for different angles. The wind becomes subsonic in the disc plane at around 0.5 Myr. Lower panel: the pressure of the shocked wind is isotropic and it rises rapidly once the gas can no longer cool efficiently. We further illustrate the outflow dynamics in Fig. 6. Depending on the angle, the outflow has typical velocities of 100–10 000 km s−1 in our fiducial model. The bottom panel nicely shows how the outflow is re-accelerated in the perpendicular direction beyond ∼5 pc due to the exponentially decreasing ISM density. We show the time evolution of the shock structure in Fig. 7. The outflow starts spherically symmetric, but develops an elongated shape, perpendicular to the disc plane after >10 kyr. After ∼1 Myr, the outflow has stopped in the disc plane in our fiducial model and keeps propagating along the poles. Figure 6. View largeDownload slide Velocity maps for our fiducial galaxy model with the galactic disc along the x-axis and the virial radius highlighted in cyan. The linear scaling in the top panel illustrates nicely the overall geometry, and the logarithmic scaling in the lower panel highlights different zones of interest: at small radii, the wind starts with a high velocity, but is decelerated by the increasing mass of the swept-up gas. Whereas the velocity is a decreasing function with radius in the disc plane, the shock accelerates in the direction perpendicular to the disc beyond ∼5 pc due to the strongly decreasing gas density in this direction. Figure 6. View largeDownload slide Velocity maps for our fiducial galaxy model with the galactic disc along the x-axis and the virial radius highlighted in cyan. The linear scaling in the top panel illustrates nicely the overall geometry, and the logarithmic scaling in the lower panel highlights different zones of interest: at small radii, the wind starts with a high velocity, but is decelerated by the increasing mass of the swept-up gas. Whereas the velocity is a decreasing function with radius in the disc plane, the shock accelerates in the direction perpendicular to the disc beyond ∼5 pc due to the strongly decreasing gas density in this direction. Figure 7. View largeDownload slide Time evolution of the shock front with the initial density distribution in the background. The shock front starts spherically symmetric but it quickly develops an elongated shape in the direction of least resistance, i.e. perpendicular to the disc plane. The virial radius for the galaxy is ∼1 kpc. Figure 7. View largeDownload slide Time evolution of the shock front with the initial density distribution in the background. The shock front starts spherically symmetric but it quickly develops an elongated shape in the direction of least resistance, i.e. perpendicular to the disc plane. The virial radius for the galaxy is ∼1 kpc. To analyse the influence of the individual driving and restoring forces, we plot the different contributions to the effective acceleration on a shell element in Fig. 8. The main forces are the acceleration by the AGN (either L/c or the pressure of the adiabatically expanding shell) and the inertia of the swept-up gas, which has to be accelerated. The gravity from the BH and DM as well as the self-gravity of the swept-up gas are negligibly small on all scales of interest. The total acceleration changes sign several times: within ≲0.01 pc the wind accelerates the shell until the mass of the newly swept-up material becomes significant. Within ≲5 pc, the wind decelerates due to the inertia of the swept-up gas. Perpendicular to the disc, the wind accelerates again at ≳5 pc due to the decreasing gas density and hence smaller gas inertia. In the disc plane, the velocity is too small and the shock becomes subsonic. The deceleration and subsequent acceleration perpendicular to the disc plane around ∼5 pc can also be seen in Fig. 6. Figure 8. View largeDownload slide Accelerations as a function of radius for the different contributions perpendicular to the disc (top) and in the disc plane (bottom) for our fiducial galaxy model. Solid lines are outward and dashed lines are inward contributions. Dotted lines represent accelerations that were calculated a posteriori and are not included self-consistently in the dynamics. The black line represents the total acceleration. The dynamics is dominated by the equilibrium between acceleration due to the AGN and the inertia of the swept-up ISM. Figure 8. View largeDownload slide Accelerations as a function of radius for the different contributions perpendicular to the disc (top) and in the disc plane (bottom) for our fiducial galaxy model. Solid lines are outward and dashed lines are inward contributions. Dotted lines represent accelerations that were calculated a posteriori and are not included self-consistently in the dynamics. The black line represents the total acceleration. The dynamics is dominated by the equilibrium between acceleration due to the AGN and the inertia of the swept-up ISM. The contribution of IR and UV account for absorption in the UV and re-emission in the IR, following Ishibashi & Fabian (2015) with the IR and UV opacities κIR = 5 cm2 g−1 and κUV = 103 cm2 g−1. Even for solar metallicity, these effects are subdominant compared to other contributions, but most importantly they are only relevant in the momentum-driven regime, where the momentum of the photons couples directly to the gas. See also Bieri et al. (2017) for a more detailed discussion of the contribution of different photon groups to AGN feedback. 3.2 Parameter study The advantage of our analytical model is that we can easily explore large sets of parameters to study the dependence of the AGN-driven wind on, e.g. the BH mass, redshift, or the AGN lifetime. In this section, we present different parameter studies and analyse their effect on the dynamics of the outflow. Some of the tested parameter combinations represent extreme scenarios and might be very rare in nature or not appear at all. We will also comment on the likelihood of the presented parameter choices. 3.2.1 AGN lifetime We vary the time for which the AGN shines with a constant luminosity, ton, before we set its luminosity to zero (see Fig. 9). For ton > 0.1 Myr, the dynamics does not change, or phrased differently; for our fiducial galaxy model, it does not make a difference whether the AGN is active for 0.1 or 10 Myr. In addition, this illustrates that the most significant part of the momentum transfer happens in the first ∼105 yr. If the AGN shines only for ton ≲ 103 yr, it cannot eject any gas out of the galaxy. The mass outflow rate increases with the lifetime as the AGN-driven wind sweeps up more ISM and reaches its peak value ∼30 kyr after the AGN starts to shine. At this time, the asymmetry of the outflow starts to develop, the pressure-driving in the disc plane is less efficient because the internal energy can escape more easily perpendicularly to the disc plane, and the shock in the disc plane slows down. For ton > 30 kyr, the mass outflow rate is independent of the lifetime. Typical AGN lifetimes are of the order 0.1–10 Myr with variations of the accretion rate on time-scales as short as 1000 yr (Park & Ricotti 2011; Eilers et al. 2017; Negri & Volonteri 2017; Sugimura et al. 2017). Figure 9. View largeDownload slide Dependence of different quantities on the lifetime of the AGN. The ability of the outflow to remove gas from the galaxy increases with ton up to maximum values above ton ≥ 0.1 Myr. The mass accretion and mass outflows rate are defined in equations (18) and (19). Figure 9. View largeDownload slide Dependence of different quantities on the lifetime of the AGN. The ability of the outflow to remove gas from the galaxy increases with ton up to maximum values above ton ≥ 0.1 Myr. The mass accretion and mass outflows rate are defined in equations (18) and (19). 3.2.2 Redshift In our fiducial model, we focus on low-mass galaxies at z = 6, but our 2D approach can also be extrapolated to galaxies at other redshift and we show the redshift dependence of the AGN-driven outflow in Fig. 10. The AGN-driven outflow is generally less efficient at higher redshift due to the deeper gravitational potential; i.e. a halo at higher redshift has a smaller virial radius for the same Mh, due to the higher mean cosmic density. This trend would be amplified by a smaller gas mass fraction at lower redshift, which we do not include in our model: with cosmic time, gas is converted into stars and the remaining ISM exerts a weaker resistance against the outflow, which therefore can propagate further. Figure 10. View largeDownload slide Dependence of different quantities on redshift. The decreasing efficiency with increasing redshift is due to the deeper gravitational potential at high redshift, because of the smaller virial radius. Figure 10. View largeDownload slide Dependence of different quantities on redshift. The decreasing efficiency with increasing redshift is due to the deeper gravitational potential at high redshift, because of the smaller virial radius. 3.2.3 Initial wind velocity The initial wind velocity is crucial, because it is the least well constrained parameter and directly influences the inner wind velocity, hence the Compton cooling time, and consequently the transition from momentum- to energy-driving (Fig. 11). For $$v$$in ≤ 0.1c, the choice of the initial wind velocity does not affect the long-term evolution and the dynamics on galactic scales. However, for values slightly larger than our fiducial value $$v$$in = 0.1c, the transition to energy-driving in the disc occurs at significantly later times (see Section 3.3 for a more detailed discussion). Figure 11. View largeDownload slide Time evolution of the AGN-driven outflow for different values of the initial wind velocity for the direction θ = 45°. For values above our fiducial parameter of $$v$$in = 0.1c, the transition to energy-driving occurs at much larger radii (red, dark blue). Also, the evolution of the asymmetry of the outflow is remarkable: while the shock is more spherical in the momentum-driven phase, it becomes more asymmetric after transition to energy-driving than the outflow that is energy-driven from early on (see the text). Figure 11. View largeDownload slide Time evolution of the AGN-driven outflow for different values of the initial wind velocity for the direction θ = 45°. For values above our fiducial parameter of $$v$$in = 0.1c, the transition to energy-driving occurs at much larger radii (red, dark blue). Also, the evolution of the asymmetry of the outflow is remarkable: while the shock is more spherical in the momentum-driven phase, it becomes more asymmetric after transition to energy-driving than the outflow that is energy-driven from early on (see the text). We can observe another interesting feature: the later the transition from momentum- to energy-driving occurs, the more asymmetric the outflow will be in the late phase (Rperp/Rdisc at t ≳ 106 yr). Except for the tiny contribution of the shell self-gravity, the momentum-driven expansion in one direction does not depend on the outflow dynamics in another direction. This is in contrast to the energy-driven case, where the accelerating pressure depends on the volume enclosed by the shock front, which directly couples the outflow dynamics in different directions. Consequently, the energy-driven outflow develops an asymmetry once it encounters a path of least resistance, which makes the pressure fall and the acceleration in the disc plane smaller. In contrast, the momentum-driven outflow starts more spherical, even at larger radii, where the asymmetric density profile allows the outflow to propagate into a low-density directions. However, when the outflow now becomes energy-driven, with a higher momentum boost, the pressure does not continue to accelerate the outflow in the disc plane, because there is already a very well-defined path of least resistance in the perpendicular direction. In addition, the volume enclosed by the shock front in the moment of the late transition is smaller than the volume of an outflow with an early transition to energy-driving at the same time. Hence, the pressure rises above the value of the pressure of the outflow with an early transition, compare P(t ≈ 105 yr). At the same time, we can see from the time evolution of the momentum boost that this higher pressure (at ∼105 yr) does not drive the acceleration in all directions to the same degree. It rather boosts the expansion perpendicular to the disc plane, for which the momentum-driven outflow has already paved a path of least resistance. We confirm that the momentum boost depends on the confinement time, i.e. the period over which the outflow is enclosed by gas of roughly the same density. 3.2.4 AGN luminosity We have verified that the outflow dynamics on all relevant scales is determined solely by the AGN luminosity, i.e. the product of Eddington ratio and BH mass and not their individual values (see Fig. 12). A small BH with a high Eddington ratio is therefore equivalent to a more massive BH with a smaller Eddington ratio, keeping in mind that high Eddington ratios are much less abundant then low or intermediate values (Habouzit et al. 2017; Weigel et al. 2017). The indicated free-fall times yield the time that gas needs at least to fall back towards the centre once the AGN is off. It scales linearly with the radius out to which the gas is pushed in the disc plane with tff = 5 Myr(R/R0). In Section 3.3, we discuss the origin of the transition to momentum-driving at higher AGN luminosities. Figure 12. View largeDownload slide We compare the outflow, quantified by tff in the disc plane, for different combinations of the Eddington ratio and the BH mass for two halo masses: 108 M⊙ (top panel) and 1010 M⊙ (bottom panel). The blue vertical lines and hashed regions indicate the typical MBH for these halo masses, based on the M–σ relation (Gültekin et al. 2009). The outflow dynamics depends only on the product fEddMBH and we clearly see two ranges of energy- and momentum-driven outflows, separated by a luminosity of ∼1043 erg s−1, as we will discuss in Section 3.3. The same AGN in a two orders of magnitude more massive halo generates fallback times of the gas that are about one order of magnitude shorter. Figure 12. View largeDownload slide We compare the outflow, quantified by tff in the disc plane, for different combinations of the Eddington ratio and the BH mass for two halo masses: 108 M⊙ (top panel) and 1010 M⊙ (bottom panel). The blue vertical lines and hashed regions indicate the typical MBH for these halo masses, based on the M–σ relation (Gültekin et al. 2009). The outflow dynamics depends only on the product fEddMBH and we clearly see two ranges of energy- and momentum-driven outflows, separated by a luminosity of ∼1043 erg s−1, as we will discuss in Section 3.3. The same AGN in a two orders of magnitude more massive halo generates fallback times of the gas that are about one order of magnitude shorter. We compare the models for different BH masses, at the same halo mass, fixing the Eddington ratio at fEdd = 0.3, in Figs 13 and 14. We also observe an interesting trend for the momentum boost, where we observe peak values of ≳100, with the lowest mass black BH generating the highest momentum boosts. Phrased differently, a lower mass BH is more efficient in converting its momentum and energy input into momentum of the outflow. For all considered BH masses, the AGN can drive outflows perpendicular to the disc beyond Rvir, but only for MBH ≳ 105 M⊙ the AGN can eject gas at θ ≥ 45°, whereas the wind powered by less massive BHs is not strong enough and it becomes subsonic before reaching Rvir in the diagonal direction. Figure 13. View largeDownload slide Time evolution of the AGN-driven outflow for different BH masses, i.e. different AGN luminosities for a fixed fEdd = 0.3 with Mh = 108 M⊙. The models with MBH ≤ 105 M⊙ experience an early transition from momentum- to energy-driving (see sudden jump in the pressure) and the model with MBH = 106 M⊙ experiences a later transition (see also Section 3.3 for more details). The radius and velocity represent the evolution at 45°. Figure 13. View largeDownload slide Time evolution of the AGN-driven outflow for different BH masses, i.e. different AGN luminosities for a fixed fEdd = 0.3 with Mh = 108 M⊙. The models with MBH ≤ 105 M⊙ experience an early transition from momentum- to energy-driving (see sudden jump in the pressure) and the model with MBH = 106 M⊙ experiences a later transition (see also Section 3.3 for more details). The radius and velocity represent the evolution at 45°. Figure 14. View largeDownload slide Dependence of the gas ejection efficiency and the mass outflow rate on the BH mass for Mh = 108 M⊙ and fEdd = 0.3. The ability to eject gas from the galaxy has a local peak at MBH ≈ 105 M⊙, because for higher BH masses, i.e. AGN luminosities, the outflow is momentum-driven in the disc (see Section 3.3). Figure 14. View largeDownload slide Dependence of the gas ejection efficiency and the mass outflow rate on the BH mass for Mh = 108 M⊙ and fEdd = 0.3. The ability to eject gas from the galaxy has a local peak at MBH ≈ 105 M⊙, because for higher BH masses, i.e. AGN luminosities, the outflow is momentum-driven in the disc (see Section 3.3). The dependence of the outflow on the Eddington ratio for different BH masses, at fixed halo mass, is shown in Fig. 12. For small BHs (MBH ≤ 103 M⊙), outflows cannot push the gas beyond disc radius even for fEdd = 1. A large range of the parameter space allows for outflows to affect gas in the whole disc (yellow and green symbols in the upper panel of Fig. 12), possibly regulating the time over which an AGN has to remain off before new gas flows in to feed a new accretion episode as well as star formation in the disc. To push gas beyond R0, which is still only $${\sim } 3\,\,\rm{per\,\,cent}$$ of the virial radius, a BH should be as massive as 1 per cent of the halo mass and have fEdd = 1, or 10 per cent of the halo mass and have fEdd > 0.1. Note that these ratios between BH and halo mass are very large. The virial velocity of a halo, i.e. the circular velocity at the virial radius, is given by Vvir = 200 km s−1(Mh/1011 M⊙)1/3 at z = 6 (Ferrarese 2002; Wyithe & Loeb 2002; Volonteri et al. 2011). Assuming that $$\sigma = V_\mathrm{vir}/\sqrt{3}$$ and the z = 0 scaling between BH mass and σ (Gültekin et al. 2009) yields \begin{eqnarray} M_\mathrm{BH} = 6000\,\mathrm{M}_{\odot }\left( \frac{M_h}{10^8\,\mathrm{M}_{\odot }} \right) ^{1.41}, \end{eqnarray} (47)for the typical BH mass in a halo of mass Mh. Simulations suggest that to correctly reproduce the luminosity function of AGN, only a small fraction ($${\lesssim } 25\,\,\rm{per\,\,cent}$$) of high-redshift BHs in galaxies with stellar mass <1010 M⊙ is expected to have Eddington ratios above 0.01, with the majority of BHs accreting at fEdd < 10−4 (Habouzit et al. 2017). For these more realistic parameter combinations, with BHs on the M−σ relation and small Eddington ratios, we expect only short fallback times of the order ≲1 Myr. Therefore, BHs can grow almost continuously and the accretion phases are interrupted by comparably short episodes of AGN feedback: a ∼104 M⊙ BH in a 108 M⊙ halo has a recovery time of ∼1 Myr, indicating that accretion is not interrupted for long. A ∼106 M⊙ BH in a 1010 M⊙ halo creates even shorter recovery times, making BH growth easier in higher mass systems. The role of AGN outflows in driving large amounts of gas out of haloes appears limited, if they are the only source of feedback. A viable possibility is that AGN feedback needs SN feedback as a precursor. Costa et al. (2014) and Prieto et al. (2017) show that AGN feedback is inefficient without the aid of SN feedback: SN winds heat the cold gas in the galaxy, creating a rarefied environment where energy injection from AGN feedback can easily accelerate the gas. In the next section, we discuss how an outflow can still regulate the duty cycle of BHs. In the rest of the paper, we do not vary the BH mass and the Eddington ratio independently, but keep the Eddington ratio at fEdd = 0.3 and vary the BH mass and consequently the AGN luminosity. The results for a different Eddington ratio can be rescaled accordingly. 3.2.5 Halo mass We compare the outflow efficiency for different combinations of the halo and BH mass in Fig 15. As physically expected, the outflow becomes generally more powerful with a higher AGN luminosity and lower halo mass. However, there is a discontinuity of this trend at MBH ≈ 105 M⊙, which we will discuss in more detail in the next section. We further note that contours of equal efficiency are steeper than linear in the log(Mh)–log(MBH) plane. Even for the most extreme scenario of a ∼104 M⊙ BH accreting at 30 per cent Eddington in a 106 M⊙ halo, the outflow hardly pushes the gas beyond R0 in the disc mid-plane. This indicates that the gas recovery times after AGN-driven winds are very short, of the order ∼1 Myr. Figure 15. View largeDownload slide Efficiency of ejecting gas and preventing further accretion for different combinations of the BH and halo mass, quantified by the maximum radius the outflow reaches in the disc plane. The black dashed lines represent $$M_\mathrm{BH} \propto M_h^{1.6}$$, illustrating thresholds of constant maximal radii. The blue line and shaded area show the expected BH masses and their error margins (Gültekin et al. 2009). Note that this relation is valid for local and more massive galaxies and should rather guide the eye and indicate typical host environments for a given BH mass than be used for a quantitative comparison. Figure 15. View largeDownload slide Efficiency of ejecting gas and preventing further accretion for different combinations of the BH and halo mass, quantified by the maximum radius the outflow reaches in the disc plane. The black dashed lines represent $$M_\mathrm{BH} \propto M_h^{1.6}$$, illustrating thresholds of constant maximal radii. The blue line and shaded area show the expected BH masses and their error margins (Gültekin et al. 2009). Note that this relation is valid for local and more massive galaxies and should rather guide the eye and indicate typical host environments for a given BH mass than be used for a quantitative comparison. 3.3 Momentum-driven outflows for AGN luminosities > 1043 erg s−1 We see in Fig. 16 that the transition from an initially momentum-driven to an energy-driven outflow occurs at different radii or not at all depending on the halo and BH mass for a fixed Eddington ratio. The nature of this transition can be better seen in Fig. 17, where we show the time evolution of the pressure in a halo with Mh = 2 × 109 M⊙ and an outflow powered by BHs of different mass. For MBH ≤ 105 M⊙ the transition to energy-driving occurs early (note that this is not our fiducial galaxy model, for which the transition occurs already at t < 0.1 yr). For BH masses of the order 105–106 M⊙ the outflow is always momentum-driven and for even higher BH masses of ≳106 M⊙ we see a late transition to energy-driving (small peak in the red curve around t ≈ 106 yr). This trend can also be seen from a comparison of the relevant time-scales in Fig. 18. The three different regimes can be understood as follows: for a low BH mass (i.e. AGN luminosity), we see an early transition to energy-driving, because of the initial smaller acceleration by the AGN (green, light-blue). For significantly higher AGN luminosities (dark-blue, purple), the outflows remain momentum-driven and reaccelerate beyond ∼5 pc (compare Fig. 6), which decrease the flow-time and make the outflow energy-driven. The intermediate cases (orange, red) also remain momentum-driven, but become subsonic before they reach ∼5 pc. Figure 16. View largeDownload slide Parameter study of the transition radius in the disc from momentum- to energy-driven outflows at $$30\,\,\rm{per\,\,cent}$$ Eddington. We can identify three different regimes: the transition occurs at very small radii close to Rmin (green), the transition occurs at larger radii ≳R0 (yellow), or there is no transition at all and the outflow is always momentum-driven (orange). The last case also includes scenarios where the transition occurs at >Rvir, which is irrelevant for our study. Figure 16. View largeDownload slide Parameter study of the transition radius in the disc from momentum- to energy-driven outflows at $$30\,\,\rm{per\,\,cent}$$ Eddington. We can identify three different regimes: the transition occurs at very small radii close to Rmin (green), the transition occurs at larger radii ≳R0 (yellow), or there is no transition at all and the outflow is always momentum-driven (orange). The last case also includes scenarios where the transition occurs at >Rvir, which is irrelevant for our study. Figure 17. View largeDownload slide Time evolution of the pressure in a halo with a mass of 2 × 109 M⊙ for BH masses spanning two orders of magnitude. This represents a vertical slice in Fig. 16. The transition to an energy-driven outflow occurs at different times or not at all, depending on the BH mass, i.e. on the AGN luminosity. Figure 17. View largeDownload slide Time evolution of the pressure in a halo with a mass of 2 × 109 M⊙ for BH masses spanning two orders of magnitude. This represents a vertical slice in Fig. 16. The transition to an energy-driven outflow occurs at different times or not at all, depending on the BH mass, i.e. on the AGN luminosity. Figure 18. View largeDownload slide For a fixed halo mass of Mh = 2 × 109 M⊙, we plot the ratio of the flow and the cooling time of the AGN-driven outflow as a function of radius (position of the shell) for different BH masses in the disc plane. All outflows start momentum-driven, but if the flow time becomes shorter than the cooling time, Compton cooling can no longer cool the shock efficiently and the shock becomes energy-driven. However, above MBH ≳ 105 M⊙ the cooling time is always shorter than the flow time in the disc and both have the same scaling with radius of t ∝ r2. For even higher AGN luminosities, a late transition to energy-driving occurs on galactic scales. Figure 18. View largeDownload slide For a fixed halo mass of Mh = 2 × 109 M⊙, we plot the ratio of the flow and the cooling time of the AGN-driven outflow as a function of radius (position of the shell) for different BH masses in the disc plane. All outflows start momentum-driven, but if the flow time becomes shorter than the cooling time, Compton cooling can no longer cool the shock efficiently and the shock becomes energy-driven. However, above MBH ≳ 105 M⊙ the cooling time is always shorter than the flow time in the disc and both have the same scaling with radius of t ∝ r2. For even higher AGN luminosities, a late transition to energy-driving occurs on galactic scales. We can understand these transitions at systematically different radii based on simple analytical arguments. The Compton cooling time as a function of radius r is given by equation (41) or for typical parameters: \begin{eqnarray} t_\mathrm{compton} &=& 1.4\, \mathrm{kyr} \left( \frac{f_\mathrm{Edd}}{0.3} \right) ^{-1} \left( \frac{M_\mathrm{BH}}{10^6\,\mathrm{M}_{\odot }} \right) ^{-1} \left( \frac{v_\mathrm{in}}{0.1 c} \right) ^{-2} \left( \frac{r}{1\, \mathrm{pc}} \right) ^{2}. \nonumber \\ \end{eqnarray} (48)To derive the relevant flow time, we have to make several simplifying assumptions. If the transition to energy-driving occurs, it occurs at <1 pc (compare Fig. 18), which is within the sphere of influence and the density in this inner part can be assumed to be constant. Furthermore, we neglect the gravity of the BH. This assumption might seem to be in contradiction to being within the BH's sphere of influence, but the gravitational force is proportional to the mass of the swept-up ISM, which is negligibly small at the radii of interest (compare Fig. 8). The equation of motion is then given by \begin{eqnarray} \frac{\mathrm{d}}{\mathrm{d}t}(M_\mathrm{shell}v)=\frac{L}{c}, \end{eqnarray} (49)with the shell's velocity $$v=\dot{r}$$ and the shell mass $$M_\mathrm{shell} = 4/3 \pi r^3 \rho _0$$. The central gas density ρ0 can be derived from the exponential disc profile as \begin{eqnarray} \rho _0 = \frac{\Sigma _0}{2 H_0} = \frac{G}{2\pi c_s^2}\frac{(m_d M_\mathrm{halo})^2}{(\lambda R_\mathrm{vir})^4}, \end{eqnarray} (50)where H0 = H(R = 0) is the disc scale height in the centre. The differential equation can be solved by a function of the form \begin{eqnarray*} r(t)&=&r_0 \left( \frac{t}{t_0} \right) ^{0.5},\\ v(t)&=&\frac{r_0}{2} (t_0t)^{-0.5}, \end{eqnarray*} which yields directly the same scaling with radius for the flow time as for the cooling time: \begin{eqnarray} t_\mathrm{flow} = \frac{r}{v} = 2 t \propto r^2. \end{eqnarray} (51)We note again that this result is only valid for the inner part of the disc, where the density can be assumed to be constant and the mass of the swept-up gas is still negligible for the gravity. For the flow time as a function of radius, we obtain \begin{eqnarray} t_\mathrm{flow} = \sqrt{\frac{8 \pi c \rho _0}{3L}} r^2 \propto m_d \lambda ^{-2} (1+z)^{3/2} f_\mathrm{Edd}^{-1/2} M_\mathrm{BH}^{-1/2} R^2. \end{eqnarray} (52)Or expressed with fiducial values and fixed md = λ = 0.05 at z = 6: \begin{eqnarray} t_\mathrm{flow} = 2.2\, \mathrm{kyr} \left( \frac{f_\mathrm{Edd}}{0.3} \right)^{-1/2} \left( \frac{M_\mathrm{BH}}{10^6\,\mathrm{M}_{\odot }} \right)^{-1/2} \left( \frac{r}{1\, \mathrm{pc}} \right)^{2}. \end{eqnarray} (53)Setting it equal to the cooling time, we can conclude that the outflow is always momentum-driven if the BH and Eddington ratio fulfil \begin{eqnarray} \frac{f_\mathrm{Edd}}{0.3} \frac{M_\mathrm{BH,crit}}{4 \times 10^5\,\mathrm{M}_{\odot }} > \left( \frac{v_\mathrm{in}}{0.1 c} \right) ^{-4} \left( \frac{1+z}{7} \right) ^{-3}, \end{eqnarray} (54)or expressed via the AGN luminosity \begin{eqnarray} L_\mathrm{AGN} > 10^{43}\, \mathrm{erg}\, \mathrm{s}^{-1} \left( \frac{v_\mathrm{in}}{0.1 c} \right) ^{-4} \left( \frac{1+z}{7} \right) ^{-3}. \end{eqnarray} (55)Note the independence on the halo mass and the strong dependence on the not well-constrained initial velocity of the inner wind $$v$$in. Uncertainties in this parameter can change the resulting transition luminosity by up to three orders of magnitude. The analytically derived BH transition mass, 4 × 105 M⊙, is close to the one obtained in the 2D simulation of ∼2 × 105 M⊙. The small difference is related to the necessary approximations in order to analytically solve the equation of motion. The strong redshift dependence indicates that the transition mass above which a BH can no longer drive an efficient energy-driven wind in the disc is higher at low redshift, caused by the redshift dependence of the central density. The parameter dependences can be explained as follows: a higher Eddington ratio fEdd yields more photons by the AGN, which boost the cooling via inverse Compton scattering. At lower redshift, the gravitational potential is shallower for the same halo mass, which makes the shell velocity higher, the flow time smaller, and hence requires more efficient cooling to sustain a momentum-driven outflow. The initial wind velocity is directly proportional to the post-shock temperature (equation 28) and hence defines the internal energy of the shocked wind. Momentum-driven outflows around the threshold luminosity LAGN ≈ 1043 erg s−1 or slightly above generate a smaller momentum boost (Fig. 13), are less powerful in ejecting gas out of the galaxy (Fig. 14), and drive the gas in the disc plane to smaller maximal radii (Figs 12 and 15). More luminous AGNs with MBH ≈ 107 M⊙ and fEdd = 1.0 in 108 M⊙ haloes are the only systems that create recovery times of the gas in excess of 5 Myr (Fig. 12) via momentum-driving. However, such a combination of Eddington ratio, BH, and halo mass is very unlikely. 3.4 Gas ejection perpendicular to the disc So far, we have focused on the gas dynamics in the disc plane and have seen that AGN-driven winds in 2D push the gas in the disc plane at most to about the scale radius R0. If we instead focus on the gas dynamics and outflow velocities at higher inclination, we can study the wind on galactic scales and relate it to observations of gas dynamics in high-redshift galaxies. In Fig. 19, we compare the outflow velocities at the virial radius for different angles with respect to the disc. For better comparison, we normalize the velocities to the virial velocity of the corresponding haloes, \begin{eqnarray} V_\mathrm{vir} = 20\, \mathrm{km}\, \mathrm{s}^{-1} \left( \frac{M_h}{10^8\,\mathrm{M}_{\odot }} \right)^{1/3}. \end{eqnarray} (56)A low halo mass and a high AGN luminosity generate outflows with higher velocities, as naively expected. As before, this general trend is interrupted by the transition between energy- and momentum-driving around a BH mass of 105 M⊙. Above this threshold, the outflow remains energy-driven at high inclination, but turns momentum-driven in the disc and a part of the post-shock thermal energy is radiated away. This reduces the pressure and therefore the adiabatic acceleration in all directions. Figure 19. View largeDownload slide Outflow velocities at the virial radius, normalized to the virial velocity for different angles with respect to the disc normal. All models are with fEdd = 0.3 and the general trend is that more massive BHs in less massive galaxies create higher outflow velocities, close to the disc normal even in excess of 104 km s−1. However, there is also a wide range of models that cannot eject gas at θ > 45° beyond the virial radius. Figure 19. View largeDownload slide Outflow velocities at the virial radius, normalized to the virial velocity for different angles with respect to the disc normal. All models are with fEdd = 0.3 and the general trend is that more massive BHs in less massive galaxies create higher outflow velocities, close to the disc normal even in excess of 104 km s−1. However, there is also a wide range of models that cannot eject gas at θ > 45° beyond the virial radius. There is a large parameter range, for which gas at θ ≥ 45° cannot escape the halo (red diamonds). In contrast, close to the disc normal, gas is energy-driven and ejected with high velocities of >104 km s−1. AGN-driven outflows eject gas out of the galaxy, and the direction-dependent velocities at Rvir range from ≲100 km s−1 to mildly relativistic velocities close to the disc normal in agreement with observations (Chartas et al. 2002; Pounds et al. 2003; Cappi 2006; Feruglio et al. 2010; Rupke & Veilleux 2011; Aalto et al. 2012; Gofford et al. 2013; Cicone et al. 2014). The subset of models that cannot even eject gas at θ = 10° (red diamonds in the right panel) are the models where the outflow is always momentum-driven in the disc (compare Fig. 16). However, one should keep in mind that these results were derived for an Eddington ratio of $$30\,\,\rm{per\,\,cent}$$ and higher values increase gas ejection (only the product of Eddington ratio and BH mass enters in the equation of motion and results can be rescaled accordingly). 3.5 Comparison to spherical case In this section, we highlight the importance of a 2D approach by comparing it directly to the solutions obtained in a 1D model. We assume that the gas in the 1D model follows the distribution of the DM, which we describe by a Hernquist profile (equation 14). We first compare the general dynamics in our fiducial model in Fig. 20. The 1D outflow is momentum-driven for a longer time and therefore starts with smaller velocities. However, once the 1D outflow becomes energy-driven, it ejects all the gas out of the virial radius, whereas the 2D model can only eject a smaller fraction of the total gas. We further quantify this effect in Fig. 21. In a spherical model, the ejected mass is either 0 per cent for ton ≤ 105 yr or 100 per cent for ton > 105 yr, whereas in our 2D model the ejected mass increases gradually with ton and reaches saturation at $$11\,\,\rm{per\,\,cent}$$ for ton ≳ 105 yr. In this regime, the 1D model overestimates the efficiency of an AGN-driven outflow by one order of magnitude in our fiducial galaxy model. It is difficult to further extrapolate these results to 3D. If we assume the same disc-like geometry in 3D, the ejection fractions and fallback times may not change much. However, a 3D model allows for a more accurate treatment of the inhomogeneous ISM, which could increases or decreases the overall outflow efficiency, as we discuss in Section 4.5. In Fig. 22, we plot the efficiency of the 1D model in ejecting gas out to a certain radius. To quantify the efficiency of the outflow in the 1D Hernquist profile, we apply the same scaling radius, R0, as for 2D exponential disc (equation 5) to check if the outflow can reach at least this radius. The enclosed gas mass within R0 in both scenarios is roughly the same with less than a factor 2 difference. There is a large parameter range for which 1D models can eject all the gas out of a galaxy. This range is of special interest, because it encloses halo masses of ∼107–109 M⊙, which are typical masses of dwarf galaxies at low-redshift and of high-redshift galaxies, which might host the first massive BH seeds. Above a certain halo mass of ∼109 M⊙, which corresponds to a typical BH mass of 105 M⊙, the AGN (at 30 per cent Eddington) is not powerful enough to eject gas out of the galaxy. Figure 20. View largeDownload slide Same as Fig. 5, but with the additional solution of a 1D model (orange, dotted), where we distribute the disc gas mass following a Hernquist profile. The transition to momentum-driving occurs at later time and larger radius, which is partially caused by the different radial density profiles for the Hernquist profile in 1D and an exponential disc in 2D. In the spherical model, all gas is ejected beyond Rvir. Figure 20. View largeDownload slide Same as Fig. 5, but with the additional solution of a 1D model (orange, dotted), where we distribute the disc gas mass following a Hernquist profile. The transition to momentum-driving occurs at later time and larger radius, which is partially caused by the different radial density profiles for the Hernquist profile in 1D and an exponential disc in 2D. In the spherical model, all gas is ejected beyond Rvir. Figure 21. View largeDownload slide Fraction of the disc mass that is ejected by the AGN-driven wind as a function of the time over which the AGN shines with a constant luminosity (fiducial model with fEdd = 0.3). In the 1D model either no or all the gas is ejected out of the virial radius and this transition occurs at ton ≈ 105 yr. In our 2D model, we see a gradual rise of the ejected mass until it reaches the final value of $${\sim } 10\,\,\rm{per\,\,cent}$$ at ton ≳ 105 yr. This demonstrates the strength and importance of a 2D treatment of AGN-driven feedback, compared to a simplistic 1D model. Figure 21. View largeDownload slide Fraction of the disc mass that is ejected by the AGN-driven wind as a function of the time over which the AGN shines with a constant luminosity (fiducial model with fEdd = 0.3). In the 1D model either no or all the gas is ejected out of the virial radius and this transition occurs at ton ≈ 105 yr. In our 2D model, we see a gradual rise of the ejected mass until it reaches the final value of $${\sim } 10\,\,\rm{per\,\,cent}$$ at ton ≳ 105 yr. This demonstrates the strength and importance of a 2D treatment of AGN-driven feedback, compared to a simplistic 1D model. Figure 22. View largeDownload slide Maximum radius the AGN-driven wind reaches in spherical symmetry for different combinations of the BH and halo mass (compare to the 2D results in Fig. 15). The purple pentagons in this plot indicate complete gas ejection out of the galaxy. Figure 22. View largeDownload slide Maximum radius the AGN-driven wind reaches in spherical symmetry for different combinations of the BH and halo mass (compare to the 2D results in Fig. 15). The purple pentagons in this plot indicate complete gas ejection out of the galaxy. We compare these results to Park et al. (2016) and Pacucci, Volonteri & Ferrara (2015), who use 1D radiation-hydrodynamic simulations to study the effect of radiative feedback on the growth of BHs in small haloes. They show that light BHs of ∼100 M⊙ cannot grow efficiently via radiation-regulated accretion. Pacucci et al. (2015) further derive an analytical expression for this critical BH mass, which is in the range of 10–106 M⊙, depending on the accretion scenario and the host properties. In our 1D comparison model, we find a comparable mass range of MBH ≈ 102–107 M⊙, where the BH can completely eject the gas out of a host halo (Mh ≲ 109 M⊙) and hence prevent rapid mass growth of the BH. For a fixed halo mass and Eddington ratio in two dimensions, a higher mass BH generates stronger outflows and ejects the gas to larger radii. This makes lower mass BHs less efficient in stopping their own gas supply and hence more susceptible for mass growth at high duty cycles. This is the inverse trend of Pacucci et al. (2015) and Park et al. (2016), who find low-mass BHs to be more efficient in stopping their own gas supply. These differences are likely caused by different assumptions on what drives the outflow: Pacucci et al. (2015) and Park et al. (2016) study the direct influence of radiation pressure from the BH on the accretion process, while we model AGN-driven winds. For instance, the critical mass for radiation pressure to significantly change the mass accretion rate (Pacucci et al. 2015) and reduce growth is related to the significance of the Eddington limit. The inflow rate is determined by the halo properties. If a given halo provides the same amount of gas inflow to a small BH or a large BH, it will be the smaller BH that will reach the Eddington limit first, having its growth stunted. We note that the differences between the 1D and 2D model presented in this section might partially depend on the different scaling of the density with radius and hence on the different radii where the outflow becomes energy-driven. However, the qualitative differences between the 1D and 2D model remain still valid beyond these differences in the density profile. 4 DISCUSSION We have developed a new 2D analytical model of AGN-driven outflows and present the results for various galaxies and a range of typical conditions. The additional dimension with respect to previously existing 1D models allows us to account for the direction-dependent density profiles. This has significant consequences on the maximally ejected gas mass fraction, on the momentum boost, and on the condition to stop further accretion on the central BH. 4.1 Advantage of 2D With a more realistic disc profile, we confirm the results of 3D hydrodynamic simulations that the outflow propagates preferentially towards low-density regions. In the momentum-driven regime, this is simply due to the lower column densities perpendicular to the disc and hence the lower gas masses that have to be accelerated. In the energy-driven regime, however, the evolution of the entire confined volume is relevant for the acceleration of the shock. Once the energy-driven shock has created a chimney perpendicular to the disc, the pressure and internal energy can escape and the driving in the disc plane is significantly reduced. Naturally, we find that an AGN drives the gas out to different radii for different angles, before the shock front becomes subsonic or crosses the virial radius. In 1D models, there is only one threshold AGN luminosity above which all the ISM is ejected and below which the AGN wind is not efficient enough. For none of our models, the AGN can eject all the gas out of the galaxy. If we define the criteria to stop further accretion on to the SMBH via the minimum time that the swept-up gas needs at least to fall back to the centre, or equivalently that gas is swept out at least to a given radius in the disc plane, we find that this condition results in a proportionality of $$M_\mathrm{BH} f_\mathrm{Edd} \propto M_h^{1.6}$$ (Fig. 15). The power of this scaling is independent of the exact maximal radius or minimal time chosen. However, the absolute value of this threshold has to be quite small (Rmax ≲ 0.01Rvir), if we want to reproduce the observed M–σ relation with our 2D analytical model of AGN feedback (the blue line in Fig. 15). For $$\sigma \approx V_\mathrm{vir} \propto M_\mathrm{h}^{1/3}$$ (Ferrarese 2002; Loeb 2010; Volonteri et al. 2011), we find a scaling between the BH mass and the halo velocity dispersion of \begin{eqnarray} M_\mathrm{BH} \propto \sigma ^{4.8}, \end{eqnarray} (57)which is close to the proposed slope of MBH ∝ σ5 for an energy-driven outflow in 1D models (Haehnelt et al. 1998; Silk & Rees 1998; Fabian 2012; McQuillin & McLaughlin 2013). We do not include SN feedback in our model, which is also expected to change this relation between the galaxy mass and the efficiency to eject gas. We compare SN-driven outflows to AGN-driven winds in more detail in Dashyan et al. (2018). Moreover, our properties of the galactic disc are intrinsically linked to the DM halo by assuming that a constant fraction of the gas with a certain fraction of the angular momentum settles into the disc (Mo et al. 1998). Therefore, we cannot disentangle the effect of our 2D density distribution and the underlying DM halo on the M–σ scaling with our analytical model. Our new 2D model further allows for new interpretations of existing observations. Cicone et al. (2015) present C ii and FIR continuum observations with the Bure Interferometer of a quasar and its host galaxy at z = 6.4189. The C ii reveals velocities up to 1400 km s−1 and outflow with a complex morphology out to ∼30 kpc. They find no evidence for a regular rotation pattern and identify 48 individual clumps with dynamical times in the range 106.6−8.0 yr. They interpret the spread in outflow times as a non-constant AGN luminosity causing various outbursts. However, it could also be related to projection effects or simply different flow times in different directions with respect to the disc plane (compare Fig. 2). This supports the importance of modelling AGN outflows in 2D to allow for an outflow with different velocities in different directions. 4.2 Outflow efficiency as a function of the AGN luminosity The efficiency of an outflow can be defined in different ways: either by its ability to convert input energy into momentum of the outflow (momentum boost), or by its ability to evacuate the ISM within a certain radius and preventing further gas accretion for a certain time. In this section, we use these two quantifications to compare AGN-driven outflows in the momentum- and in the energy-driven regime. BHs of lower mass or accreting at lower Eddington ratio have a higher momentum boost; that is, they have a stronger push relative to the input momentum. The momentum loading of the outflow increases with decreasing AGN luminosity, because if the luminosity is too high, the outflow rapidly propagates towards lower density regions and paves a path of least resistance perpendicular to the disc plane. Consequently, the shocked wind can adiabatically expand in this direction without transferring its momentum to the denser gas in the disc plane. An outflow driven by a lower luminosity AGN remains confined for longer and therefore has more time to transfer its momentum to the gas. However, the strength of the outflow, is higher the more massive and highly accreting the BH is (Figs 12 and 15). High-luminosity AGN are more efficient in driving an outflow out to large galactic radii and hence suppress their own gas supply for a significant amount of time. For instance, Fig. 12 shows that in a halo of 108 M⊙ an AGN with luminosity below <1040 erg s−1 shuts off accretion for only <1 Myr, an AGN with luminosity 1040 erg s−1 < LAGN ≲ 1043 erg s−1 shuts off accretion for ∼1 Myr, and an AGN with luminosity above >1044 erg s−1 shuts off accretion for >10 Myr. These AGN luminosities of log(LAGN/erg s−1) = 40, 43, 44 correspond to BH masses of roughly 102, 105, 106 M⊙ at fEdd = 1.0. In a more massive halo, these time-scales are even shorter, since they depend on the virial radius. For AGN luminosities ≳1043 erg s−1, the outflow remains momentum-driven in the disc out to galactic radii. For AGN luminosities around this threshold, the momentum-driven outflow is less powerful in ejecting gas out of the galactic potential and in preventing further mass accretion in the disc plane, compared to an energy-driven outflow in the same galaxy. However, at even higher AGN luminosities, the AGN is more efficient and can prevent further gas infall for >10 Myr, whereas the energy-driven outflow can only push gas in the disc plane out to distances, which correspond to a recovery time of ∼1 Myr. The main reason for this difference is that the direct momentum input acts more isotropically and works still in the higher-density disc plane, even if there is already a path of least resistance in the perpendicular direction. In contrast to the energy-driven outflow, where all the pressure can escape perpendicular and does not push the gas in the disc plane to larger radii. 4.3 Driving mechanism For low AGN luminosities, the transition from momentum- to energy-driving in the disc plane occurs on very small, i.e. sub-pc, scales in our models. For higher AGN luminosities beyond LAGN ≳ 1043 erg s−1, we observe that the outflow remains momentum-driven in the disc. In all models, we see a mildly relativistic, energy-driven outflow close to the disc normal. This outflow develops due to the low column density in this direction, which causes very high accelerations. The ejected mass in this cone is negligible compared to the gas mass in the disc. Previous 1D models find that the transition from momentum- to energy-driving occurs on scales of ∼1 kpc (Ciotti & Ostriker 1997; King 2003; King, Zubovas & Power 2011; Zubovas & King 2012), several orders of magnitude larger than what we find in our 2D model. In our 1D comparison model, this transition occurs around 1 pc (Fig. 20), but our halo mass is lower than in other 1D studies. Bourne & Nayakshin (2013) model the expected X-ray signature of inverse Compton cooling in a one temperature medium and with thermally decoupled electrons in the post-shock wind. They argue that current observations of AGN do not show evidence of Compton cooling from a one temperature medium and weak constraints on a possible cooling from a two-temperature medium. This observation supports the theory of energy-driven winds on galactic scales that do not radiate away their thermal internal energy. This suggests that most AGN-driven outflows are energy-driven. We derive a characteristic AGN luminosity above which the wind remains momentum-driven in the galactic disc (equation 54). For our fiducial parameters, this characteristic luminosity corresponds to a BH mass of the order ∼105 M⊙, which is independent of the halo mass, but it strongly depends on the redshift and the initial wind velocity. In this paper, we focused on galactic outflows driven by the inner disc wind of the AGN. However, the ISM could also be accelerated directly by the radiation pressure from the AGN (Ishibashi & Fabian 2014, 2015; Thompson et al. 2015; Costa et al. 2018). Although this mechanism can be efficient close to the SMBH, radiation pressure can accelerate gas only up to an effective optical depth of order unity. For electron scattering this transparency radius, where τ ≈ 1, is (King & Pounds 2014, 2015) \begin{eqnarray} R_\mathrm{tr} \approx 50\, \mathrm{pc} \left( \frac{f_\mathrm{gas}}{0.17} \right) \left( \frac{\sigma }{200\, \mathrm{km}\, \mathrm{s}^{-1}} \right) ^{2}, \end{eqnarray} (58)beyond which direct driving by radiation pressure from the AGN cannot accelerate galactic outflows. The inclusion of dust increases the efficiency of direct radiation pressure on the ISM, due to the higher cross-section of dust and the consequently larger transparency radius of several kpc (Ishibashi & Fabian 2012, 2015; King & Pounds 2015). The momentum input of such a radiation pressure-driven wind can be comparable to that of a momentum-driven wind, but the main difference is the frequency-dependent cross-section of dust (with a peak in the UV). Consequently, not all gas might experience the radiation pressure due to its locally low opacity or due to self-shielding by higher opacity regions. Even hybrid models might be realistic, where the inner wind does not directly emerge from the disc, but is accelerated by radiation pressure. In any case, if radiation pressure on electrons and dust drives galactic gas outflows, the accelerated ISM prefers the path of least resistance, which requires a 2D treatment. 4.4 Inferences on quasar outflows The AGN-driven outflows of our fiducial model eject gas out of the galaxy and the direction-dependent velocities at Rvir range from ≲100 km s−1 to mildly relativistic velocities close to the disc normal. Most observations, however, target higher-luminosity, lower-redshift AGN. Our model is not intended to simulate quasars in high-mass galaxies at relatively low redshift, where the gas content is expected to be lower than the universal baryon fraction and a large stellar component is expected. We can, however, speculate on the qualitative behaviour of the outflows by studying the properties of a luminous AGN with luminosity 1047 erg s−1 in a halo with mass ∼1012.5 M⊙ at z ∼ 2 as a reference case (e.g. Prochaska et al. 2013) with 10 per cent and 100 per cent gas fraction. Also, in these models, velocities range from ≲100 km s−1 to almost 104 km s−1, in agreement with observations. Gas is able to escape from ∼60 per cent to 10 per cent of the total solid angle for 10 per cent and 100 per cent gas fractions, respectively. This is the region in the halo within which the outflow is energy-driven (see Section 3.4 for details). The amount of gas escaping from a given direction is, however, not a monotonic quantity. A small amount of gas is pushed at small opening angles around the normal to the disc, as the subtended gas mass is small; the gas mass increases with increasing distance from the normal, but eventually it drops, when the outflow in that direction is unable to escape, analogously to what we see in Fig. 7 for our fiducial set of parameters. This appears in nice agreement with the observational results of Hennawi & Prochaska (2007) and Prochaska et al. (2013), who find a low incidence of gas absorbers along the line of sight to quasars (‘down-the-barrel’), and they argue that gas in this direction is photoionized, while the incidence of absorbers in the transverse direction is larger towards smaller scales, implying an increasing density and covering fraction of neutral gas towards the quasar. The outflow properties in terms of velocity, asymmetry, and energetics are also in good agreement with observations (Prochaska & Hennawi 2009). One important difference is the gas temperature. The observed gas is cool, ∼104 K, and cannot have been shock heated recently, which disfavours AGN-driven outflows (especially energy-driven ones) as an interpretation. However, depending on the time-scales, the gas could already have cooled down, after being pushed out by the AGN wind. 4.5 Caveats We treat the ISM as homogeneous and assume that the AGN-driven shock sweeps up all the enclosed gas. Real galaxies have a multiphase ISM with H ii regions and dense molecular clumps, which alter this simplistic treatment. Theoretical models of the impact of AGN winds (Bieri et al. 2017) and jets (Wagner, Bicknell & Umemura 2012; Wagner, Umemura & Bicknell 2013; Cielo et al. 2017) on a clumpy ISM show that dense clumps can be dissolved by IR photons that can penetrate into these clumps. The efficiency of the AGN-driven wind depends on the size and filling factor of these dense clumps. A galaxy with many small isolated clouds experiences efficient cloud dispersion compared to a galaxy with fewer but bigger cloud complexes. By neglecting the realistic ISM structure, we overestimate the outflow efficiency, because a fractal ISM predefines already paths of least resistance and the momentum transfer to dense clumps is smaller compared to a homogeneous ISM. On the other hand, the low-density gas between the clumps is accelerated and ejected more easily. The shock can also form dense cold clumps via dynamical instabilities as shown by Costa et al. (2015), which we do not include in our model. The contact discontinuity is strongly Rayleigh-Taylor unstable and might give rise to the formation of clumps that form behind the shock front and decouple from the outflow. These denser entities could also be the constituents of the high-speed molecular outflows. Recently, Ferrara & Scannapieco (2016) model the formation of molecular clumps in AGN-driven outflows with 3D simulations. They find that clumps might form at the transition from momentum- to energy-driving, but they get rapidly dissolved by the hot shock gas flowing past them (see also Bieri et al. 2017; Richings & Faucher-Giguère 2018). Our analytical model focuses on isolated galaxies to better study the influence of individual parameters, independently from the cosmological environment. Therefore, we are not able to take into account large-scale effects, such as the clustering of galaxies (Gaspari et al. 2011a,b, 2014), the accretion of cold gas along cosmic filaments, or galaxy mergers. We assume the AGN luminosity to be constant during the lifetime of the AGN. This is a necessary assumption to clearly disentangle the influence of different input parameters. Different groups study the radiation-regulated accretion on to the BH and demonstrate the importance of multidimensional hydrodynamical simulation to self-consistently follow the accretion and resulting luminosity (Park & Ricotti 2011; Negri & Volonteri 2017; Sugimura et al. 2017). The accretion rate and hence the AGN luminosity oscillate by up to two orders of magnitude on time-scales of several thousand years shorter than the lifetime of the AGN in our model. Gilli et al. (2017) investigate galactic outflows driven by an AGN with exponentially increasing luminosity; that is, they self-consistently account for the mass growth of the central BH at constant Eddington ratio. They find that the late time expansion of the radius in their 1D model is also exponential, irrespective of the details of the driving mechanism. We only distinguish between the momentum- and energy-driven phase, whereas Faucher-Giguère & Quataert (2012) also account for the ‘intermediate partially radiative bubble stage’. In this phase, the cooling time is shorter than the flow time but longer than the crossing time of the shock and the shock cools only partially. Although we agree that a more realistic transition from one driving mechanism to another is necessary, this does not affect our fiducial model since the transition occurs already very early and is very sharp. The weak collisional coupling between protons and electrons in the shocked wind can increase the Compton cooling time by two orders of magnitude (Section 2.4.2). We do not implicitly include this effect, but as pointed out by Faucher-Giguère & Quataert (2012), the more realistic treatment of the shock as a two-temperature medium significantly decreases the efficiency of inverse Compton cooling. This does not change anything in our fiducial model, because the transition to energy-driven occurs already at very small radii. However, the AGN transition luminosity, above which the outflow is always momentum-driven in the disc (equation 54), will be higher as it scales $$L_\mathrm{AGN,crit}\propto t_\mathrm{compton}^2$$. An additional increase of the Compton cooling time can be achieved by the inclusion of Compton heating. Sazonov, Ostriker & Sunyaev (2004) determine the equilibrium Compton temperature of gas exposed to a characteristic AGN spectrum to be TC = 2 × 107 K, which is insensitive to obscuration effects. The equilibrium Compton temperature is equal to the mean photon energy averaged over the AGN spectra. They conclude that the characteristic Compton heating and cooling rates per particle should be the same within a factor of ∼2. Most of the Compton cooling will be provided by the IR component, whereas the Compton heating is dominated by the high-energy component. The maximum radius out to which an AGN can heat low density, highly ionized gas by Compton heating is \begin{eqnarray} r_C = 0.4\, \mathrm{kpc}\left( \frac{f_\mathrm{Edd}}{0.1} \right) ^{1/2} \left( \frac{M_\mathrm{BH}}{10^8 \,\mathrm{M}_{\odot }} \right) ^{1/2}. \end{eqnarray} (59)For ionization parameters below ξ ≲ 105, the gas cannot be heated to TC due to other more efficient cooling mechanisms. A more detailed discussion of the effect of Compton heating is given in Sazonov et al. (2004, 2005). Hydrodynamical simulations of a geometrically thin and optically thick accretion disc show that the disc wind might not be isotropic, as assumed in our model, but has a covering factor of $${\sim } 20\,\,\rm{per\,\,cent}$$ (Proga, Stone & Kallman 2000). Moreover, the initial wind velocity itself might depend on the angle to the disc with higher velocities in the perpendicular direction, and a possible self-shadowing effect by the outer part of the disc amplifies this anisotropy (Sugimura et al. 2017). 5 SUMMARY AND CONCLUSION We present a new 2D analytical model for AGN-driven outflows and demonstrate the importance of a more realistic gaseous disc profile for the outflow dynamics. In contrast to simplistic 1D models, we predict smaller gas ejection fractions and shorter fallback times, both facilitating an efficient mass growth of BHs via feeding from galactic scales. This is related to the energy-driven nature of the outflows, which pave a path of least resistance perpendicular to the disc and hence prevent efficient driving in the disc plane. We confirm earlier results that AGN-driven winds can transport energy and momentum from sub-pc scales around the BH to galactic scales and thereby regulate the co-evolution of a galaxy and its central BH. Our main results are as follows: For typical high-redshift galaxies, the AGN can eject at most $${\sim } 10\,\,\rm{per\,\,cent}$$ of the ISM out of the halo, whereas 1D models predict complete gas ejection under similar conditions. At high redshift, the ejected gas mass fractions are lower due to the deeper gravitational potential, compared to similar galaxies at low redshift. The characteristic time for which the AGN can suppress further gas supply is remarkably short, of the order a few million years. We find AGNs with a low luminosity to be more efficient in converting their input energy into outflow momentum, because the swept-up medium is confined for longer by the shock front and has consequently more time to store up internal energy. We also find a systematic transition in the outflow nature: for AGN luminosities below 1043 erg s−1, the outflow is energy-driven, independent of the halo properties. At higher luminosities, the outflow remains momentum-driven in the disc plane. Independently of the exact criterion to suppress further gas accretion, we find a slope of MBH ∝ σ4.8 for the M–σ relation. Our new model highlights the importance of a realistic 2D density profile to predict the ejected gas mass fraction and fallback time. We can reproduce results of 3D hydrodynamical simulations with much less computational effort, although more realistic models are still necessary to correctly account for the substructure of the ISM and for environmental effects. Our results can be used as an improved sub-grid model in cosmological simulations or semi-analytical models of galaxy formation to predict the efficiency of AGN feedback. ACKNOWLEDGEMENTS We thank the reviewer for her/his suggestions and careful reading of the manuscript. We are grateful to Alex Wagner and Rebekka Bieri for valuable discussions and helpful comments on our model. The authors acknowledge funding under the European Community's Seventh Framework Programme (FP7/2007-2013) via the European Research Council Grants ‘BLACK’ under the project number 614199. REFERENCES Aalto S., Muller S., Sakamoto K., Gallagher J. S., Martín S., Costagliola F., 2012, A&A , 546, A68 CrossRef Search ADS Baldassare V. F., Reines A. E., Gallo E., Greene J. E., 2015, ApJ , 809, L14 https://doi.org/10.1088/2041-8205/809/1/L14 CrossRef Search ADS Barai P., Gallerani S., Pallottini A., Ferrara A., Marconi A., Cicone C., Maiolino R., Carniani S., 2018, MNRAS , 473, 4003 https://doi.org/10.1093/mnras/stx2563 CrossRef Search ADS Behroozi P. S., Wechsler R. H., Conroy C., 2013, ApJ , 770, 57 https://doi.org/10.1088/0004-637X/770/1/57 CrossRef Search ADS Bieri R., Dubois Y., Rosdahl J., Wagner A., Silk J., Mamon G. A., 2017, MNRAS , 464, 1854 https://doi.org/10.1093/mnras/stw2380 CrossRef Search ADS Bourne M. A., Nayakshin S., 2013, MNRAS , 436, 2346 https://doi.org/10.1093/mnras/stt1739 CrossRef Search ADS Bower R. G., Benson A. J., Malbon R., Helly J. C., Frenk C. S., Baugh C. M., Cole S., Lacey C. G., 2006, MNRAS , 370, 645 https://doi.org/10.1111/j.1365-2966.2006.10519.x CrossRef Search ADS Cappi M., 2006, Astron. Nachr. , 327, 1012 https://doi.org/10.1002/asna.200610639 CrossRef Search ADS Chartas G., Brandt W. N., Gallagher S. C., Garmire G. P., 2002, ApJ , 579, 169 https://doi.org/10.1086/342744 CrossRef Search ADS Cicone C. et al. , 2014, A&A , 562, A21 CrossRef Search ADS Cicone C. et al. , 2015, A&A , 574, A14 CrossRef Search ADS Cielo S., Bieri R., Volonteri M., Wagner A., Dubois Y., 2017, preprint (arXiv:1712.03955) Ciotti L., Ostriker J. P., 1997, ApJ , 487, L105 https://doi.org/10.1086/310902 CrossRef Search ADS Costa T., Sijacki D., Trenti M., Haehnelt M. G., 2014, MNRAS , 439, 2146 https://doi.org/10.1093/mnras/stu101 CrossRef Search ADS Costa T., Sijacki D., Haehnelt M. G., 2014, MNRAS , 444, 2355 https://doi.org/10.1093/mnras/stu1632 CrossRef Search ADS Costa T., Sijacki D., Haehnelt M. G., 2015, MNRAS , 448, L30 https://doi.org/10.1093/mnrasl/slu193 CrossRef Search ADS Costa T., Rosdahl J., Sijacki D., Haehnelt M. G., 2018, MNRAS , 473, 4197 https://doi.org/10.1093/mnras/stx2598 CrossRef Search ADS Courteau S., de Jong R. S., Broeils A. H., 1996, ApJ , 457, L73 https://doi.org/10.1086/309906 CrossRef Search ADS Croton D. J. et al. , 2006, MNRAS , 365, 11 https://doi.org/10.1111/j.1365-2966.2005.09675.x CrossRef Search ADS Croton D. J. et al. , 2016, ApJS , 222, 22 https://doi.org/10.3847/0067-0049/222/2/22 CrossRef Search ADS Dashyan G., Silk J., Mamon G. A., Dubois Y., Hartwig T., 2018, MNRAS , 473, 5698 https://doi.org/10.1093/mnras/stx2716 CrossRef Search ADS Decarli R. et al. , 2018, ApJ , 854, 97 Dotti M., Colpi M., Maraschi L., Perego A., Volonteri M., 2010, in Maraschi L., Ghisellini G., Della Ceca R., Tavecchio F., eds, Vol. 427 of ASP Conf. Ser. Vol. 427, Accretion and Ejection in AGN: a Global View . Astron. Soc. Pac., San Francisco, p. 19 Eilers A.-C., Davies F. B., Hennawi J. F., Prochaska J. X., Lukić Z., Mazzucchelli C., 2017, ApJ , 840, 24 https://doi.org/10.3847/1538-4357/aa6c60 CrossRef Search ADS Fabian A. C., 1999, MNRAS , 308, L39 https://doi.org/10.1046/j.1365-8711.1999.03017.x CrossRef Search ADS Fabian A. C., 2012, ARA&A , 50, 455 CrossRef Search ADS Faucher-Giguère C.-A., Quataert E., 2012, MNRAS , 425, 605 https://doi.org/10.1111/j.1365-2966.2012.21512.x CrossRef Search ADS Ferrara A., Scannapieco E., 2016, ApJ , 833, 46 https://doi.org/10.3847/1538-4357/833/1/46 CrossRef Search ADS Ferrarese L., 2002, ApJ , 578, 90 https://doi.org/10.1086/342308 CrossRef Search ADS Ferrarese L., Merritt D., 2000, ApJ , 539, L9 https://doi.org/10.1086/312838 CrossRef Search ADS Feruglio C., Maiolino R., Piconcelli E., Menci N., Aussel H., Lamastra A., Fiore F., 2010, A&A , 518, L155 CrossRef Search ADS Förster Schreiber N. M. et al. , 2009, ApJ , 706, 1364 https://doi.org/10.1088/0004-637X/706/2/1364 CrossRef Search ADS Freeman K. C., 1970, ApJ , 160, 811 https://doi.org/10.1086/150474 CrossRef Search ADS Gabor J. M., Bournaud F., 2013, MNRAS , 434, 606 https://doi.org/10.1093/mnras/stt1046 CrossRef Search ADS Gabor J. M., Bournaud F., 2014, MNRAS , 441, 1615 https://doi.org/10.1093/mnras/stu677 CrossRef Search ADS Gaspari M., Melioli C., Brighenti F., D'Ercole A., 2011a, MNRAS , 411, 349 https://doi.org/10.1111/j.1365-2966.2010.17688.x CrossRef Search ADS Gaspari M., Brighenti F., D'Ercole A., Melioli C., 2011b, MNRAS , 415, 1549 https://doi.org/10.1111/j.1365-2966.2011.18806.x CrossRef Search ADS Gaspari M., Brighenti F., Ruszkowski M., 2013, Astronomische Nachrichten , 334, 394 https://doi.org/10.1002/asna.201211865 CrossRef Search ADS Gaspari M., Brighenti F., Temi P., Ettori S., 2014, ApJ , 783, L10 https://doi.org/10.1088/2041-8205/783/1/L10 CrossRef Search ADS Gebhardt K. et al. , 2000, ApJ , 539, L13 https://doi.org/10.1086/312840 CrossRef Search ADS Genzel R. et al. , 2008, ApJ , 687, 59 https://doi.org/10.1086/591840 CrossRef Search ADS Genzel R. et al. , 2017, Nature , 543, 397 https://doi.org/10.1038/nature21685 CrossRef Search ADS PubMed Gilli R., Calura F., D'Ercole A., Norman C., 2017, A&A , 603, A69 CrossRef Search ADS Gofford J., Reeves J. N., Tombesi F., Braito V., Turner T. J., Miller L., Cappi M., 2013, MNRAS , 430, 60 https://doi.org/10.1093/mnras/sts481 CrossRef Search ADS Goodman J., Tan J. C., 2004, ApJ , 608, 108 https://doi.org/10.1086/386360 CrossRef Search ADS Greene J. E., Ho L. C., Barth A. J., 2008, ApJ , 688, 159 https://doi.org/10.1086/592078 CrossRef Search ADS Gültekin K. et al. , 2009, ApJ , 698, 198 https://doi.org/10.1088/0004-637X/698/1/198 CrossRef Search ADS Habouzit M., Volonteri M., Dubois Y., 2017, MNRAS , 468, 3935 https://doi.org/10.1093/mnras/stx666 CrossRef Search ADS Haehnelt M. G., Natarajan P., Rees M. J., 1998, MNRAS , 300, 817 https://doi.org/10.1111/j.1365-8711.1998.t01-1-01951.x CrossRef Search ADS Häring N., Rix H.-W., 2004, ApJ , 604, L89 https://doi.org/10.1086/383567 CrossRef Search ADS Heckman T. M., Kauffmann G., 2011, Science , 333, 182 https://doi.org/10.1126/science.1200504 CrossRef Search ADS PubMed Hennawi J. F., Prochaska J. X., 2007, ApJ , 655, 735 https://doi.org/10.1086/509770 CrossRef Search ADS Hernquist L., 1990, ApJ , 356, 359 https://doi.org/10.1086/168845 CrossRef Search ADS Hirschmann M., Dolag K., Saro A., Bachmann L., Borgani S., Burkert A., 2014, MNRAS , 442, 2304 https://doi.org/10.1093/mnras/stu1023 CrossRef Search ADS Ishibashi W., Fabian A. C., 2012, MNRAS , 427, 2998 https://doi.org/10.1111/j.1365-2966.2012.22074.x CrossRef Search ADS Ishibashi W., Fabian A. C., 2014, MNRAS , 441, 1474 https://doi.org/10.1093/mnras/stu672 CrossRef Search ADS Ishibashi W., Fabian A. C., 2015, MNRAS , 451, 93 https://doi.org/10.1093/mnras/stv944 CrossRef Search ADS King A., 2003, ApJ , 596, L27 https://doi.org/10.1086/379143 CrossRef Search ADS King A., 2005, ApJ , 635, L121 https://doi.org/10.1086/499430 CrossRef Search ADS King A. R., 2010, MNRAS , 402, 1516 https://doi.org/10.1111/j.1365-2966.2009.16013.x CrossRef Search ADS King A. R., Pounds K. A., 2014, MNRAS , 437, L81 https://doi.org/10.1093/mnrasl/slt144 CrossRef Search ADS King A., Pounds K., 2015, ARA&A , 53, 115 CrossRef Search ADS King A. R., Pringle J. E., Hofmann J. A., 2008, MNRAS , 385, 1621 https://doi.org/10.1111/j.1365-2966.2008.12943.x CrossRef Search ADS King A. R., Zubovas K., Power C., 2011, MNRAS , 415, L6 https://doi.org/10.1111/j.1745-3933.2011.01067.x CrossRef Search ADS Kormendy J., Ho L. C., 2013, ARA&A , 51, 511 CrossRef Search ADS Loeb A., 2010, How Did the First Stars and Galaxies Form? Princeton Univ. Press, Princeton, NJ Google Scholar CrossRef Search ADS Magorrian J. et al. , 1998, AJ , 115, 2285 https://doi.org/10.1086/300353 CrossRef Search ADS Maiolino R. et al. , 2017, Nature , 544, 202 https://doi.org/10.1038/nature21677 CrossRef Search ADS PubMed Marconi A., Hunt L. K., 2003, ApJ , 589, L21 https://doi.org/10.1086/375804 CrossRef Search ADS Martin C. L., 2005, ApJ , 621, 227 https://doi.org/10.1086/427277 CrossRef Search ADS McConnell N. J., Ma C.-P., 2013, ApJ , 764, 184 https://doi.org/10.1088/0004-637X/764/2/184 CrossRef Search ADS McQuillin R. C., McLaughlin D. E., 2012, MNRAS , 423, 2162 https://doi.org/10.1111/j.1365-2966.2012.21028.x CrossRef Search ADS McQuillin R. C., McLaughlin D. E., 2013, MNRAS , 434, 1332 https://doi.org/10.1093/mnras/stt1101 CrossRef Search ADS Mo H. J., Mao S., White S. D. M., 1998, MNRAS , 295, 319 https://doi.org/10.1046/j.1365-8711.1998.01227.x CrossRef Search ADS Murray N., Quataert E., Thompson T. A., 2005, ApJ , 618, 569 https://doi.org/10.1086/426067 CrossRef Search ADS Negri A., Volonteri M., 2017, MNRAS , 467, 3475 https://doi.org/10.1093/mnras/stx362 CrossRef Search ADS O'Shea B. W., Wise J. H., Xu H., Norman M. L., 2015, ApJ , 807, L12 https://doi.org/10.1088/2041-8205/807/1/L12 CrossRef Search ADS Ostriker J. P., Choi E., Ciotti L., Novak G. S., Proga D., 2010, ApJ , 722, 642 https://doi.org/10.1088/0004-637X/722/1/642 CrossRef Search ADS Pacucci F., Volonteri M., Ferrara A., 2015, MNRAS , 452, 1922 https://doi.org/10.1093/mnras/stv1465 CrossRef Search ADS Park K., Ricotti M., 2011, ApJ , 739, 2 https://doi.org/10.1088/0004-637X/739/1/2 CrossRef Search ADS Park K., Ricotti M., Natarajan P., Bogdanović T., Wise J. H., 2016, ApJ , 818, 184 https://doi.org/10.3847/0004-637X/818/2/184 CrossRef Search ADS Pawlik A. H., Milosavljević M., Bromm V., 2011, ApJ , 731, 54 https://doi.org/10.1088/0004-637X/731/1/54 CrossRef Search ADS Penny S. J. et al. , 2018, MNRAS , preprint (arXiv:1710.07568) Pounds K. A., Reeves J. N., King A. R., Page K. L., O'Brien P. T., Turner M. J. L., 2003, MNRAS , 345, 705 https://doi.org/10.1046/j.1365-8711.2003.07006.x CrossRef Search ADS Prieto J., Escala A., Volonteri M., Dubois Y., 2017, ApJ , 836, 216 https://doi.org/10.3847/1538-4357/aa5be5 CrossRef Search ADS Prochaska J. X., Hennawi J. F., 2009, ApJ , 690, 1558 https://doi.org/10.1088/0004-637X/690/2/1558 CrossRef Search ADS Prochaska J. X. et al. , 2013, ApJ , 776, 136 https://doi.org/10.1088/0004-637X/776/2/136 CrossRef Search ADS Proga D., Stone J. M., Kallman T. R., 2000, ApJ , 543, 686 https://doi.org/10.1086/317154 CrossRef Search ADS Reeves J. N., O'Brien P. T., Ward M. J., 2003, ApJ , 593, L65 https://doi.org/10.1086/378218 CrossRef Search ADS Reines A. E., Volonteri M., 2015, ApJ , 813, 82 https://doi.org/10.1088/0004-637X/813/2/82 CrossRef Search ADS Richings A. J., Faucher-Giguère C.-A., 2018, MNRAS , 474, 3673 https://doi.org/10.1093/mnras/stx3014 CrossRef Search ADS Roos O., Juneau S., Bournaud F., Gabor J. M., 2015, ApJ , 800, 19 https://doi.org/10.1088/0004-637X/800/1/19 CrossRef Search ADS Rupke D. S. N., Veilleux S., 2011, ApJ , 729, L27 https://doi.org/10.1088/2041-8205/729/2/L27 CrossRef Search ADS Safranek-Shrader C., Agarwal M., Federrath C., Dubey A., Milosavljević M., Bromm V., 2012, MNRAS , 426, 1159 https://doi.org/10.1111/j.1365-2966.2012.21852.x CrossRef Search ADS Sazonov S. Y., Ostriker J. P., Sunyaev R. A., 2004, MNRAS , 347, 144 https://doi.org/10.1111/j.1365-2966.2004.07184.x CrossRef Search ADS Sazonov S. Y., Ostriker J. P., Ciotti L., Sunyaev R. A., 2005, MNRAS , 358, 168 https://doi.org/10.1111/j.1365-2966.2005.08763.x CrossRef Search ADS Sérsic J. L., 1963, Boletin de la Asociacion Argentina de Astronomia La Plata Argentina , 6, 41 Shakura N. I., Sunyaev R. A., 1973, A&A , 24, 337 Sharma M., Nath B. B., 2013, ApJ , 763, 17 https://doi.org/10.1088/0004-637X/763/1/17 CrossRef Search ADS Shull J. M., Draine B. T., 1987, in Hollenbach D. J., Thronson H. A. Jr, eds, Interstellar Processes Vol. 134 of Astrophysics and Space Science Library , The Physics of Interstellar Shock Waves. Springer-Verlag, Heidelberg, p. 283 Silk J., Nusser A., 2010, ApJ , 725, 556 https://doi.org/10.1088/0004-637X/725/1/556 CrossRef Search ADS Silk J., Rees M. J., 1998, A&A , 331, L1 Smit R. et al. , 2018, Nature , 553, 178 CrossRef Search ADS PubMed Somerville R. S., Hopkins P. F., Cox T. J., Robertson B. E., Hernquist L., 2008, MNRAS , 391, 481 https://doi.org/10.1111/j.1365-2966.2008.13805.x CrossRef Search ADS Sturm E. et al. , 2011, ApJ , 733, L16 https://doi.org/10.1088/2041-8205/733/1/L16 CrossRef Search ADS Sugimura K., Hosokawa T., Yajima H., Omukai K., 2017, MNRAS , 469, 62 https://doi.org/10.1093/mnras/stx769 CrossRef Search ADS Sutherland R. S., Dopita M. A., 1993, ApJS , 88, 253 https://doi.org/10.1086/191823 CrossRef Search ADS Thompson T. A., Fabian A. C., Quataert E., Murray N., 2015, MNRAS , 449, 147 https://doi.org/10.1093/mnras/stv246 CrossRef Search ADS Tombesi F., 2016, Astron. Nachr. , 337, 410 https://doi.org/10.1002/asna.201612322 CrossRef Search ADS Tombesi F., Cappi M., Reeves J. N., Palumbo G. G. C., Yaqoob T., Braito V., Dadina M., 2010, A&A , 521, A57 CrossRef Search ADS Tozzi P., Norman C., 2001, ApJ , 546, 63 https://doi.org/10.1086/318237 CrossRef Search ADS Trakhtenbrot B., Lira P., Netzer H., Cicone C., Maiolino R., Shemmer O., 2017, ApJ , 836, 8 https://doi.org/10.3847/1538-4357/836/1/8 CrossRef Search ADS Trebitsch M., Volonteri M., Dubois Y., Madau P., 2017, MNRAS , preprint (arXiv:1712.05804) Tremaine S. et al. , 2002, ApJ , 574, 740 https://doi.org/10.1086/341002 CrossRef Search ADS Volonteri M., Haardt F., Ghisellini G., Della Ceca R., 2011, MNRAS , 416, 216 Wagner A. Y., Bicknell G. V., Umemura M., 2012, ApJ , 757, 136 https://doi.org/10.1088/0004-637X/757/2/136 CrossRef Search ADS Wagner A. Y., Umemura M., Bicknell G. V., 2013, ApJ , 763, L18 https://doi.org/10.1088/2041-8205/763/1/L18 CrossRef Search ADS Wang R. et al. , 2013, ApJ , 773, 44 https://doi.org/10.1088/0004-637X/773/1/44 CrossRef Search ADS Weigel A. K., Schawinski K., Caplar N., Wong O. I., Treister E., Trakhtenbrot B., 2017, ApJ , 845, 134 https://doi.org/10.3847/1538-4357/aa803b CrossRef Search ADS Wyithe J. S. B., Loeb A., 2002, ApJ , 581, 886 https://doi.org/10.1086/344249 CrossRef Search ADS Zubovas K., King A., 2012, ApJ , 745, L34 https://doi.org/10.1088/2041-8205/745/2/L34 CrossRef Search ADS © 2018 The Author(s) Published by Oxford University Press on behalf of the Royal Astronomical Society

Monthly Notices of the Royal Astronomical Society – Oxford University Press

**Published: ** May 1, 2018

Loading...

personal research library

It’s your single place to instantly

**discover** and **read** the research

that matters to you.

Enjoy **affordable access** to

over 18 million articles from more than

**15,000 peer-reviewed journals**.

All for just $49/month

Query the DeepDyve database, plus search all of PubMed and Google Scholar seamlessly

Save any article or search result from DeepDyve, PubMed, and Google Scholar... all in one place.

Get unlimited, online access to over 18 million full-text articles from more than 15,000 scientific journals.

Read from thousands of the leading scholarly journals from *SpringerNature*, *Elsevier*, *Wiley-Blackwell*, *Oxford University Press* and more.

All the latest content is available, no embargo periods.

## “Hi guys, I cannot tell you how much I love this resource. Incredible. I really believe you've hit the nail on the head with this site in regards to solving the research-purchase issue.”

Daniel C.

## “Whoa! It’s like Spotify but for academic articles.”

@Phil_Robichaud

## “I must say, @deepdyve is a fabulous solution to the independent researcher's problem of #access to #information.”

@deepthiw

## “My last article couldn't be possible without the platform @deepdyve that makes journal papers cheaper.”

@JoseServera

DeepDyve ## Freelancer | DeepDyve ## Pro | |
---|---|---|

Price | FREE | $49/month |

Save searches from | ||

Create lists to | ||

Export lists, citations | ||

Read DeepDyve articles | Abstract access only | Unlimited access to over |

20 pages / month | ||

PDF Discount | 20% off | |

Read and print from thousands of top scholarly journals.

System error. Please try again!

or

By signing up, you agree to DeepDyve’s Terms of Service and Privacy Policy.

Already have an account? Log in

Bookmark this article. You can see your Bookmarks on your DeepDyve Library.

To save an article, **log in** first, or **sign up** for a DeepDyve account if you don’t already have one.

All DeepDyve websites use cookies to improve your online experience. They were placed on your computer when you launched this website. You can change your cookie settings through your browser.

ok to continue