# Momentum-driven Winds from Radiatively Efficient Black Hole Accretion and Their Impact on Galaxies

Ryan Brennan, Ena Choi, Rachel S. Somerville, Michaela Hirschmann,

Thorsten Naab, Jeremiah P. Ostriker

► To cite this version:

Ryan Brennan, Ena Choi, Rachel S. Somerville, Michaela Hirschmann, Thorsten Naab, et al. Momentum-driven Winds from Radiatively Efficient Black Hole Accretion and Their Impact on Galaxies. The Astrophysical Journal, 2018, 860, 10.3847/1538-4357/aac2c4. insu-03747697

HAL Id: insu-03747697

https://insu.hal.science/insu-03747697v1

Submitted on 20 Aug 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_2.png" /></figure>

# Momentum-driven Winds from Radiatively Efficient Black Hole Accretion and Their Impact on Galaxies

Ryan Brennan¹, Ena Choi¹,²,<fnref n="7" /> Ⓓ , Rachel S. Somerville¹,³, Michaela Hirschmann⁴ Ⓓ , Thorsten Naab⁵ Ⓓ , and Jeremiah P. Ostriker²,⁶ Ⓓ

¹Department of Physics and Astronomy, Rutgers, The State University of New Jersey, 136 Frelinghuysen Road, Piscataway, NJ 08854, USA
ena.choi@columbia.edu

²Department of Astronomy, Columbia University, 550 W 120th Street, New York, NY 10027, USA

³Center for Computational Astrophysics, Flatiron Institute, 162 5th Avenue, New York, NY 10010, USA

⁴Sorbonne Universités, UPMC-CNRS, UMR7095, Institut d'Astrophysique de Paris, F-75014 Paris, France

⁵Max-Planck-Institute for Astrophysics, Karl-Schwarzschild-Str. 1, D-85748 Garching, Germany

⁶Department of Astrophysical Sciences, Princeton University, Princeton, NJ 08544, USA

Received 2018 January 29; revised 2018 April 23; accepted 2018 May 2; published 2018 June 7

## Abstract

We explore the effect of momentum-driven winds representing radiation-pressure-driven outflows from accretion onto supermassive black holes in a set of numerical hydrodynamical simulations. We explore two matched sets of cosmological zoom-in runs of 24 halos with masses ~10¹².⁰–10¹³.⁴ $M_\odot$ run with two different feedback models. Our “NoAGN” model includes stellar feedback via UV heating, stellar winds and supernovae, photoelectric heating, and cosmic X-ray background heating from a metagalactic background. Our fiducial “MrAGN” model is identical except that it also includes a model for black hole seeding and accretion, as well as heating and momentum injection associated with the radiation from black hole accretion. Our MrAGN model launches galactic outflows, which result in both “ejective” feedback—the outflows themselves that drive gas out of galaxies—and “preventative” feedback, which suppresses the inflow of new and recycling gas. As much as 80% of outflowing galactic gas can be expelled, and accretion can be suppressed by as much as a factor of 30 in the MrAGN runs when compared with the NoAGN runs. The histories of NoAGN galaxies are recycling dominated, with ~70% of material that leaves the galaxy eventually returning, and the majority of outflowing gas reaccretes on 1 Gyr timescales without AGN feedback. Outflowing gas in the MrAGN runs has a higher characteristic velocity (500–1000 km s⁻¹ versus 100–300 km s⁻¹ for outflowing NoAGN gas) and travels as far as a few megaparsecs. Only ~10% of ejected material is reaccreted in the MrAGN galaxies.

**Key words:** galaxies: active – galaxies: evolution – galaxies: interactions – quasars: general

## 1. Introduction

It has been estimated by studying absorption lines associated with active galactic nuclei (AGNs) across a large range of luminosities that upwards of 60% of AGNs exhibit outflows, meaning that outflows are an important part of the process of supermassive black hole (SMBH) growth and accretion (Ganguly & Brotherton 2008). This AGN-outflow connection could even be universal if the typical covering fraction of outflows were ~60%. Evidence of these outflows being potentially very powerful comes in the form of observations of broad absorption line quasars, which make up ~20% of the observed QSO population (Knigge et al. 2008). These systems show evidence of outflows being launched at velocities as fast as 10,000 km s⁻¹ or more in the vicinity of the SMBH (Crenshaw et al. 2003; Moe et al. 2009; Dunn et al. 2010; Tombesi et al. 2010, 2013, 2015). There is also evidence for winds driven by AGN in low-redshift systems (Cicone et al. 2014; Cheung et al. 2016). However, it is quite difficult to constrain the size and extent of outflows, and therefore the mass outflow rate, based on absorption line observations (Arav et al. 2012; Maiolino et al. 2012; Arav et al. 2013; Chamberlain & Arav 2015). Recently, integral field unit (IFU) observations have been used to more fully map the gas around AGNs and gain insight into the kinematics of AGN-driven winds. These studies have revealed that outflows can act on large scales and entrain large quantities of gas on their way out, although it is

difficult to say for sure how much of a role AGNs play in driving these winds as opposed to stellar feedback processes (Tremonti et al. 2007; Prochaska & Hennawi 2009; Alatalo et al. 2011; Rupke & Veilleux 2011, 2013; Sturm et al. 2011; Zakamska & Greene 2014; Perna et al. 2015).

The overall effect that these outflows have on their host galaxies, specifically on the gas cycle and future star formation rate, is unclear, although in general outflows are important for the regulation of galaxy growth (Lilly et al. 2013). Some observational studies have investigated the relationship between AGN activity and star formation in galaxies, but the results have been controversial. Some studies conclude that winds can remove enough gas to cause negative feedback on star formation (Rupke & Veilleux 2013), although it is not clear if the gas will remain outside of the galaxy and whether the presence of a strong AGN and a suppression of star formation are correlated (Page et al. 2012). However, other studies find that AGNs preferentially live in galaxies on or even above the star-forming main sequence (Santini et al. 2012; Rosario et al. 2013a, 2013b), where they can drive powerful outflows (Genzel et al. 2014), or even that AGNs are stimulating star formation (positive feedback; Feain et al. 2007; Zinn et al. 2013). The disagreement between different observational studies originates from the different selection criteria used by various works and the different timescales on which galaxies and their central SMBHs operate (see Harrison 2017 for details).

On the theoretical side, AGN feedback is often invoked to explain many observed properties of galaxy populations. It is

<footnote n="7">⁷ Corresponding author.</footnote>
believed that almost all galaxies with a bulge component host an SMBH at their center (Magorrian et al. 1998). The mass of these central SMBHs is also known to be correlated with several properties of their host galaxies, specifically the mass, velocity dispersion, and luminosity of their bulge components (Kormendy & Richstone 1995; Ferrarese & Merritt 2000; Gebhardt et al. 2000; Tremaine et al. 2002; Marconi & Hunt 2003; Häring & Rix 2004). Meanwhile, there is a bimodality in the colors of the galaxy population (Baldry et al. 2004; Bell et al. 2004), and the growth of these two populations suggests that blue, star-forming galaxies are being transformed into red, quiescent galaxies, with the dearth of galaxies between the two peaks in color or star formation rate space indicating a fast transition timescale for a significant fraction of transitioning galaxies (Bell et al. 2004; Faber et al. 2007; Pandya et al. 2017).

With these phenomena in mind (and the large amounts of energy that can be released by an accreting SMBH; Lynden-Bell 1969), AGN feedback becomes an appealing explanation for this galaxy “quenching,” as well as for explaining the origin of the various BH scaling relations mentioned above. It has also frequently been invoked to solve the overcooling problem in massive galaxies, suppressing the late formation of stellar mass and bringing theoretical predictions into alignment with the observed high-mass ends of the stellar mass function and stellar-mass-halo-mass relation (Somerville & Davé 2015; Naab & Ostriker 2017, and references therein).

Theorists generally categorize AGN feedback into two broad types: “radiative” mode (sometimes called “quasar” mode or “bright” mode) and “jet” mode (Heckman & Best 2014). Radiative mode is associated with radiatively efficient accretion, relatively high accretion rates (above a few hundredths of the Eddington rate), and is thought to be fueled by a classical optically thick, geometrically thin accretion disk (Shakura & Sunyaev 1973). The hard radiation field emanating from the accretion disk can Compton- and photoheat as well as photoionize and photodissociate gas. In addition, radiation pressure on dust and free electrons can drive outflows, and this is likely the origin of the broad absorption line winds discussed above (Proga & Kurosawa 2009; Higginbottom et al. 2014; Gaskell et al. 2016). When accretion rates drop lower than about ~0.01 of the Eddington rate, the accretion becomes radiatively inefficient, with most of the energy instead emerging as highly collimated relativistic jets. These jets are often observed at radio frequencies, giving rise to the term “radio mode” feedback. The jets appear to be able to heat the diffuse hot halo gas via giant bubbles (seen in X-ray observations), sound waves, and weak shocks (see Fabian 2012; Heckman & Best 2014 for reviews).

Phenomenologically, theorists often speak of “ejective” feedback, in which star formation is quenched due to the removal of the ISM from the galaxy, and “preventative” feedback, in which star formation is eventually choked off by the lack of fuel, as the inflow of fresh or recycled gas is suppressed. As noted by Peng et al. (2015) and others, the implications and effects of these two types of mechanisms for various aspects of galaxy evolution will differ in important ways. It is often assumed that radiative mode feedback works in a solely ejective manner, while jet mode is solely preventative. (As a result, jet mode feedback is sometimes also referred to as “maintenance mode”). Indeed, these assumptions are built into the implementations of AGN feedback in most semi-analytic models. However, the work

presented here will call into question the first assumption, while observations of powerful outflows in systems with giant radio jets seem to challenge the second (Torresi et al. 2012).

Many previous simulations have included prescriptions for AGN feedback, although usually the radiative mode of feedback is implemented via deposition of thermal energy, while mechanical feedback is reserved for lower Eddington ratios and is associated with the jet mode of feedback (Springel et al. 2005; Di Matteo et al. 2012; Dubois et al. 2013; Hirschmann et al. 2014; Vogelsberger et al. 2014; Khandai et al. 2015; Schaye et al. 2015; Steinborn et al. 2015; Barai et al. 2016; Weinberger et al. 2017). Although such simulations can be successful at reproducing galaxy properties, it is unclear whether just depositing thermal energy associated with the radiative mode of BH accretion can drive winds similar to those seen in observations. Moreover, thermal energy input leads to halo X-ray luminosities in disagreement with observations. If no winds are launched, the thermal energy that is radiated away is concentrated in the center of the galaxy, resulting in halo X-ray luminosities that are too low (Bogdán et al. 2015), while if there are only weak winds, the radiating hot gas that gets pushed into the halo results in X-ray luminosities that are too high when compared with observations (Choi et al. 2012, 2014; Schaye et al. 2015).

The motivation behind the work presented here is that outflowing hot, shocked gas on scales below those that can be explicitly resolved in cosmological simulations can impart momentum to gas, which cannot be radiated away. The physics is similar to the momentum boost occurring at the end of the Sedov-Taylor phase in a supernova (SN) explosion. Of the simulations that have modeled AGN feedback, relatively few have included a momentum-driven prescription associated with radiatively efficient BH accretion. Momentum-driven wind scalings associated with AGNs were adopted in the “Santa Cruz” semi-analytic model (Somerville et al. 2008). Ostriker et al. (2010) examined the effect of momentum-driving on the feeding of AGNs in a series of one- and two-dimension simulations. Debuhr et al. (2011) investigated momentum-driven winds from radiation pressure in hydrodynamic simulations of idealized galaxy mergers. Eisenreich et al. (2017) investigated the effect of the AGN prescription used in this work on the metal content of a series of idealized elliptical galaxies. Hopkins et al. (2016) implemented BH feedback via broad absorption line winds in a similar manner to the mechanical AGN feedback model presented in Choi et al. (2012) and in the current paper. They investigated the interplay between stellar and AGN feedback in isolated galaxies and found that powerful outflows driven by momentum flux injected via AGN feedback strongly suppressed star formation by removing gas. None of these hydrodynamical works, however, were in a cosmological context. In terms of cosmological simulations, Anglés-Alcázar et al. (2017a) implemented momentum-driven outflows from AGNs in cosmological zoom simulations and found that they had a strong effect on galaxy and BH growth. In addition, the kinetic feedback model of Weinberger et al. (2017) was shown to bring the baryon content of galaxies in the Illustris simulation into better agreement with observations (Pillepich et al. 2018), although their mechanical feedback prescription is associated only with low black hole accretion rates (BHARs; jet mode).

In the cosmological zoom simulations that we study here, in addition to thermal energy, we also deposit momentum
associated with radiatively efficient accretion into the gas particles surrounding the BH, and, further, model photoionization and photoheating by this radiation. The prescriptions for momentum and radiation feedback from AGNs (MrAGN) presented in Choi et al. (2012, 2014, 2015) have been implemented into SPH GAL (Hu et al. 2014), an updated version of the SPH code GADGET-3 (Springel 2005). The code has also recently been modified to include updated treatments of many physical processes, including an improved treatment of stellar and SN feedback (Núñez et al. 2017), chemical enrichment and metal line cooling (Aumer et al. 2013), photoelectric heating, and X-ray heating by the metagalactic X-ray background. These updates are described in detail in Choi et al. (2017).

We address several questions: (1) What are the histories of inflow and outflow for these galaxies? (2) How much of the material ejected by winds is removed permanently, and how much comes back, and on what timescales? (3) How do these winds affect the host galaxy? We compare simulations including both AGN feedback and feedback due to stars and SNe (MrAGN) to a matched set of simulations that include only stellar and SN feedback (NoAGN). In this way, we can attempt to isolate the effect of the AGN-driven winds on the gas cycle in these halos, which have virial masses of ~$10^{12}-10^{13.4} M_\odot$, and consider how the answers to the questions posed above differ in the MrAGN versus NoAGN cases.

In Section 2, we describe the simulations used as well as our analysis methods, and in Section 3, we present our results. In Section 4, we discuss our results. We summarize in Section 5.

## 2. Simulation and Methods

In this work, we present simulated massive galaxies from cosmological zoom-in runs performed with a version of the parallel smoothed particle hydrodynamics (SPH) code GAGDET-3 (Springel 2005). Specifically, we use the modified code, SPH GAL, which mitigates problems that previous SPH codes have traditionally had with fluid mixing (Hu et al. 2014). A detailed description of the code can be found in Hu et al. (2014) and relevant updates can be found in Choi et al. (2017), but below we give a brief overview of the physics relevant to our study.

### 2.1. Code Basics and Setup

The code employs the pressure-entropy SPH formulation of Hopkins (2013) and also has improved force accuracy due to the use of the Wendland $C^4$ kernel with 200 neighboring particles (Dehnen & Aly 2012). We also include the improved artificial viscosity implementation presented by Cullen & Dehnen (2010) and an artificial thermal conductivity according to Read & Hayfield (2012) in order to reduce the noise in pressure estimates in the presence of strong shocks. Finally, a time-step limiter is employed according to Saitoh & Makino (2009) and Durier & Dalla Vecchia (2012) to ensure that neighboring particles have similar time steps and that ambient particles do not remain inactive when a shock is approaching.

### 2.2. Star Formation and Stellar/Supernova Feedback

Star formation and chemical evolution are modeled as described in Aumer et al. (2013); chemical enrichment is achieved via winds driven by Type I and II SNe and asymptotic

giant branch (AGB) stars. Eleven species of metals are tracked explicitly, and cooling rates, as well as the temperature and density of the gas, are calculated based on abundances. Redshift-dependent metagalactic UV/X-ray and cosmic microwave backgrounds are included with a modified Haardt & Madau (2012) spectrum.

Stars are formed stochastically where the density is greater than the density threshold for star formation. This threshold is given as $n_{\text{th}} \equiv n_0(T_{\text{gas}}/T_0)^3 (M_0/M_{\text{gas}})^2$, where $n_0 = 2.0 \text{ cm}^{-3}$ and $T_0 = 12,000$ K, with $M_0$ being the gas particle mass. This corresponds to the density value for the Jeans gravitational instability of a mass $M_{\text{gas}}$ at temperature $T_{\text{gas}}$. The star formation rate is calculated as $d\rho_*/dt = \eta\rho_{\text{gas}}/t_{\text{dyn}}$, where $\rho_*$, $\rho_{\text{gas}}$, and $t_{\text{dyn}}$ are the stellar density, gas density, and local dynamical time for the gas particle, respectively. $\eta$ is the star formation efficiency and is set to 0.025.

Stellar feedback is included in the form of stellar winds and heating by ionizing radiation from young massive stars. Momentum from stellar winds is added to the surrounding gas particles, while cold gas within the Strömgren radius of hot stars is heated to $T = 10^4$ K.

In the SN feedback model, SN energy and momentum are distributed to the surrounding ISM from the SN event. Depending on the distance from the SN events, we assume that each nearby gas particle is affected by one of three successive phases of SN remnant (SNR) evolution: the ejecta-dominated free expansion (FE) phase, the energy-conserving Sedov-Taylor blast-wave SNR phase, and the momentum-conserving snowplow phase. SN energy is transferred by conserving the ejecta momentum for gas particles within the radius of the FE phase. For gas particles lying outside the FE radius but within the Sedov-Taylor phase, the SN energy is transferred as 30% kinetic and 70% thermal. Finally, at larger radii in the snowplow phase, a fraction of the original SN energy is transferred as radiative cooling becomes significant. An allowance is also made for more efficient propagation of the remnant in a multiphase interstellar medium where appropriate. A detailed description of the implementation of the early stellar and SN feedback prescription can be found in Núñez et al. (2017).

Feedback from low- and intermediate-mass stars is also included in the form of slow winds as from AGB stars. Energy and momentum are transferred to surrounding gas particles such that momentum is conserved. The initial outflowing wind velocities are assumed to be $v_{\text{out,AGB}} = 10 \text{ km s}^{-1}$, which is typical for AGB stars (Nyman et al. 1992). Metal-enriched gas from all of these prescriptions is continuously added to the ISM. We also include metal diffusion, which allows for the mixing and spreading of metals in the enriched gas (see Aumer et al. 2013 for details).

### 2.3. Black Hole Growth and Feedback

Seed BHs are treated as collisionless sink particles and are placed in the centers of halos that reach a mass threshold ($10^{5.15} M_\odot$ BHs are placed in $10^{11.15} M_\odot$ halos). These BHs can then grow by gas accretion or by merging with another BH, as soon as the two BHs fall within each other's local SPH smoothing lengths and their relative velocities are smaller than the local sound speed. In the case of gas accretion, infall onto the BH is governed by the Bondi-Hoyle-Littleton rate (Bondi 1952). The size of the SPH gas particles is taken into account so that full accretion is allowed when the entire volume
of a particle is within the Bondi radius. Gas particles that fall
only partially inside the Bondi radius are given a probability of
being absorbed by the BH based on the volume within the
Bondi radius (Springel et al. 2005; Choi et al. 2012, 2014).

Massive disk winds driven from the nuclear region surrounding
the central SMBH are thought to produce many observed AGN
spectral features, including the broad absorption lines (e.g.,
Krigge et al. 2008). These wind-driving processes on scales
below our resolution limit, such as radiation pressure acting on
spectral lines, are treated as subgrid processes and assumed to
impart momentum to the gas. AGN-driven winds are launched
from the central region around the BH, with a fixed wind velocity
of $v_{\text{outf,AGN}} = 10,000$ km s$^{-1}$, and the number of particles selected
to receive a kick due to this wind is determined by a parameter for
feedback efficiency. The total energy flux carried by the wind is
$\dot{E}_w \equiv \epsilon_w \dot{M}_{\text{acc}} c^2$, where the efficiency parameter $\epsilon_w$ is set to
0.005 (Ostriker et al. 2010; Choi et al. 2017). The mass flux
and momentum flux carried by the wind are $\dot{M}_{\text{outf}} =
2 \dot{M}_{\text{acc}} \epsilon_w c^2 / v_{\text{outf,AGN}}^2$ and $\dot{p}_{\text{outf}} =
2 \epsilon_w \dot{M}_{\text{acc}} c^2 / v_{\text{outf,AGN}}$, respec-
tively, and therefore, for our selected wind velocity and feedback
efficiency, we have $\dot{M}_{\text{outf}} = 9 \dot{M}_{\text{acc}}$ and $\dot{p}_{\text{outf}} =
30 L_{\text{BH}}/c$. Thus,
90\% of the inflowing mass entering the central region is expelled
while 10\% is accreted onto the BH.

The selected gas particles receive the wind kick in a direction
parallel or antiparallel to their angular momentum vectors. The
emitted wind particle shares its momentum with its two nearest
neighbors to reproduce the shock-heated momentum-driven
flows. The momentum is split between the gas particles and is
conserved, but the kinetic energy of the outflow is not; the
excess energy is deposited into the gas particles as thermal
energy. This gives 2:1 divisions into thermal and kinetic
energies, which is similar to that in the Sedov–Taylor blast
wave, making the shocked winds approach the Sedov solution
faster, thus making us less sensitive to the resolution.

Radiation feedback from Compton and photoionization
heating due to X-ray radiation from the accreting BH, radiation
pressure associated with the heating, and the Eddington force is
also included. We utilize the AGN spectrum and metal line
heating prescription of Sazonov et al. (2004). X-ray radiation is
coupled to the surrounding gas using an approximation from
Sazonov et al. (2005). The radiation pressure on each gas
element is also calculated. Accretion is not artificially capped at
the Eddington rate, but the Eddington force acting on electrons
is included such that super-Eddington accretion can occasion-
ally occur but naturally reduces inflow while stimulating
outflows.

Also included are metallicity-dependent heating prescriptions due to photoelectric emission and metal line absorption. The emissivity of background AGNs is calculated as $\epsilon(z) = \epsilon \frac{d\rho_{\text{BH}}(z)}{dt}c^2$, where the radiative efficiency $\epsilon$ is set to 0.1. From this, the heating by the cosmic X-ray background is derived. More details about these and all of our feedback prescriptions can be found in Choi et al. (2017).

We note here that these feedback prescriptions have been
shown to produce galaxies with fairly realistic observable
properties (Choi et al. 2017). Without AGN feedback, our
model produces galaxies with stellar to halo mass ratios that are
too high by a factor of three and continues to form stars until
$z=0$, in conflict with observations of galaxies of this mass,
which are predominantly quenched. This results in compact
stellar cores, effective radii that fall a factor of five below the
observed size–mass relation, and high velocity dispersions

(Choi et al. 2017). The addition of AGN feedback alleviates all
of these problems, resulting in quiescent galaxies with
observationally consistent sizes and velocity dispersions at
$z=0$. Although AGN feedback alleviates the overcooling
problem and results in galaxy sizes in good agreement with
observations, our lower-mass galaxies, and high-redshift
galaxies, tend to have sizes that are too small for their mass
compared with observations (Choi et al. 2017). While the cold
gas fractions of our galaxies are lower than observed at low
redshift due to shortcomings in our treatment of the ISM and
stellar feedback, AGN feedback effectively removes hot gas
from halos, resulting in lower, more physical X-ray luminos-
ities and total gas mass fractions (Choi et al. 2017).
Postprocessing the MrAGN simulations with newly developed
nebular emission line models additionally shows that the
evolution of the optical nebular emission line ratios of massive
galaxies is widely consistent with observations (Hirschmann
et al. 2017). Satisfied that our model produces galaxies that
compare favorably with observed high-mass elliptical galaxies,
in this work we turn to analyze the wind properties and gas
cycle in these systems.

## 2.4. Zoom Simulations

The “zoom-in” initial conditions that we use are described in
detail in Oser et al. (2010, 2012). We run our zooms on dark
matter halos picked from a dark matter only simulation that
employed the parameters from WMAP3 (Spergel et al. 2007)
and assumed a flat cosmology: $h = 0.72$, $\Omega_b = 0.044$,
$\Omega_{\text{dm}} = 0.216$, $\Omega_\Lambda = 0.74$, $\sigma_8 = 0.77$, and an initial power
spectrum slope $n_{\text{rms}} = 0.95$. We trace back the dark matter
particles close to the halos of interest in each snapshot and then
replace those dark matter particles with high-resolution dark
matter and gas particles. The high-resolution zoom is then
evolved from $z=43$ to today.

The selected dark matter halos have final virial masses between
$1.4 \times 10^{12} M_{\odot}$ and $2.3 \times 10^{13} M_{\odot}$ and are made up of dark
matter particles with a mass of $m_{\text{dm}} = 3.57 \times 10^7 M_{\odot}$. The final
central galaxy masses in these halos (for our runs with AGNs) are
between $8.2 \times 10^{10} M_{\odot}$ and $1.5 \times 10^{12} M_{\odot}$, with gas and star
particles both having a mass of $m_{\star,\text{gas}} = 6.0 \times 10^6 M_{\odot}$. We use
comoving gravitational softening lengths of $\epsilon_{\text{gas,star}} = 0.556$ kpc
and $\epsilon_{\text{halo}} = 1.236$ kpc for gas/star particles and dark matter
particles, respectively.

In what follows, we will examine zoom regions run with two
different models, MrAGN and NoAGN. MrAGN is the fiducial
model, which includes all of the physics described above,
including the different ways in which AGNs can affect
feedback on their surrounding galaxies. NoAGN is the same
in every way, except that it contains no BHs and thus no AGN
feedback. In this way, the effects of the different feedback
mechanisms can be isolated. Out of our initial sample of 30
pairs of zoom runs, we focus on the 24 that we can use to make
a direct comparison between the two models, as described
below.

## 2.5. Methods

### 2.5.1. *Galaxy Matching*

In order to compare MrAGN galaxies with their NoAGN
counterparts, and also to examine inflowing and outflowing gas
in our simulations, we must first find the center of our halos of
interest and make sure that we are tracking the same progenitor
back with redshift in both runs. We utilize two different codes to track the halo centers: ROCKSTAR (Behroozi et al. 2013) and GTRACE, a tool built specifically for analyzing GADGET snapshots. In 24 of our 30 pairs, at least one of these tools found the center of the same progenitor in both runs by $z=1$. Going back to $z=2$, we have 23 matched pairs, and by $z=3$, we have 21 matched pairs. There are brief periods when even the center of these matched pairs is lost, such as when a satellite of comparable mass approaches the main progenitor, but this is generally only for a short time and can be easily identified when a relatively stable quantity, such as stellar mass, suddenly dips, as can be seen in Figure 1 below. We exclude the six galaxies for which the found centers do not correspond to the same progenitors by $z=1$.

### 2.5.2. Particle Tracking

With the centers found, we track the flow of particles across two shells: one at 10% of the virial radius and one at the virial radius. We refer to the former as the “galaxy radius” $r_g$ and the latter as the “halo radius” $r_h$. The tracking begins when there is a reliable center found for both the MrAGN and NoAGN runs; for our three case study galaxies, for which it is important to track over the same period (as we display cumulative quantities), this tracking begins at $z \sim 4$, give or take one time step. Almost all of the rest of our galaxies are tracked by $z=3$, as mentioned above. At each time step, the gas particles that are inside both of these radii are catalogued. In the next time step, any particle that was inside (outside) the radius of interest in the last time step and is now outside (inside) that radius and has a positive (negative) radial velocity is considered outflowing (inflowing). Note that the typical time step between snapshots at around $z=0$ is 0.13 Gyr.

The inflowing and outflowing mass at each time step is stored, and the inflow and outflow rate can be calculated by using the time step between snapshots. We also keep track of whether a gas particle is accreting or outflowing for the first time, or if it has done so previously. Finally, we keep track of how long it takes for individual gas particles to be recycled, as well as the maximum radial displacement experienced by the particle during each recycling event (see Übler et al. 2014 for a similar treatment of gas particles in disk galaxies). Unless otherwise specified, stellar and gas masses of the galaxy are calculated within 10% of the virial radius.

## 3. Results

Here we present the results of our analysis, first examining the detailed histories of a few representative galaxies, then looking at broader trends in galaxy and inflow/outflow properties for the entire set of galaxies. We note again that outflows are driven by both stars and SNe and AGNs. Our analysis also captures inflow and outflow due simply to thermal motions of gas. This is more apparent in our NoAGN galaxies, as will be mentioned below.

### 3.1. Case Studies

When examining the histories of individual galaxies in terms of their basic properties and the properties of their inflows and outflows, some broad classes emerge, mainly as a result of galaxy mass. Below we present three galaxies representative of their mass bins: “high mass” m0163 ($M_{*,final} \sim 10^{11.4} M_\odot$ and $M_{h,final} \sim 10^{13.1} M_\odot$), “intermediate

mass” m0329 ($M_{*,final} \sim 10^{11.3} M_\odot$ and $M_{h,final} \sim 10^{12.7} M_\odot$), and “low mass” m0501 ($M_{*,final} \sim 10^{11.2} M_\odot$ and $M_{h,final} \sim 10^{12.5} M_\odot$). We study the baryon cycles in these halos in detail and use them as exemplars of more general trends that we discuss later.

#### 3.1.1. Galaxy Property and Gas Flow Histories

Halo m0163 has a history characteristic of the more massive end of our population. It has a final halo mass of $log(M_h/M_\odot) \sim 13.1$ and a final stellar mass of $log(M_*/M_\odot) \sim 11.4$ in the MrAGN run. Figure 1 shows the evolutionary histories of different galaxy, BH, inflow, and outflow properties. In panels (a) and (e), we see that the halo mass (shown multiplied by the universal baryon fraction, 0.1658) is largely unaffected by AGN feedback, while the final stellar mass is reduced by roughly 0.4 dex. While the cold gas mass ($T < 2 \times 10^4$ K) decreases only mildly in the NoAGN run, all of the cold gas within $r_g$ is removed or heated by AGN feedback in the MrAGN run. The hot gas mass of the MrAGN galaxy is only slightly smaller than that for the NoAGN galaxy, due to the removal of some of the hot gas from the halo. We can also see the growth of the BH in the MrAGN run, which reaches a final mass of $log(M_*/M_\odot) = 9.2$. The green dashed lines denote mergers for which the halo ratio is greater than 1:10.

Panels (b) and (f) show that as the cold gas mass decreases, the star formation rate falls just as drastically, which is partially responsible for the smaller final stellar mass found in the MrAGN run. Star formation continues at a nearly constant rate in the NoAGN run (several tens of solar masses per year). We also see the BH accretion rate (BHAR), which gradually increases from ~4 Gyr until about 12 Gyr after the start of the simulation. After this point, BH accretion itself is quenched.

Panels (c), (d), (g), and (h) show the inflow and outflow rates for both runs on galaxy scales $r_g$ and halo scales $r_h$. In the MrAGN run, the inflow rate at $r_g$ is suppressed relative to the NoAGN run value of $\sim 50 M_\odot \text{yr}^{-1}$ to only $1-2 M_\odot \text{yr}^{-1}$, and overcome by the outflow rate after the two early merger events at $\sim 3-5 \text{ Gyr}$ (just after $z \sim 2$). The outflow rate continues to mirror the inflow rate thereafter, at $\sim 10-30 M_\odot \text{yr}^{-1}$. At $r_h$, this same spike in outflow rate occurs after a slight delay, due to material that was caught up in the initial outflow at $r_g$ and eventually crosses the virial radius, pushing out even more material on its way. This results in outflow rates of as much as $\sim 300 M_\odot \text{yr}^{-1}$. The outflow rate at $r_h$, however, rarely overtakes the inflow rate, even in the MrAGN run. In the NoAGN run, inflow dominates at almost all times at both radii. The outflow rate at $r_g$ in the NoAGN case steadily rises, eventually reaching an equilibrium with the inflow rate at a value higher than is seen in the MrAGN case. This is due to the fact that we have more gas in the central region in the NoAGN case, which is cycling in and out of the galaxy due to thermal motions. The inflow rate at $r_h$ is largely unaffected by AGN feedback, while the outflow rate in the MrAGN run is enhanced relative to the NoAGN case. This means that the AGN feedback acts mostly in an ejective way at this halo mass ($\sim 10^{13.1} M_\odot$) at halo scales, rather than preventatively. See Table 1 for the cumulative inflow and outflow gas masses that cross both shells over the duration of the simulation. Already apparent from Table 1 is the importance of preventative feedback, as represented by the cumulative mass in baryons that accreted onto the halo divided by the final halo mass times the universal baryon fraction.
<figure><img src="image_3.png" /></figure>

**Figure 1.** Baryon cycle history of the main progenitor galaxy in halo m0163. Left column: results for the MrAGN run. Right column: results for the NoAGN run. Top row: evolution of galaxy properties—baryon fraction times halo mass, stellar mass, total black hole mass (for MrAGN run), cold gas mass within the galaxy radius ($T < 2 \times 10^4$ K), and hot gas mass within the halo radius ($T > 2 \times 10^4$ K). Note that we show the total black hole mass within the galaxy radius $r_g$. Second row: evolution of star formation rate and black hole accretion rate (for the MrAGN run; black hole accretion rate is scaled up by a factor of 1000). Third row: inflow and outflow rates at the galaxy radius. Fourth row: inflow and outflow rates at the halo radius. Green vertical dashed lines indicate halo merger events for which the halo mass ratio is 1:10 or greater. The final stellar mass in the MrAGN run is smaller due to the decrease in cold gas supply and thus star formation rate. AGN feedback does not affect the inflow of gas into the halo but does mildly affect the inflow rate of gas at galactic scales. It also enhances outflows on both scales at early times.
<table><caption>**Table 1**
Final Properties of Three Example Galaxies</caption><thead><tr><th>Galaxy</th><th>Halo Mass</th><th>Stellar Mass</th><th>BH Mass</th><th>Hot Gas Mass</th><th>Inflow<sub>g</sub></th><th>Outflow<sub>g</sub></th><th>Inflow<sub>h</sub></th><th>Outflow<sub>h</sub></th><th>M<sub>in,h</sub>/f<sub>h</sub> M<sub>h,0</sub></th></tr></thead><tbody><tr><td>m0163 (MrAGN)</td><td>13.1</td><td>11.4</td><td>9.2</td><td>11.9</td><td>11.2</td><td>11.2</td><td>12.2</td><td>11.9</td><td>0.90</td></tr><tr><td>m0163 (NoAGN)</td><td>13.1</td><td>11.8</td><td>N/A</td><td>12.1</td><td>11.8</td><td>11.5</td><td>12.3</td><td>11.6</td><td>1.00</td></tr><tr><td>m0329 (MrAGN)</td><td>12.7</td><td>11.3</td><td>9.2</td><td>11.2</td><td>10.9</td><td>10.7</td><td>11.9</td><td>11.7</td><td>0.85</td></tr><tr><td>m0329 (NoAGN)</td><td>12.8</td><td>11.7</td><td>N/A</td><td>11.7</td><td>11.5</td><td>10.9</td><td>12.0</td><td>11.1</td><td>1.05</td></tr><tr><td>m0501 (MrAGN)</td><td>12.5</td><td>11.2</td><td>9.3</td><td>7.3</td><td>11.0</td><td>10.8</td><td>11.5</td><td>11.4</td><td>0.57</td></tr><tr><td>m0501 (NoAGN)</td><td>12.6</td><td>11.5</td><td>N/A</td><td>11.4</td><td>11.4</td><td>10.8</td><td>11.8</td><td>10.9</td><td>1.00</td></tr></tbody></table>

Note. The final halo, stellar, black hole, and hot gas masses; the cumulative gas inflow and outflow masses that crossed shells at the galaxy radius and the halo radius; and the cumulative inflow mass at the halo radius divided by the final halo mass times the universal baryon fraction for our three example galaxies, both for our MrAGN run and our NoAGN run. All masses are given in units of log solar masses.

Figure 2 is the same as Figure 1, but now for galaxy m0329. This galaxy falls near the middle of our mass range, with a final halo mass of $log(M_h/M_\odot) \sim 12.7$ and a final stellar mass of $log(M_*/M_\odot) \sim 11.3$ in the MrAGN run. We see again the difference between the stellar mass (again $\sim 0.4$ dex) and cold gas mass in the MrAGN run versus those in the NoAGN run. The galaxy in the MrAGN run has a steeply decreasing cold gas mass, resulting in a steeply decreasing SFR. In the case of m0329, this also corresponds to an increasing BHAR. The suppression of inflow at both radii of interest is more pronounced than for m0163; the inflow rate at $r_g$ is significantly (and permanently) decreased from $>10 M_\odot \text{yr}^{-1}$ to $1-2 M_\odot \text{yr}^{-1}$ following a strong outflow event at $\sim 4-6$ Gyr, concluding at $z \sim 1$. The inflow rate at $r_h$ is also more noticeably decreased after this outflow event than it was for the more massive halo m0163; the final inflow rate at $r_h$ is $\sim 1.7$ times larger in the NoAGN run than in the MrAGN run. This suggests that as we look at smaller halo masses and therefore shallower potential wells, preventative feedback, where AGN feedback not only removes material but also prevents new material from accreting, becomes more important (see the rightmost column of Table 1). The outflow rates at both radii match or exceed the inflow rates at almost all times after the initial outflow event in the MrAGN run. In the NoAGN case, inflows and outflows behave much as they did for m0163, with inflow always dominating.

Finally, Figure 3 is the same as Figures 1 and 2, but now for galaxy m0501. This galaxy is less massive than m0329, with a final halo mass of $log(M_h/M_\odot) \sim 12.52$ and a final stellar mass of $log(M_*/M_\odot) \sim 11.22$ in the MrAGN run, and has a history more characteristic of the low-mass galaxies in our suite. In the top four panels, the galaxy properties of m0501 evolve very similarly to those already examined for m0163 and m0329, except for hot gas in the MrAGN run, which is much more strongly affected. When focusing on the inflow and outflow properties, m0501 is very different. In the NoAGN case, things are much the same, with inflow dominating outflow at both radii at almost all times. In the MrAGN case, several bursts of outflow of between $10-30 M_\odot \text{yr}^{-1}$ dominate inflow at $r_g$ at around 4 Gyr. This outflow is powerful enough to halt inflow, after which the outflow stops as well because the gas within $r_g$ has been completely depleted. At $r_h$, we see this outflow, having swept up gas in the halo and now removing $\sim 100 M_\odot \text{yr}^{-1}$, peak at a slightly later time and once again

halt inflow across the virial radius. This outflow then also tapers off as it has cleared most of the galaxy's halo gas as well. This is an extreme case where AGN feedback acts in both an ejective and a preventative way on dark matter halo scales. This is fairly rare, occurring only in 4 of our 30 MrAGN galaxies. These galaxies all have relatively small halo masses for our sample (the most massive has $M_h \sim 10^{12.5} M_\odot$), and all contain a relatively large BH mass for their halo mass, although they still sit on the $M-\sigma$ relation.

### 3.1.2. Gas Morphology and Velocity Structure

In Figure 4, we show vector maps of the gas velocities in m0163 overlaid onto color maps of the gas temperature for both runs at four different redshifts. The luminance of the gas in the images corresponds to its density. At $z \sim 3$ and $z \sim 2$, the images of the two runs look very similar to each other, with the main velocity features being filaments of gas funneling into the forming galaxy. At $z \sim 1$, however, the two runs look very different. The galaxy in the MrAGN run is going through a large bout of outflow, which can be seen in panels (c) and (d) of Figure 1 at $\sim 5$ Gyr. Meanwhile, in the NoAGN run, the bulk of the velocity features are still due to inflowing material onto the galaxy. At $z \sim 0$, while the two runs look very similar, the MrAGN run exhibits a more diffuse hot gas halo. While gas has been heated by stellar and SN feedback, the NoAGN run still exhibits somewhat ordered motions and inflows as opposed to the MrAGN run, in which we see material being pushed out of and away from the galaxy.

In Figure 5, the galaxy m0329 in the MrAGN run is undergoing a major outflow by $z \sim 2$, in contrast to its NoAGN counterpart, which is dominated by filamentary accretion. At $z \sim 1$, while the NoAGN galaxy also appears to be undergoing feedback as demonstrated by the bulk of hot gas, there is no strong outflow signature as we see again in the MrAGN case. Even in the MrAGN run, however, there is still a strong filamentary inflow feature outside of the outflow sphere of influence. When this inflowing gas reaches the outflow, some of it is heated and/or turned around. At $z \sim 0$, we again end up with two galaxies with similar (to the eye) gas contents, although we know from Figure 2 that there is more to the story.

Figure 6 shows that the two runs of m0501 are very similar until $z \sim 1$, at which point a large AGN-driven outflow occurs in the MrAGN case. This outflow goes on to clear the galaxy of gas and eventually destroy the filament supplying the galaxy
<figure><img src="image_4.png" /><figcaption>**Figure 2.** Same as Figure 1, but for halo m0329. Inflow in the MrAGN case is more noticeably suppressed at late times at both radii of interest than in the more massive halo m0163. Following an initial outflow event at ~4–6 Gyr, outflow dominates inflow at the galaxy radius by a small amount, but consistently. In the NoAGN case, inflow dominates outflow at all times at both radii.</figcaption></figure>
<figure><img src="image_5.png" /><figcaption>**Figure 3.** Same as Figures 1 and 2, but for halo m0501. While the NoAGN case is similar to the previous galaxies, the MrAGN case is exceptional in that AGN feedback completely clears the galaxy of gas and prevents any new gas from accreting for several billion years, even on the scale of the halo.</figcaption></figure>
<figure><img src="image_6.png" /><figcaption>**Figure 4.** Gas temperature and velocity vector maps for galaxy m0163. Images created using `pygad` (Röttgers 2017). In the gas image, the color corresponds to gas temperature, while the luminance corresponds to gas density. The velocity vectors are calculated for a slice with a thickness of the halo radius at the redshift of interest, with contributing gas velocities weighted by their densities. The arrow lengths are normalized by the number of bins, and then normalized by the average arrow length. Top row: results for the MrAGN run. Bottom row: results for the NoAGN run. The velocity vectors are overlaid onto gas temperature maps at four redshifts. White circles denote the current galaxy radius and the halo radius. At *z* ~ 3, the two runs are nearly identical and dominated by accretion along filaments. At *z* ~ 2, the two runs are still very similar, with the bulk velocity flow due to gas inflowing onto the central galaxy. At *z* ~ 1, the MrAGN run is undergoing a bout of outflow, while the NoAGN run is still steadily accreting gas. At *z* ~ 0, though the remnants appear similar, the gas in the MrAGN run is more diffuse.</figcaption></figure>

<figure><img src="image_7.png" /><figcaption>**Figure 5.** Same as Figure 4, but for halo m0329. At *z* ~ 3 and *z* ~ 0, the two galaxies look very similar, while at *z* ~ 2 and *z* ~ 1, we see the MrAGN galaxy at the height of its AGN-driven outflow activity in stark contrast to the NoAGN case.</figcaption></figure>

with gas, halting accretion even at the virial radius. At *z* ~ 0, there is basically no gas left to track in the MrAGN run, whereas in the NoAGN run, the gas has actually settled into a cold disk on galactic scales. Our mechanical AGN feedback is very efficient at removing gas from the lower-mass galaxies of our sample.

### 3.1.3. Recycling Fractions

In Figure 7, two different measures of the instantaneous recycling fraction are plotted. The blue curves represent the fraction of inflowing material across $r_g$ that was previously within the galaxy, while the green curves represent the future recycling fraction, i.e., the fraction of currently
<figure><img src="image_8.png" /><figcaption>**Figure 6.** Same as Figures 4 and 5, but for halo m0501. Between $z \sim 1$ and $z \sim 0$, AGN-driven outflows drive the gas completely out of the galaxy in the MrAGN run, even destroying the gas filament that was supplying the galaxy with new gas.</figcaption></figure>

outflowing material across $r_g$ that will return to the galaxy at some future time. The recycling fraction necessarily falls toward zero at late times, as outflowing material runs out of time to fall back in before the end of the simulation. From top to bottom, we see both the MrAGN and NoAGN runs for m0163, m0329, and m0501. In m0163, much of the outflowing material at early times in the MrAGN run (50%–70%) is destined to come back until the first big bout of outflow at ~4 Gyr, after which there are only a few isolated incidents of outflowing material that will eventually fall back in. This decrease in the future recycling fraction also corresponds with the BHAR starting to increase, as can be seen in panel (c) of Figure 1.

As a result, the fraction of inflowing material that is recycled never rises much above ~20%. The future recycling fraction in the NoAGN run remains very large, between 60% and 90% for a large portion of the galaxy’s history. As a result, the inflowing recycling fraction by the end of the simulation is almost 70%. In the MrAGN run, of the total $10^{11.2} M_{\odot}$ of gas that accreted onto the galaxy throughout its history, $10^{10.1} M_{\odot}$, or ~8%, was contributed by recycled material. By comparison, in the NoAGN run, 31% of the total $10^{11.8} M_{\odot}$ was recycled material. Conversely, the fraction of the cumulative outflow mass that was destined to come back when it left is 8.1% for the MrAGN run and 67.8% for the NoAGN run.

The galaxy m0329 behaves similarly to m0163, but there are some differences. The drop in future recycling fraction is more extreme than in m0163, with the inflowing recycling fraction rarely rising above a few percent. In the NoAGN galaxy, the future recycling fraction still reaches heights of 90% and greater, but the inflowing recycling fraction does not exceed 40%. This is partially because less material is being recycled than in m0163, but also because outflowing material falls back

in on shorter timescales, leading to a more consistent inflowing recycling fraction. The total contribution of recycled material to the cumulative inflow mass in the MrAGN galaxy is ~3%; in the NoAGN galaxy, it is ~16%. The fraction of the cumulative outflow mass contributing to the future recycling fraction is 5.1% for the MrAGN run and 78.0% for the NoAGN run.

The galaxy m0501 shows similar trends, except that here there is no inflow or outflow in the MrAGN case after ~6 Gyr. We can still see that the inflowing recycling fraction by $z=0$ in the NoAGN case is ~50%. The final total contribution of recycled material to the cumulative inflow mass is ~9% for the MrAGN galaxy and ~17% for the NoAGN galaxy. The fraction of cumulative outflow mass contributing to the future recycling fraction is 14.2% in the MrAGN run and 69.3% in the NoAGN run.

### 3.1.4. Accretion of Gas and Stars

The fraction of the $z=0$ total galactic baryonic mass ($M_{\text{gas}}+M_*$ within $r_g$) of m0163 that is made up of accreted stars (formed in separate halos that have merged with the main progenitor, as opposed to being formed in situ) is ~63% in the MrAGN run and ~33% in the NoAGN run. These fractions are shown in the far-left bars of the top row of Figure 8. Accreted stars make up a larger fraction of the final baryonic mass of the MrAGN galaxy than of the NoAGN galaxy because of the reduction of in situ star formation in the former. Also plotted is the total gas mass that has accreted since $z \sim 4$ divided by the total baryonic mass at $z=0$, in bins of the number of times the gas particle accreted onto the galaxy. This is further subdivided into the fates of the accreted gas particles. Gas in the NoAGN run is more likely to be accreted several times and is also more likely to remain in the galaxy in the form of gas or stars. In the MrAGN galaxy, gas tends to accrete fewer times, and the gas
<figure><img src="image_9.png" /><figcaption>**Figure 7.** Two different measures of instantaneous recycling fraction for our case study galaxies. The green curves show the future recycling fraction, i.e., the fraction of ejected material that will return at some future time. The blue curves show the inflowing recycling fraction, i.e., the fraction of currently inflowing material that was previously inside the galaxy. Left panel: MrAGN. Right panel: NoAGN. Top panels: m0163. Middle panels: m0329. Bottom panels: m0501. In all cases, both recycling fractions are much higher in the NoAGN runs.</figcaption></figure>

that is accreted is more likely to be outside of the galaxy by $z=0$. The fraction of total accreted gas that is ejected from the galaxy by $z=0$ in the MrAGN run is 67.1% as compared with 31.4% in the NoAGN run.

For m0329, the fraction of the $z=0$ baryonic mass that was accreted as stars is again larger for the MrAGN galaxy than for the NoAGN galaxy (~53% versus ~24%). Also apparent is that practically none of the accreted gas remains as gas within the MrAGN galaxy at $z=0$, whereas some of this gas remains in the NoAGN run, while a larger fraction is turned into stars. In the MrAGN run, 32.6% of all accreted gas is ejected by $z=0$, while only 10.7% is ejected in the NoAGN run. These trends continue for m0501, which again has a larger accreted stellar fraction in the MrAGN run (~39% versus ~24%). There is a slight preference for more recycling events in the NoAGN galaxy, and again none of the accreted gas in the MrAGN galaxy remains within the galaxy by $z=0$. Of the accreted gas, 20.8% is ejected in the MrAGN run, while the rest is converted into stars. This figure is relatively small when compared with the 14.4% in the NoAGN run, but this is because m0501 never

<figure><img src="image_10.png" /><figcaption>**Figure 8.** Total mass of accreted stars and gas that accreted onto the galaxy N times divided by the total baryonic mass of the galaxy at $z=0$. The baryonic mass is given in each panel. Left panels: MrAGN. Right panels: NoAGN. Top panels: m0163. Middle panels: m0329. Bottom panels: m0501. The far-left bar indicates the mass of accreted stars, which only ever accrete once, divided by the final total baryonic mass within the galaxy radius. The rest of the bars indicate gas particles that accreted N times. These bars are subdivided into the gas particle's fate at $z=0$: blue if it is still gas within the galaxy, gray if it is gas outside of the galaxy, and red if it is now a star particle. Accreted stars make up a larger fraction of the final baryonic mass of MrAGN galaxies than of NoAGN galaxies. Gas in the NoAGN run is more likely to accrete several times, and accreted gas is more likely to remain in the galaxy or form stars than in the MrAGN run, in which a significant fraction of accreted gas ends up being ejected.</figcaption></figure>

has an opportunity to accrete any more gas that would be affected by future AGN feedback.

### 3.1.5. Recycling and Ejection Timescales

Figure 9 shows the distributions of the timescales of gas recycling and gas ejection events for our case study galaxies. We define the recycling timescale as the amount of time between the ejection of a gas particle in a wind and reaccretion into a galaxy. Every time a gas particle crosses a shell at $r_g$, the time is recorded. The time is recorded again if it flows back into the galaxy. This may happen multiple times for a single gas particle, each of which is recorded as a separate recycling event and is included in the distribution. If the gas particle outflows but never comes back, the ejection timescale is recorded as the age of the universe minus the time of ejection. The top row shows us that in the MrAGN run of m0163, there is far less recycling on short timescales than in the NoAGN counterpart.
<figure><img src="image_11.png" /><figcaption>**Figure 9.** Distribution of timescales of recycling and ejection events for our case study galaxies. Top row: m0163. Middle row: m0329. Bottom row: m0501. Left panels: MrAGN runs. Right panels: NoAGN runs. The blue histograms show the distribution of the recycling timescale, the time spent outside the galaxy by gas that has been ejected past the galaxy radius until it later came back. The red histograms show the distributions of the time since ejection for gas particles that have exited the galaxy at the galaxy radius and have not yet returned. There are more recycling events in the NoAGN runs. Recycled gas particles in higher-mass NoAGN galaxies like m0163 tend to have longer recycling timescales than those in lower-mass galaxies. Ejected gas particles in lower-mass MrAGN galaxies like m0501 are more likely to have been ejected at early times by a powerful bout of AGN feedback than those in higher-mass galaxies.</figcaption></figure>

This appears to be because of a few very long-term (seemingly permanent) ejection events that occurred 8–10 Gyr ago. Those gas particles were then not available to be recycled by stellar and SN feedback processes as in the NoAGN run. At later times (ejection timescales $\lesssim 2$ Gyr), the ejection distributions look very similar. This is most likely due to random motions of gas particles in both runs that have not had time to return to the galaxy.

The NoAGN run of m0329 also has many more recycling events than its MrAGN counterpart. This time, however, there are consistently more ejection events of longer timescales in the MrAGN case, and these are not only gas particles outflowing at early times as for m0163. The AGN activity more consistently removes material throughout this galaxy’s history, again removing gas that would have contributed to the short-timescale recycling events due to stellar and SN feedback. We also see the preference for shorter-timescale recycling in

the NoAGN galaxy as compared with m0163, which is responsible for the relatively constant inflowing recycling fraction seen in the middle row of Figure 7.

A cumulative gas mass of $\sim 4 \times 10^{10} M_{\odot}$ in the MrAGN run of m0501 is ejected in a series of outflows between 7 and 12 Gyr ago. The NoAGN run is dominated by recycling events, with far fewer gas particles being ejected ($<10^{10} M_{\odot}$) relative to being recycled on short timescales ($>4 \times 10^{10} M_{\odot}$).

### 3.1.6. Recycling and Ejection Displacements

Figure 10 shows the distribution of maximum displacements of gas particles, again both for completed recycling events and ejection events in which the particle is still outside of the galaxy at the end of the simulation. When a particle is first tagged as outflowing, its current radius is recorded as the initial radius (usually just outside of $r_g$). We keep track of its subsequent radius and store the maximum radius reached during a given recycling or ejection event. We then subtract the initial radius to get the maximum displacement experienced by a given gas particle. For m0163, we see once again that there is far more mass ejected ($\sim 10^{11} M_{\odot}$) than recycled ($<10^{10} M_{\odot}$) in the MrAGN run, due to a larger number of gas particles, which are part of outflows at early times and never reaccrete onto the galaxy. It is also evident that the outflowing gas particles removed by AGN feedback travel much larger distances than those removed by stellar and SN feedback. In some cases, gas particles travel several megaparsecs, while in the NoAGN run only the very tail of the distribution reaches $\sim 1$ Mpc.

The distribution of maximum displacements for gas particles in the MrAGN run of m0329 is even more extreme than for m0163, with gas particles being driven out beyond the zoom region of our simulation. While we keep track of these distances, here we cap the particles’ displacements at the boundaries of our zoom region, which ranges from 5 to 10 Mpc for our 24 halos. Outside of this zoom-in region, only the low-resolution dark matter particles are populated to approximate the long-distance tidal force (Oser et al. 2010), and therefore the wind particles pushed outside of the zoom-in region can travel even farther. However, outflows can be stalled due to intergalactic medium pressure and then fall back to the galaxy halo in reality. This is the caveat that arises from the limited size of the zoom-in region; however, its effect on our infall rate measurement is believed to be small, given that the displacement of the recorded recycling events is two orders of magnitude smaller than the scale of the zoom-in region. Ejected particles in the NoAGN case can travel several hundred kiloparsecs, but only the extreme tail of the distribution exhibits radii larger than this. In both cases, the recycled material is restricted to these smaller displacements as well.

In m0501, with the halo depleted and the surrounding filament destroyed, there is nothing to stop ejected gas from traveling tens of megaparsecs in the MrAGN run. In contrast, ejected particles in the NoAGN run mostly only travel 100–200 kpc. In both runs, recycled gas particles are restricted to smaller radii, as is expected.

### 3.2. Broad Trends in Galaxy and Outflow Properties

In the last section, we focused on three individual example galaxies and found that the strength of gas outflows and inflows depends on halo mass. We now turn to the broad trends and scaling relations for galaxy and outflow properties that can be
<figure><img src="image_12.png" /><figcaption>**Figure 10.** Distribution of displacements for recycling and ejection events in our case study galaxies. Displacement is measured from the particle's initial radius outside of the galaxy radius, recorded when the particle is first tagged as outflowing. Top row: m0163. Middle row: m0329. Bottom row: m0501. Left panels: MrAGN runs. Right panels: NoAGN runs. The blue histograms show the distribution of the maximum distance traveled by gas particles that have been ejected past the galaxy radius and later come back. The red histograms show the distributions of the maximum distance for gas particles that have exited the galaxy past the galaxy radius and have not yet returned. Recycled particles that only ever travel back inwards from their initially recorded radius outside the galaxy are counted in our lowest bin. Gas particles ejected by AGN feedback travel much farther than those ejected by stellar and supernova feedback, and this trend becomes more extreme toward lower halo masses.</figcaption></figure>

extracted from our complete suites of galaxies. In the following
plots, each point represents a single galaxy at the specified
redshift. Some galaxies do not have a matched pair or reliable
center by *z* = 3, or in some cases *z* = 2, and so these are not
plotted at high redshift.

3.2.1. Galaxy Properties

Figure 11 shows the gas and stellar mass of all of our galaxies versus halo mass for four different redshifts. The full cold gas mass and stellar mass histories of our three case study galaxies can be found in Figures 1–3. At z = 3, before most
<figure><img src="image_13.png" /><figcaption>**Figure 11.** Gas and stellar content of all galaxies at $z=0$ vs. halo mass. Top row: total mass of gas within 10% of the virial radius, $r_g$. Middle row: cold gas mass ($T < 2 \times 10^4$ K) within $r_g$. Bottom row: stellar mass. The solid black curves are the abundance matching estimates from Moster et al. (2013), and the dashed black curve is the $z=0$ abundance matching estimate of Kravtsov et al. (2014). Black circles denote galaxies in the NoAGN run, while red triangles denote galaxies in the MrAGN run. Quantities plotted at $10^7 M_\odot$ with downward arrows are galaxies for which the y-axis value is 0. By $z=0$, all MrAGN galaxies have significantly less gas (and in some cases no gas), and all have been cleared of cold gas completely. All MrAGN galaxies also have lower stellar masses, as was seen in the individual galaxy histories shown above.</figcaption></figure>

major AGN activities (when SMBHs are still fairly small), galaxies in both runs have very similar gas content and fall on a fairly tight relation between gas mass and halo mass. At $z=2$, AGN feedback begins to kick in, and MrAGN galaxies start to drop below the NoAGN relation. By $z=1$, much of the gas affected by star formation and stellar feedback has fallen back into the galaxies, while in the MrAGN runs, much of this gas is lost. At $z=0$, galaxies in the two runs fall on distinct relations, with some of the lowest halo mass galaxies having no gas left within $r_g$. The middle row tells a similar story about the cold gas content (defined as gas with $T < 2 \times 10^4$ K), except for two important differences. The relationship between cold gas and halo mass for NoAGN galaxies has more scatter than the total gas mass to halo mass relation. Second, by $z=0$, none of the MrAGN galaxies has any cold gas left at all. Finally, the bottom panels show the stellar mass of the central galaxy versus the halo mass and reveal a more gradual separation of the two relations due to the permanent removal and/or heating of cold gas causing the buildup of stars in the MrAGN runs to lag behind that of the NoAGN runs. The solid black curves are

an estimate of the stellar-mass-halo-mass relation from the abundance matching analysis of Moster et al. (2013). The dashed black curve is the $z=0$ abundance matching estimate of Kravtsov et al. (2014). Our MrAGN galaxies are a very good match to the Kravtsov et al. (2014) estimate, which used improved photometric techniques and took into account intracluster light to measure stellar masses.

### 3.2.2. Inflow and Outflow Properties

Next, we turn to inflow and outflow properties. Again, the full inflow and outflow histories of our case study galaxies have been presented above; here, we present snapshots of these quantities at the specified redshifts for our entire sample of galaxies. Figure 12 illustrates that at $z=3$, the inflow rates for MrAGN and NoAGN galaxies are very similar and dependent on halo mass. By $z=2$, the higher-mass MrAGN galaxies begin to have inflow suppressed by $\gtrsim 1$ dex as their AGNs turn on. By $z=1$, this phenomenon is widespread, and the inflow rate is again broadly dependent on halo mass. At $z=0$, the highest-mass galaxies are again accreting at a rate comparable
<figure><img src="image_14.png" /><figcaption>**Figure 12.** Inflow rate of gas at the specified redshift. Top row: inflow across the galaxy radius. Bottom row: inflow across the halo radius. Black circles denote galaxies in the NoAGN runs, while red triangles denote galaxies in the MrAGN runs. The inflow rate at the galaxy radius for MrAGN galaxies is decreased by as much as 1.5 dex when $z=0$. We see evidence of downsizing in the initial decrease of inflow for high-mass galaxies at $z=2$ when their AGNs begin to switch on. At the halo radius, the inflow rate is only suppressed for lower-mass halos.</figcaption></figure>

<figure><img src="image_15.png" /><figcaption>**Figure 13.** Outflow rate of gas at the specified redshift. Top row: outflow across the galaxy radius. Bottom row: outflow across the halo radius. Black circles denote galaxies in the NoAGN runs, while red triangles denote galaxies in the MrAGN runs. The outflow rate at the galaxy radius for MrAGN galaxies is not appreciably different from that of NoAGN galaxies. At the halo radius, we see an elevated outflow rate for MrAGN galaxies, mainly at $z=2$ and $z=1$. We again see evidence for downsizing, with the most massive galaxies experiencing enhanced outflows at the halo radius before lower-mass galaxies.</figcaption></figure>

to their NoAGN counterparts, while the rest of the MrAGN galaxies continue to have their inflow rates suppressed by up to 1.5 dex. This again shows the importance of preventative feedback, especially in lower-mass galaxies.

At $r_h$, the evolution is very similar, although the suppression of inflow at this larger scale is much less pronounced. At $z=1$

and $z=0$, it is the lower-mass MrAGN galaxies that have their inflow suppressed (although now by ~0.5 dex), while more massive galaxies have inflow rates closer to their NoAGN counterparts. However, as shown in Table 1, the cumulative inflowing mass at $r_h$ expressed as a fraction of the universal baryon fraction times the halo mass at $z=0$ is smaller in
<figure><img src="image_16.png" /><figcaption>**Figure 14.** Outflow rate of gas across the galaxy radius $r_g$ vs. the bolometric luminosity of the galaxy's AGN. MrAGN galaxies are shown as red triangles. The dashed blue and green lines are scaling relations found for AGN-driven molecular and ionized outflows, respectively, and presented in Fiore et al. (2017). These scaling relations were derived for outflows observed at a number of different radii, and there is a large uncertainty on some of the values used to derive the relations, but our results appear to be in qualitative agreement with the observed scaling relations found in the literature and compiled in Fiore et al. (2017).</figcaption></figure>

MrAGN galaxies than in NoAGN galaxies, even at larger masses. Still, the effect is more pronounced at lower masses. At both scales, there is not much evolution in the inflow rates of NoAGN galaxies. At both radii, we see evidence of downsizing: AGN feedback occurs in more massive galaxies first, as evidenced by the initial decrease in inflow rate for the most massive galaxies at $z \sim 2$, and more strongly affects lower-mass galaxies at later times.

Figure 13 is the same as Figure 12, but now for the outflow rate. At $r_g$, the outflow rates of MrAGN galaxies at all redshifts are comparable to the outflow rates of NoAGN galaxies. At $r_h$, however, we see elevated outflow rates, first for more massive galaxies at $z=2$, then for most galaxies by $z=1$. The outflow rates of MrAGN galaxies remain slightly elevated by $z=0$ for many galaxies. This confirms what was shown earlier for the case studies: outflowing gas driven by AGN feedback is less likely to fall back into the galaxy and is more likely to travel out past the virial radius.

Figure 14 again shows the outflow rates of our MrAGN galaxies at four different redshifts, this time versus the bolometric luminosity of their AGNs. Higher-luminosity AGNs tend to be correlated with larger outflow rates. The blue and green dashed lines represent the scaling relations presented in Fiore et al. (2017) for molecular and ionized winds, respectively, for a compilation of observed AGN galaxies with detectable outflows. While the measurements of outflow velocity probe different radii in the galaxies used to define these relationships, we are encouraged that our MrAGN galaxies occupy a realistic portion of the OFR–$L_{bol}$ plane, and thus qualitatively agree with observed winds. Yesuf et al. (2017) have also observed winds in elliptical galaxies with AGNs at $z \sim 0.1$ and found velocities as large as $\sim 680 \text{ km s}^{-1}$, also in qualitative agreement with our results.

Figure 15 further quantifies the fraction of gas that has crossed a shell at $r_g$ that eventually crosses a shell at $r_h$ as a function of halo mass. Whenever a gas particle is tagged as outflowing past $r_g$, its radius is tracked, and if it later crosses $r_h$, it is considered expelled. At each redshift shown, we have divided the cumulative mass that has been expelled by the cumulative mass of gas that has been considered outflowing

past $r_g$ up until that redshift. These cumulative masses may include certain gas particles multiple times if they are recycled. At $z=3$, when galaxies are small, this ratio is similar for MrAGN galaxies and NoAGN galaxies. At $z=2$, we see the effect of AGN turning on in higher-mass galaxies as in Figures 12 and 13. At later redshifts, once AGN feedback kicks in for all galaxies, this ratio is much larger for MrAGN galaxies than for NoAGN galaxies, as gas particles ejected by AGN feedback tend to travel much farther. Interestingly, there does not appear to be much of a consistent trend with halo mass for this quantity.

In Figure 16, we examine different types of loading factors by plotting outflow rates versus inflow rates, both measured at $r_g$ (top row), SFR (middle row), and BHAR (bottom row). In the top row, we consider outflow rate versus inflow rate. At $z=3$, the two runs sit on top of each other at loading factors less than unity. By $z=2$, the loading factors of MrAGN galaxies tend to be larger by as much as 10 times, and by $z=1$, as much as a hundred times, with some galaxies having outflow rates larger than their inflow rates. This is due to a combination of both higher outflow rates and lower inflow rates for MrAGN galaxies. By $z=1$, almost every MrAGN galaxy has a loading factor greater than unity, while the NoAGN galaxies all sit below this line. Now, however, this difference is almost entirely due to the suppressed inflow rates of MrAGN galaxies. Finally, at $z=0$, it appears that outflow and inflow are regulating each other in the case of the MrAGN galaxies, while the NoAGN galaxies hit a floor in inflow rate, as stellar and SN feedback and gravitational heating do not suppress inflow as efficiently as AGN feedback does. Turning to outflow rate versus star formation rate, at $z=3, 2$, and $1$, this loading factor follows a similar pattern to the outflow rate/inflow rate. At $z=2$, the outflow rates of MrAGN galaxies are higher than those of NoAGN galaxies, while their SFRs are lower. At $z=1$, the difference in loading factor is due almost entirely to the decreased SFRs of MrAGN galaxies. At $z=0$, the one MrAGN galaxy still forming stars has a much lower SFR than any of the NoAGN galaxies, but also a very low outflow rate. Finally, the bottom panels show that the outflow rates of MrAGN galaxies are 100 times or more greater than their
<figure><img src="image_17.png" /><figcaption>**Figure 15.** Fraction of the cumulative mass which has crossed a shell at the galaxy radius and which has subsequently crossed a shell at the halo radius by the specified redshift. Black circles denote galaxies in the NoAGN runs, while red triangles denote galaxies in the MrAGN run. We see that a much larger fraction of gas that outflows past the galaxy radius eventually crosses the halo radius in the MrAGN runs.</figcaption></figure>

BHARs at $z=3$, although this value falls slightly with redshift. There are also fewer galaxies with appreciable BHARs as we go from $z=2$ to $z=0$.

In Figures 17 and 18, we examine, very broadly, the kinematics of outflowing material. Figure 17 shows the kinetic energy outflow rate for particles that have crossed a shell at $r_g$ since the last time step versus stellar mass. The kinetic energy outflow rate is calculated by summing the total kinetic energy of all of the outflowing particles and dividing by the time between snapshots. At $z=3$, the MrAGN galaxies and NoAGN galaxies are very similar, but by $z=2$, some of the MrAGN galaxies have much more kinetic energy in outflowing particles. By $z=1$, this trend is even clearer, and the MrAGN galaxies have started to become displaced to lower stellar masses as well. At $z=0$, the kinetic energy in outflowing particles is once again similar in the two runs, but the outflows have left their mark on the MrAGN galaxies, which have smaller stellar masses.

Figure 18 instead illustrates the average radial velocity of the same outflowing gas considered in Figure 17. At $z=3, 2$, and $1$, the trends are quite similar to those seen in the last figure. Large outflow velocities driven by AGN winds begin to be seen at $z=2$ and at $z=1$, and almost all MrAGN galaxies have outflowing material with significant radial velocities. At $z=0$, whereas the total kinetic energies of outflowing particles for MrAGN and NoAGN galaxies were very similar, the radial velocities of those outflowing gas particles are characteristically larger in the MrAGN run. The similarity in kinetic energies is due to the fact that there is less mass in outflowing gas in the MrAGN galaxies; much more gas has already been removed permanently by AGN feedback. In the NoAGN run, the same gas is allowed to outflow and reaccrete. The higher radial velocities of MrAGN galaxies at $z=0$ are mostly due to the higher gas temperature where outflows are roughly balanced by gas inflows with loading factor $\eta \sim 1$ (see Figure 16), except for a few galaxies with $\eta > 1$, which had AGN outburst in the recent past.

### 3.2.3. Gas Recycling

Figure 19 depicts the ratio of the inflow rate of recycled material at $r_g$ at the specified redshift to the inflow rate of new material. This is related to the inflowing recycling fraction depicted for our case studies. Down to $z=1$, galaxies in both runs are more likely to accrete a larger portion of new gas regardless of halo mass. This remains true for galaxies in the MrAGN run to $z=0$, but NoAGN galaxies experience more

recycled accretion, with higher-mass halos being dominated by recycling.

In Figure 20, we have taken the histograms like those shown in Figure 9 for all of our galaxies and turned them into cumulative timescale distributions. While we lose information on the relative number of recycling and ejection events, those properties are in line with those found for the example galaxies above. The cumulative distributions are color-coded by halo mass, with purple corresponding to the least massive and yellow to the most massive. MrAGN galaxies have characteristically longer recycling timescales due to the more efficient and energetic feedback. This is especially true for more massive galaxies; while the feedback is strong enough to lengthen recycling timescales, it is harder to remove the gas completely from a very massive halo, resulting in gas that is gone for a long time but nevertheless comes back. Although MrAGN galaxies show characteristically longer recycling timescales compared to NoAGN galaxies overall, we note that some low-mass MrAGN galaxies show cumulative distributions toward shorter recycling timescales compared to NoAGN galaxies. This is because gas recycling happens only at early epochs with short timescales, and the vast majority of gas is ejected for these galaxies. We also find that ejection timescales are characteristically longer for MrAGN galaxies; this time, the longest timescales correspond to the less-massive halos, which experienced strong early bouts of AGN feedback that removed gas permanently.

## 4. Discussion

In this section, we will first compare our results with those from other simulations, then discuss our results in the context of several questions asked in the Introduction.

### 4.1. Comparison with Other Works

Other studies have examined the baryon cycle in simulations, although usually these have been focused on lower-mass galaxies than the ones we study here. Despite this, we now discuss how our results compare with some of these studies.

Oppenheimer et al. (2010) examined a set of galaxies run with GADGET-2 (Springel 2005) with various prescriptions for stellar-driven winds, tied to the galaxies’ star formation rate. They find that when galactic winds are included, the contribution of recycled material to accretion at late times (especially $z<1$) is very significant. These outflow models result in too many galaxies with $M_* \gtrsim 2 \times 10^{12} M_\odot$ because
<figure><img src="image_18.png" /><figcaption>**Figure 16.** Outflow rate at the galaxy radius vs. the inflow rate at the galaxy radius (top row), star formation rate (middle row), and black hole accretion rate (bottom row, for MrAGN runs). Black circles denote galaxies in the NoAGN runs, while red triangles denote galaxies in the MrAGN run. Black lines specify constant mass-loading factors, $\eta$. Top row: $\eta$ = outflow rate/inflow rate. Middle row: $\eta$ = outflow rate/star formation rate. Bottom row: $\eta$ = outflow rate/black hole accretion rate. At $z=2$ and $z=1$, the MrAGN galaxies tend to have higher loading factors than the NoAGN galaxies, although by $z=0$, there is only one MrAGN galaxy that still exhibits star formation or black hole accretion.</figcaption></figure>

the material that is removed is reaccreted, and they posit that AGN feedback may suppress this reaccretion. These results are consistent with our large recycling fractions in the NoAGN runs; even if our sample is not cosmologically representative, we confirm their speculation that AGN feedback can suppress late reaccretion of gas.

Hirschmann et al. (2013) examine cosmological zoom simulations using the same stellar-driven wind model as Oppenheimer et al. (2010) and find the same large inflow rate at late times, leading to too much in situ star formation and galaxies that are too massive. They also found outflow rate versus SFR mass-loading factors at high redshift ($z \sim 2-3$) as large as $\sim 100$, which is higher than those we find in either of our MrAGN or NoAGN runs.

Übler et al. (2014) examined a series of cosmological zoom simulations, also run with GADGET, with masses slightly below our mass range and overlapping with our lower-mass galaxies. Much like this work, they compared two feedback models, one with weak feedback from massive stars, the other with strong feedback. They examined the baryon cycle around these galaxies and found that stronger stellar feedback produced disk galaxies with properties comparable to observations. Their strong feedback model is similar to our NoAGN model, and they find similar gas behavior. They find that at early times, some gas is ejected from the galaxies permanently, due to their shallower potential wells, while later the galaxies are dominated by the accretion of recycled gas, like ours. Their Figure 6 is the inspiration for our Figure 8.
<figure><img src="image_19.png" /><figcaption>**Figure 17.** Kinetic energy outflow rate of all gas particles that have crossed a shell at the galaxy radius since the last time step vs. stellar mass. The kinetic energy outflow rate is calculated by summing the total kinetic energy of all of the outflowing particles and dividing by the time between snapshots. Black circles denote galaxies in the NoAGN runs, while red triangles denote galaxies in the MrAGN run. At $z=3$ and $z=0$, the total kinetic energy of outflowing particles is very similar between the two runs, although the stellar masses of MrAGN galaxies are characteristically lower. At $z=2$ and $z=1$, MrAGN galaxies have more kinetic energy in outflows by as much as two orders of magnitude.</figcaption></figure>

<figure><img src="image_20.png" /><figcaption>**Figure 18.** Average radial velocity of all gas particles just crossing a shell at the galaxy radius at the specified redshift vs. stellar mass. Black circles denote galaxies in the NoAGN runs, while red triangles denote galaxies in the MrAGN run. For redshifts of 2 and less, the average radial velocity of gas outflowing past the galaxy radius in MrAGN galaxies is much higher than that of gas in NoAGN galaxies.</figcaption></figure>

<figure><img src="image_21.png" /><figcaption>**Figure 19.** Instantaneous inflow rate at the galaxy radius composed of recycled accretion divided by the inflow rate of new gas at the specified redshift. Black circles denote galaxies in the NoAGN runs, while red triangles denote galaxies in the MrAGN run. While MrAGN galaxies are almost always dominated by new accretion, NoAGN galaxies at $z=0$ experience more recycled accretion and are even dominated by recycled accretion in the large halo mass regime.</figcaption></figure>

Christensen et al. (2016) carried out a study of outflows in 20 field galaxies run with the SPH code GASOLINE (Wadsley et al. 2004). They looked at galaxies with halo masses of $10^{9.5}-10^{12} M_{\odot}$. They find that 50% of gas that leaves the galaxies is later reaccreted regardless of mass. Our results are in agreement with theirs in that we find that the future recycling fraction for our NoAGN galaxies is very similar for all galaxies, also regardless of mass. However, we also find that the cumulative inflowing mass of higher-mass galaxies is more likely to be recycled (~25%-35%) than that of lower-mass galaxies (as low as ~7%). Our

results also suggest that a very common recycling timescale for gas in our NoAGN galaxies is ~1 Gyr, in agreement with their results.

Barai et al. (2016) implemented mechanical AGN feedback similar to ours in a cluster-sized halo using GADGET-3. Although our analysis is quite different, they come to several similar conclusions. They find bipolar bubble-like outflows of heated gas out to hundreds of kiloparsecs, which is very similar to the structure of our outflows. They also find that kinetic feedback is far more efficient at affecting the central gas and quenching star formation than thermal feedback.
<figure><img src="image_22.png" /><figcaption>**Figure 20.** Normalized cumulative distribution of recycling (top row) and ejection (bottom row) event timescales by $z = 0$, color-coded by halo mass. Purple corresponds to lower halo mass, while yellow corresponds to higher halo mass. Left panels: MrAGN galaxies. Right panels: NoAGN galaxies. MrAGN galaxies have characteristically longer recycling timescales, especially for higher-mass galaxies. MrAGN galaxies also have characteristically longer ejection timescales, often the result of large, early bouts of AGN feedback. The latter trend is strongest for lower-mass galaxies.</figcaption></figure>

Anglés-Alcázar et al. (2017b) examined a suite of galaxies in the FIRE simulation (Hopkins et al. 2014), which spans a lower range of halo masses than the ones we have looked at and also do not include AGN feedback. They find (for galaxies with stellar feedback) that galaxies with lower halo masses are more dominated by recycled accretion and become more so at low redshift, while higher-mass halos are dominated by new accretion. We find that in general, NoAGN galaxies are more dominated by recycled accretion at late times, although this is much more prevalent for the higher-mass galaxies in our halo mass range. With our AGN feedback prescription, we find that there is very little recycled accretion onto our galaxies, even at late times, regardless of galaxy mass.

Anglés-Alcázar et al. (2017b) also find that in the absence of AGN feedback, the most massive galaxies have the shortest recycling timescales, although they do find that these same galaxies have larger recycling distances than lower-mass galaxies. We find that our NoAGN galaxies are all dominated by relatively short recycling timescales, and these timescales are longer for more massive halos. Recycled gas remains outside the galaxy for longer and travels farther in MrAGN galaxies than in NoAGN galaxies, and higher-mass halos have

recycled gas that spends the longest time outside the galaxy and travels the farthest. This is because gas launched from lower-mass galaxies is much more likely to be expelled from the halo.

Finally, Weinberger et al. (2017) implemented a prescription for AGN feedback into AREPO, the moving-mesh magneto-hydrodynamics code (Springel 2010). They included thermal feedback for high-Eddington ratio BH accretion and kinetic feedback for low-Eddington ratio accretion. Pillepich et al. (2018) adopted this dual-mode feedback prescription and also implemented an improved model for galactic-scale, star-formation-driven, kinetic winds in the IllustrisTNG simulation. They found that the interplay between the new AGN feedback prescription and the new galactic wind scalings resulted in realistic elliptical galaxies, with which we are in qualitative agreement, but it is unclear if they drive the powerful outflows associated with the observed quasars.

Our NoAGN model is in broad agreement with several of the studies above, which only include stellar and SN feedback. We find similar trends in terms of inflow rate at late times and gas recycling. However, differences in our feedback prescriptions cause us not to agree on every detail. While our MrAGN model produces realistic galaxies in agreement with other studies that
have included AGN feedback, our model is still unique in that we employ mechanical feedback at all BHARs.

## 4.2. Physical Interpretation

In this work, we have examined how the cycle of gas inflow and outflow is affected by our model for strong mechanical and radiation-driven AGN feedback. We have focused on 24 massive galaxies with halo masses of $M_{vir} \sim 10^{12}-10^{13.4} M_{\odot}$ at $z=0$. For each of these galaxies, we have runs from two different models: MrAGN and NoAGN. The MrAGN model includes stellar feedback via UV heating, stellar winds and SNe, AGN feedback via momentum-driven winds and X-ray heating, photoelectric heating, and cosmic X-ray background heating from a metagalactic X-ray background. The NoAGN model is identical except that it does not include any AGN feedback. The MrAGN model has been shown to produce realistic galaxy properties for massive galaxies (Choi et al. 2015, 2017; Hirschmann et al. 2017).

We set out, in part, to answer several questions about the gas cycle in our suite of galaxies, which we will discuss here.

### 4.2.1. What are the Histories of Inflow and Outflow Like for These Galaxies?

The inflow and outflow histories for both MrAGN and NoAGN galaxies are dependent upon halo mass (Figures 1–3, 12 and 13). As demonstrated by our case studies, as well as our ensemble plots, the inflow rate at both halo and galactic scales is halo mass dependent and is naturally due to the depth of the galaxy’s potential well. The outflow rates also appear to be partially mass dependent, especially at halo scales. Outflow rarely overtakes inflow at either scale of interest, except at late times in massive halos (like m0163) where inflow and outflow reach rough equilibrium due to the high central gas density. The outflows driven by AGN activity appear to correlate with AGN luminosity in a way qualitatively similar to winds in the universe (Figure 14).

In MrAGN galaxies, halo mass governs not only the cosmological inflow rate as for NoAGN galaxies, but also the effectiveness of AGN feedback. The outflow rates at $r_g$ for MrAGN galaxies are often comparable to those of NoAGN galaxies, except for when the AGN turns on, resulting in a spike in the outflow rate. From $z \sim 2$ to $z \sim 1$, the outflow rate at $r_h$ is elevated in MrAGN galaxies relative to NoAGN galaxies, as material driven by AGN feedback travels farther than material driven by stellar and SN feedback. We also see inflow suppressed at both radii significantly and for almost all galaxies at $r_g$; after the AGN turns on, inflow and outflow tend to track each other (Figure 16, top row). The effect is more subtle and mainly significant for lower-mass galaxies at $r_h$.

Both outflow enhancement and inflow suppression are seen first in higher-mass galaxies at high redshift before manifesting in lower-mass galaxies at later times, a result that corresponds to the phenomenon of “downsizing,” or antihierarchical BH growth (Hirschmann et al. 2012, 2014). In our simulations, this is a result of the more massive halos receiving seed BHs at earlier times when the gas density in the BH’s vicinity is very high, resulting in massive BHs at higher redshift feeding voraciously. This in turn results in more massive galaxies being affected by AGN feedback earlier. This is in agreement with many observational studies of the AGN population, which conclude that the number density of the most-luminous AGN

peaks at high redshift, while less-luminous AGNs have almost constant number density and are more prevalent than high-luminosity AGNs at late times (Cristiani et al. 2004; Croom et al. 2004; Matute et al. 2006).

We also want to emphasize the increasing importance of so-called “preventative” feedback as we move to lower-mass halos. While “ejective” feedback is important for removing gas and bouts of ejective feedback are seen in halos of all masses, the inflow of gas, both new and recycled, is far more suppressed in low-mass halos than in high-mass halos. By $z=1-0$, MrAGN galaxies have inflow at galactic scales suppressed by as much as $\sim 1.6$ dex, with lower-mass halos most strongly affected. In fact, the lowest-mass halos even have inflow suppressed at halo scales by up to $\sim 0.5$ dex. The cumulative mass of inflowing material onto the galaxy is lower by 0.4–0.6 dex for our MrAGN case studies, while the cumulative mass accreted onto the halo is suppressed by 0.1–0.3 dex. This is a result of outflowing material disrupting the infall of new and recycled material by imparting momentum and energy. The infall is easier to halt in a shallower potential well. This mass-dependent preventative feedback is an interesting consequence of our mechanical feedback prescription. The strong outflows and their ability to clear out gas from the galaxy are not in and of themselves a surprise as the feedback prescription is designed to launch these winds, but their ability to sweep up gas that is on its way into the galaxy and in some cases turn it around is an interesting counterpoint to the usual theoretical method of putting in two separate feedback prescriptions, a “wind” mode and a “maintenance” mode; they are rather two sides of the same coin. Even if our AGN feedback model might not be perfect, since gas-free present-day galaxies (such as m0501) may not be realistic, the ability of strong outflows to act in both ejective and preventative ways can be considered a general result and should be kept in mind.

### 4.2.2. How Much Gas is Removed Permanently and How Much Comes Back?

The main difference between the gas cycles in MrAGN and NoAGN galaxies is that while NoAGN galaxies remain recycling dominated throughout their lives, MrAGN galaxies become ejection dominated when their BHs begin feeding and their AGNs turn on (Figure 7). While as much as 90% of the material in a given outflow episode might return to a NoAGN galaxy, this fraction is much smaller, at most 20%, for MrAGN galaxies. Even fractions this high are only seen in more massive galaxies whose potential wells are deep enough to retain some of the gas pushed out by AGN feedback. AGN-driven outflowing gas can travel longer distances than stellar and SN-driven gas because it is launched with much higher velocities and is harder to slow down and turn around (Figures 17 and 18). MrAGN galaxies can have as much as 90% of the total outflowing material that is driven out past $r_g$ cross the virial radius by $z=0$ (Figure 15). NoAGN galaxies tend to have much smaller fractions of expelled material, around 1%–2%. This results in NoAGN galaxies having much higher contributions of recycled gas to their inflow and outflow by $z=0$, with larger contributions for larger halos (Figures 7 and 19).

Since gas in NoAGN galaxies is more likely to be recycled, each gas particle is more likely to accrete two or more times than gas particles in MrAGN galaxies, and the timescales of
these recycling events tend to be shorter, and the distances traveled shorter, than those for MrAGN gas particles (Figures 8, 9, and 20). Again, this seems to be due to the higher velocities associated with outflows driven by AGN feedback, which causes the recycling history of galaxies in the two runs to diverge when the BH begins to feed. After that point, MrAGN galaxies rarely recycle outflowing material. The large distances traveled by gas particles may have interesting implications for metal enrichment in the hot gas halo around galaxies and in the IGM (R. Brennan et al. 2018, in preparation).

### 4.2.3. How are the Host Galaxies Affected?

Galaxies in our two runs tend to have different gas morphologies by $z=0$ as a result of the different effects of feedback. In higher-mass halos, the difference is mainly in the gas density and more diffuse hot gas halo around MrAGN galaxies as a result of more gas being removed from both the galaxy's inner regions and its gas halo (Figures 4 and 5). In lower-mass galaxies, however, the difference is much more stark. Lower-mass MrAGN galaxies may be left only with diffuse hot gas (or may be depleted of gas completely), while their NoAGN counterparts develop substantial disks of cold gas by $z=0$ (Figure 6). While the complete removal of gas may be too extreme, the cold gas disk, as well as the young stellar disk that accompanies it, are in conflict with observations in the mass range of our sample.

Beyond morphologies, AGN feedback is needed to quench our galaxies in order for them to resemble the observed high-mass galaxy population. Our feedback prescription accomplishes this by removing or heating the cold gas in the galaxy, which results in a steep decrease in star formation (Figures 1–3). In the case of lower-mass galaxies, as mentioned above, even the hot gas within $r_h$ is severely depleted by feedback (in m0501, this decrease is >2 dex). This results in overall reduced in situ SFRs and thus smaller gas (by ~0.5–2 dex) and stellar (by ~0.2 dex) masses in MrAGN galaxies by $z=0$ (Figure 11).

The effects of our implementation of AGN feedback bring our simulated galaxies into better agreement with observed galaxies in terms of morphology, star formation rate, and stellar mass, but it may result in cold gas fractions that are too low.

## 5. Summary

Our study of two sets of 24 cosmological zoom galaxies both with and without mechanical and radiation-driven AGN feedback has led to the following main conclusions:

1. Our model for radiation-pressure driven AGN feedback enhances galaxy-scale outflows and acts in an ejective manner. Outflow rates at $r_g$ are comparable between the two runs, but outflow rates at $r_h$ can be enhanced by up to 1 dex in MrAGN galaxies at $z \sim 1-2$. This is the result of higher outflowing gas velocities in the MrAGN runs, which cause a higher fraction of outflowing material at $r_g$ to escape the galaxy's potential well and cross $r_h$. This fraction is as high as 80% in MrAGN galaxies, compared with ~5%–10% in NoAGN galaxies.

2. Our feedback model also suppresses inflow in MrAGN galaxies relative to NoAGN galaxies, affecting feedback in a preventative way as well, especially at the low-mass end of our sample. At $r_g$, the inflow rate can be

suppressed by as much as 1.5 dex, while at $r_h$, inflow rates can be suppressed by up to 0.5 dex in lower-mass galaxies. This results in an overall smaller cumulative inflowing mass relative to the final halo masses of MrAGN galaxies versus their NoAGN counterparts (see Table 1 for these numbers for our case study halos).

3. Our NoAGN galaxies are recycling dominated throughout their lives, such that most of the material removed from the galaxy returns (~68%–78% for our case studies). By $z=0$, ~50% of material inflowing into NoAGN galaxies is returning recycled material. Once their BHs begin to feed, MrAGN galaxies become ejection dominated, with the majority of outflowing gas never returning to the galaxy (only about 8%–15% returns in our case studies). The recycled inflowing fraction of MrAGN galaxies at $z=0$ can be as low as a few percent.

4. Gas that is recycled in MrAGN galaxies tends to remain outside of the galaxy for longer (~2–3 Gyr versus 1 Gyr) and to travel farther than recycled gas in NoAGN galaxies (up to several megaparsecs). Accreted gas is more likely to undergo several recycling events in NoAGN galaxies than in MrAGN galaxies and is more likely to remain in the galaxy (either as gas or stars) until $z=0$. Between ~70% and 90% of accreted gas remains in our NoAGN case studies at $z=0$ compared with ~30%–50% in our two MrAGN case studies that have entire histories.

5. Our model for AGN feedback associated with radiatively efficient BH accretion succeeds in quenching massive galaxies over long timescales and keeping disks from re-forming.

There are several further avenues of study, which we plan to investigate in a series of future works.

One interesting point illustrated by our individual galaxy histories (Figures 1–3) is the weak correlation between AGN activity/outflow events and galaxy mergers. This is different from what has been seen in some other studies of simulations with AGN feedback, such as Tremmel et al. (2017). While it does sometimes appear that mergers and outflow events may be correlated, such as for m0163, it is clear that a major merger is not a prerequisite for quenching in our simulations. It is worth noting that the green dashed lines in Figures 1–3 represent 1:10 or greater halo mass ratio mergers. The stellar or baryonic mass ratios and the time of the eventual galaxy mergers will be somewhat different, and this could contribute to the lack of correlation seen. It is also possible, especially in the early universe, to have several smaller satellites, none of which representing a 1:10 merger, interact with the central galaxy, the collective effect of which is disruption on par with a more major interaction. AGN activity may also be triggered by a very rich gas supply along filaments that does not take the form of galaxy interactions. This lack of a strong relationship between merger events and AGN activity is also being seen in other simulations of statistically complete galaxy populations, such as the Magneticum simulations (Steinborn et al. 2018). In any event, a much closer look at how this AGN activity correlates with merger history is reserved for a future work (E. Choi et al. 2018, in preparation).

We also plan to carry out a detailed comparison with observations on several fronts. We plan to do a more in depth study of the kinematics of these outflows beyond the very
superficial discussion here, as well as the kinematics of the gas within the galaxy. This will include a theoretical study of the angular momentum, as well as a study of mock absorption lines to compare the kinematics of our simulated winds with those observed in the universe. We also plan to use these mock observations to compare with observed mass-loading factors (again, only briefly touched on here), in order to see how realistic our outflows appear (especially in the case of our more extreme blowout events in lower-mass galaxies) and how our interpretation of these events might change from different viewing angles and at different times.

The very large distances traveled by ejected gas (and even recycled gas in the MrAGN runs) has interesting implications for outflows enriching the surrounding halo and the IGM, as well as enriching other galaxies at very large distances, something which we have not studied here. We will present a detailed study of the impact of AGN-driven winds on metal enrichment in R. Brennan et al. (2018, in preparation). Finally, we also plan to explore the impact of AGN-driven outflows on nebular emission and absorption in different regions of a galaxy, producing spatially resolved emission and absorption line maps to improve the interpretation of modern integral field spectroscopic observations (M. Hirschmann et al. 2018, in preparation).

In this paper, we have presented a first theoretical look at the inflow and outflow properties in a series of cosmological zoom simulations in the presence and absence of mechanical and radiation-driven AGN feedback. This is the first step in a series of studies in which we will study these properties from the perspective of an observer in an effort to compare the gas cycle in these zoom simulations with that which can be inferred from real galaxies.

We are grateful to the anonymous referee for constructive comments that improved the paper. R.B. was supported in part by HST Theory grant HST-AR-13270-A. E.C. was supported by NASA through the HST Cycle 23 AR-14287 grant. R.S.S. thanks the Downsborough family for their generous support and gratefully acknowledges support from the Simons Foundation through a Simons Investigator grant. M.H. acknowledges financial support from the European Research Council (ERC) via an Advanced Grant under grant agreement No. 321323-NEOGAL. R.B. would like to thank David Lynch and Mark Frost for the new *Twin Peaks* and Damon Lindelof and Tom Perrotta for the final season of *The Leftovers*.

## ORCID iDs

Ena Choi https://orcid.org/0000-0002-8131-6378

Michaela Hirschmann https://orcid.org/0000-0002-3301-3321

Thorsten Naab https://orcid.org/0000-0002-7314-2558

Jeremiah P. Ostriker https://orcid.org/0000-0002-6405-9904

## References

* Alatalo, K., Blitz, L., Young, L. M., et al. 2011, ApJ, 735, 88

* Anglés-Alcázar, D., Davé, R., Faucher-Giguère, C.-A., Özel, F., & Hopkins, P. F. 2017a, MNRAS, 464, 2840

* Anglés-Alcázar, D., Faucher-Giguère, C.-A., Kereš, D., et al. 2017b, MNRAS, 470, 4698

* Arav, N., Borguet, B., Chamberlain, C., Edmonds, D., & Danforth, C. 2013, MNRAS, 436, 3286

* Arav, N., Edmonds, D., Borguet, B., et al. 2012, A&A, 544, A33

* Aumer, M., White, S. D. M., Naab, T., & Scannapieco, C. 2013, MNRAS, 434, 3142

* Baldry, I. K., Glazebrook, K., Brinkmann, J., et al. 2004, ApJ, 600, 681

* Barai, P., Murante, G., Borgani, S., et al. 2016, MNRAS, 461, 1548

* Behroozi, P. S., Wechsler, R. H., & Wu, H.-Y. 2013, ApJ, 762, 109

* Bell, E. F., Wolf, C., Meisenheimer, K., et al. 2004, ApJ, 608, 752

* Bogdán, Á., Vogelsberger, M., Kraft, R. P., et al. 2015, ApJ, 804, 72

* Bondi, H. 1952, MNRAS, 112, 195

* Chamberlain, C., & Arav, N. 2015, MNRAS, 454, 675

* Cheung, E., Bundy, K., Cappellari, M., et al. 2016, Natur, 533, 504

* Choi, E., Naab, T., Ostriker, J. P., Johansson, P. H., & Moster, B. P. 2014, MNRAS, 442, 440

* Choi, E., Ostriker, J. P., Naab, T., et al. 2017, ApJ, 844, 31

* Choi, E., Ostriker, J. P., Naab, T., & Johansson, P. H. 2012, ApJ, 754, 125

* Choi, E., Ostriker, J. P., Naab, T., Oser, L., & Moster, B. P. 2015, MNRAS, 449, 4105

* Christensen, C. R., Davé, R., Governato, F., et al. 2016, ApJ, 824, 57

* Cicone, C., Maiolino, R., Sturm, E., et al. 2014, A&A, 562, A21

* Crenshaw, D. M., Kraemer, S. B., & George, I. M. 2003, ARA&A, 41, 117

* Cristiani, S., Alexander, D. M., Bauer, F., et al. 2004, ApJL, 600, L119

* Croom, S. M., Smith, R. J., Boyle, B. J., et al. 2004, MNRAS, 349, 1397

* Cullen, L., & Dehnen, W. 2010, MNRAS, 408, 669

* Debuhr, J., Quataert, E., & Ma, C.-P. 2011, MNRAS, 412, 1341

* Dehnen, W., & Aly, H. 2012, MNRAS, 425, 1068

* Di Matteo, T., Khandai, N., DeGraf, C., et al. 2012, ApJL, 745, L29

* Dubois, Y., Gavazzi, R., Peirani, S., & Silk, J. 2013, MNRAS, 433, 3297

* Dunn, J. P., Bautista, M., Arav, N., et al. 2010, ApJ, 709, 611

* Durier, F., & Dalla Vecchia, C. 2012, MNRAS, 419, 465

* Eisenreich, M., Naab, T., Choi, E., Ostriker, J. P., & Emsellem, E. 2017, MNRAS, 468, 751

* Faber, S. M., Willmer, C. N. A., Wolf, C., et al. 2007, ApJ, 665, 265

* Fabian, A. C. 2012, ARA&A, 50, 455

* Feain, I. J., Papadopoulos, P. P., Ekers, R. D., & Middelberg, E. 2007, ApJ, 662, 872

* Ferrarese, L., & Merritt, D. 2000, ApJL, 539, L9

* Fiore, F., Feruglio, C., Shankar, F., et al. 2017, A&A, 601, A143

* Ganguly, R., & Brotherton, M. S. 2008, ApJ, 672, 102

* Gaskell, C. M., Gill, J. J. M., & Singh, J. 2016, arXiv:1611.03733

* Gebhardt, K., Bender, R., Bower, G., et al. 2000, ApJL, 539, L13

* Genzel, R., Förster Schreiber, N. M., Rosario, D., et al. 2014, ApJ, 796, 7

* Haardt, F., & Madau, P. 2012, ApJ, 746, 125

* Häring, N., & Rix, H.-W. 2004, ApJL, 604, L89

* Harrison, C. M. 2017, NatAs, 1, 0165

* Heckman, T. M., & Best, P. N. 2014, ARA&A, 52, 589

* Higginbottom, N., Proga, D., Knigge, C., et al. 2014, ApJ, 789, 19

* Hirschmann, M., Charlot, S., Feltre, A., et al. 2017, MNRAS, 472, 2468

* Hirschmann, M., Dolag, K., Saro, A., et al. 2014, MNRAS, 442, 2304

* Hirschmann, M., Naab, T., Davé, R., et al. 2013, MNRAS, 436, 2929

* Hirschmann, M., Somerville, R. S., Naab, T., & Burkert, A. 2012, MNRAS, 426, 237

* Hopkins, P. F. 2013, MNRAS, 428, 2840

* Hopkins, P. F., Kereš, D., Oñorbe, J., et al. 2014, MNRAS, 445, 581

* Hopkins, P. F., Torrey, P., Faucher-Giguère, C.-A., Quataert, E., & Murray, N. 2016, MNRAS, 458, 816

* Hu, C.-Y., Naab, T., Walch, S., Moster, B. P., & Oser, L. 2014, MNRAS, 443, 1173

* Khandai, N., Di Matteo, T., Croft, R., et al. 2015, MNRAS, 450, 1349

* Knigge, C., Scaringi, S., Goad, M. R., & Cottis, C. E. 2008, MNRAS, 386, 1426

* Kormendy, J., & Richstone, D. 1995, ARA&A, 33, 581

* Kravtsov, A., Vikhlinin, A., & Meshscheryakov, A. 2014, arXiv:1401.7329

* Lilly, S. J., Carollo, C. M., Pipino, A., Renzini, A., & Peng, Y. 2013, ApJ, 772, 119

* Lynden-Bell, D. 1969, Natur, 223, 690

* Magorrian, J., Tremaine, S., Richstone, D., et al. 1998, AJ, 115, 2285

* Maiolino, R., Gallerani, S., Neri, R., et al. 2012, MNRAS, 425, L66

* Marconi, A., & Hunt, L. K. 2003, ApJL, 589, L21

* Matute, I., La Franca, F., Pozzi, F., et al. 2006, A&A, 451, 443

* Moe, M., Arav, N., Bautista, M. A., & Korista, K. T. 2009, ApJ, 706, 525

* Moster, B. P., Naab, T., & White, S. D. M. 2013, MNRAS, 428, 3121

* Naab, T., & Ostriker, J. P. 2017, ARA&A, 55, 59

* Núñez, A., Ostriker, J. P., Naab, T., et al. 2017, ApJ, 836, 204

* Nyman, L.-A., Booth, R. S., Carlstrom, U., et al. 1992, A&AS, 93, 121

* Oppenheimer, B. D., Davé, R., Kereš, D., et al. 2010, MNRAS, 406, 2325

* Oser, L., Naab, T., Ostriker, J. P., & Johansson, P. H. 2012, ApJ, 744, 63
Oser, L., Ostriker, J. P., Naab, T., Johansson, P. H., & Burkert, A. 2010, ApJ, 725, 2312

Ostriker, J. P., Choi, E., Ciotti, L., Novak, G. S., & Proga, D. 2010, ApJ, 722, 642

Page, M. J., Symeonidis, M., Vieira, J. D., et al. 2012, Natur, 485, 213

Pandya, V., Brennan, R., Somerville, R. S., et al. 2017, MNRAS, 472, 2054

Peng, Y., Maiolino, R., & Cochrane, R. 2015, Natur, 521, 192

Perna, M., Brusa, M., Salvato, M., et al. 2015, A&A, 583, A72

Pillepich, A., Springel, V., Nelson, D., et al. 2018, MNRAS, 473, 4077

Prochaska, J. X., & Hennawi, J. F. 2009, ApJ, 690, 1558

Proga, D., & Kurosawa, R. 2009, in AIP Conf. Ser. 1171, Recent Directions in Astrophysical Quantitative Spectroscopy and Radiation Hydrodynamics: Proc. of the International Conf. in Honor of Dimitri Mihalas for His Lifetime Scientific Contributions on the Occasion of His 70th Birthday, ed. I. Hubeny et al. (Melville, NY: AIP), 295

Read, J. I., & Hayfield, T. 2012, MNRAS, 422, 3037

Rosario, D. J., Santini, P., Lutz, D., et al. 2013a, ApJ, 771, 63

Rosario, D. J., Trakhtenbrot, B., Lutz, D., et al. 2013b, A&A, 560, A72

Röttgers, B. 2017, Astrophysics Source Code Library, pygad: Analysing Gadget Simulations with Python, ascl:1569

Rupke, D. S. N., & Veilleux, S. 2011, ApJL, 729, L27

Rupke, D. S. N., & Veilleux, S. 2013, ApJ, 768, 75

Saitoh, T. R., & Makino, J. 2009, ApJL, 697, L99

Santini, P., Rosario, D. J., Shao, L., et al. 2012, A&A, 540, A109

Sazonov, S. Y., Ostriker, J. P., Ciotti, L., & Sunyaev, R. A. 2005, MNRAS, 358, 168

Sazonov, S. Y., Ostriker, J. P., & Sunyaev, R. A. 2004, MNRAS, 347, 144

Schaye, J., Crain, R. A., Bower, R. G., et al. 2015, MNRAS, 446, 521

Shakura, N. I., & Sunyaev, R. A. 1973, A&A, 24, 337

Somerville, R. S., & Davé, R. 2015, ARA&A, 53, 51

Somerville, R. S., Hopkins, P. F., Cox, T. J., Robertson, B. E., & Hernquist, L. 2008, MNRAS, 391, 481

Spergel, D. N., Bean, R., Doré, O., et al. 2007, ApJS, 170, 377

Springel, V. 2005, MNRAS, 364, 1105

Springel, V. 2010, MNRAS, 401, 791

Springel, V., Di Matteo, T., & Hernquist, L. 2005, MNRAS, 361, 776

Steinborn, L. K., Dolag, K., Hirschmann, M., Prieto, M. A., & Remus, R.-S. 2015, MNRAS, 448, 1504

Steinborn, L. K., Hirschmann, M., Dolag, K., et al. 2018, MNRAS, submitted (arXiv:1805.06902)

Sturm, E., González-Alfonso, E., Veilleux, S., et al. 2011, ApJL, 733, L16

Tombesi, F., Cappi, M., Reeves, J. N., et al. 2010, A&A, 521, A57

Tombesi, F., Cappi, M., Reeves, J. N., et al. 2013, MNRAS, 430, 1102

Tombesi, F., Meléndez, M., Veilleux, S., et al. 2015, Natur, 519, 436

Torresi, E., Grandi, P., Costantini, E., & Palumbo, G. G. C. 2012, in Int. Journal of Modern Physics Conf. Ser. 8 (World Scientific Publishing Co.), 396

Tremaine, S., Gebhardt, K., Bender, R., et al. 2002, ApJ, 574, 740

Tremmel, M., Karcher, M., Governato, F., et al. 2017, MNRAS, 470, 1121

Tremonti, C. A., Moustakas, J., & Diamond-Stanic, A. M. 2007, ApJL, 663, L77

Übler, H., Naab, T., Oser, L., et al. 2014, MNRAS, 443, 2092

Vogelsberger, M., Genel, S., Springel, V., et al. 2014, MNRAS, 444, 1518

Wadsley, J. W., Stadel, J., & Quinn, T. 2004, NewA, 9, 137

Weinberger, R., Springel, V., Hernquist, L., et al. 2017, MNRAS, 465, 3291

Yesuf, H. M., Faber, S. M., Koo, D. C., & Lee, L. L. 2017, arXiv:1706.09500

Zakamska, N. L., & Greene, J. E. 2014, MNRAS, 442, 784

Zinn, P.-C., Middelberg, E., Norris, R. P., & Dettmar, R.-J. 2013, ApJ, 774, 66