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

THE BARYON CYCLE AT HIGH REDSHIFTS: EFFECTS OF GALACTIC WINDS ON GALAXY EVOLUTION IN OVERDENSE AND AVERAGE REGIONS

RAPHAEL SADOUN¹,²,³, ISAAC SHLOSMAN²,⁴, JUN-HWAN CHOI²,⁵, AND EMILIO ROMANO-DÍAZ²,⁶

¹Department of Physics & Astronomy, University of Utah, Salt Lake City, UT 84112-0830, USA; raphael.sadoun@utah.edu

²Department of Physics & Astronomy, University of Kentucky, Lexington, KY 40506-0055, USA

³Institut d'Astrophysique de Paris (UMR 7095), F-75014 Paris, France

⁴Theoretical Astrophysics, Department of Earth & Space Science, Osaka University, Osaka 560-0043, Japan

⁵Department of Astronomy, University of Texas, Austin, TX 78712, USA

⁶Argelander Institut fuer Astronomie, Auf dem Haegel 71, D-53121 Bonn, Germany

Received 2016 April 14; revised 2016 June 22; accepted 2016 July 18; published 2016 September 23

ABSTRACT

We employ high-resolution cosmological zoom-in simulations focusing on a high-sigma peak and an average cosmological field at $z \sim 6-12$ in order to investigate the influence of environment and baryonic feedback on galaxy evolution in the reionization epoch. Strong feedback, e.g., galactic winds, caused by elevated star formation rates (SFRs) is expected to play an important role in this evolution. We compare different outflow prescriptions: (i) constant wind velocity (CW), (ii) variable wind scaling with galaxy properties (VW), and (iii) no outflows (NW). The overdensity leads to accelerated evolution of dark matter and baryonic structures, absent from the “normal” region, and to shallow galaxy stellar mass functions at the low-mass end. Although CW shows little dependence on the environment, the more physically motivated VW model does exhibit this effect. In addition, VW can reproduce the observed specific SFR (sSFR) and the sSFR-stellar mass relation, which CW and NW fail to satisfy simultaneously. Winds also differ substantially in affecting the state of the intergalactic medium (IGM). The difference lies in the volume-filling factor of hot, high-metallicity gas, which is near unity for CW, while such gas remains confined in massive filaments for VW, and locked up in galaxies for NW. Such gas is nearly absent from the normal region. Although all wind models suffer from deficiencies, the VW model seems to be promising in correlating the outflow properties with those of host galaxies. Further constraints on the state of the IGM at high $z$ are needed to separate different wind models.

**Key words:** cosmology: theory – dark ages, reionization, first stars – galaxies: formation – galaxies: high-redshift – methods: numerical

1. INTRODUCTION

With the rapidly increasing number of available multi-wavelength observations as well as the development of more sophisticated theoretical and numerical tools, our understanding of galaxy formation and evolution has progressed considerably in recent years. While the majority of large galaxy surveys have focused on low and intermediate redshifts, substantial data are now available to gain insight into the early stages of galaxy formation at $z > 6$. In particular, observations of high-redshift quasars (QSOs) at $z \gtrsim 6$ have hinted at the presence of rare, overdense regions hosting massive structures less than one billion years after the Big Bang (Fan et al. 2001). Deep observations with the wide-area Suprime-Cam of the Subaru telescope have led to the discovery of $\sim 3\sigma$ overdensities of Lyman break galaxies (LBGs) on scales of a few megaparsecs around the $z \sim 6.4$ QSO CFHQS J2329-0301 (Utsumi et al. 2010). Similar results have been obtained using the Large Binocular Telescope to study the environment of other high-$z$ QSOs and also concluded that they reside in overdense regions (e.g., Morselli et al. 2014). Recent theoretical work has shown that these regions do not necessarily evolve into the most massive clusters at $z=0$ and argued that galaxy formation might proceed differently there than in fields of average density (Trenti & Stiavelli 2008; Overzier et al. 2009; Romano-Díaz et al. 2011a, 2011b).

Star formation and the activity of active galactic nuclei (AGNs) both appear to peak around $z \sim 2-3$ (e.g., Madau & Dickinson 2014), and galactic outflows, or galactic winds, are

expected to exhibit a similar behavior when averaged over a large enough volume. However, overdense regions, such as those hosting high-$z$ QSOs, have been predicted to evolve ahead of the average and underdense regions (e.g., Barkana & Loeb 2004; Romano-Díaz et al. 2011a, 2011b, 2014; Yajima et al. 2015). Although their evolution is concurrent with the reionization epoch at $6 \lesssim z \lesssim 14$, they can be reionized ahead of the normal regions in the universe because of the elevated star formation rates (SFRs) there. Moreover, these regions should collapse early, and, because the universe is also substantially denser at these redshifts, one expects mass accretion rates to be high. The question then is to what degree the associated galactic outflows from these objects affect their evolution as well as the evolution of their environment. In this paper, we aim to study the effects of these winds on galaxy evolution using high-resolution cosmological simulations of overdense regions and compare them with the evolution in regions of average density.

Massive outflows are projected to play an important role in the evolution of a large number of astrophysical objects—from single and binary stars to supermassive black holes (SMBHs) and starburst galaxies. On galactic scales, winds are expected to affect the evolution of their hosts by regulating the SFR (e.g., Oppenheimer & Davé 2006; Scannapieco et al. 2006; Schaye et al. 2010), ejecting gas and metals into the parent dark matter (DM) halos and the intergalactic medium (IGM, e.g., Tremonti et al. 2004; Davé et al. 2008; Kirby et al. 2008; Sijacki et al. 2009), and determining the physical properties (disk mass and size, chemical evolution, galaxy luminosity function, and
cuspiness of DM halos) of these objects (e.g., Dekel & Silk 1986; El-Zant et al. 2001, 2004; Romano-Díaz et al. 2008). These outflows powered by supernovae (SNe) can be supplemented by feedback from the accretion processes onto the central SMBHs, the so-called AGN winds (e.g., Blandford & Payne 1982; Shlosman et al. 1985; Emmering et al. 1992; Konigl & Kartje 1994; Murray et al. 1995; Arav et al. 1997, p. 128; Cecil et al. 2001; Proga 2003), for which there exists plentiful observational evidence (e.g., Emonts et al. 2005; Croston et al. 2008; Cano-Díaz et al. 2012; Combes et al. 2013; Cicone et al. 2014). A compelling example that has demonstrated the necessity for such strong feedback can be found in the comparison study of disk evolution in a cosmological setting, with and without feedback (e.g., Robertson et al. 2004; Schaye et al. 2015), where in the absence of feedback the gas quickly violates the Toomre criterion and fragments. Without feedback there is an overproduction of metals, especially in small galaxies. Even more revealing is the overcooling problem, which leads directly to the angular momentum catastrophe—the gas falls into the subhalos, and hitchhikes to the bottom of the potential well of the parent galaxy (e.g., Maller & Dekel 2002), possibly contributing to the overgrowth of classical bulges.

The presence of galactic outflows is supported by numerous observations in both low- and high-redshift galaxies (e.g., Osterbrock 1960; Lynds & Sandage 1963; Mathews & Baker 1971; Heckman et al. 1990; Hummel et al. 1991; Heckman 1994, p. 234; Kunth et al. 1998), including the $z \sim 3-4$ LBGs (e.g., Pettini et al. 2000; Adelberger et al. 2003; Shapley et al. 2003), gravitationally lensed galaxies at $z \sim 4-5$ (e.g., Franx et al. 1997), luminous and ultra-luminous IR galaxies (ULIRGs) (e.g., Smail et al. 2003), etc. About 75%-100% of ULIRGs are associated with such winds—their strength appears to correlate with the SFRs but saturates for the most luminous ULIRGs (e.g., Martin 2005; Veilleux et al. 2005). Galactic winds are typically characterized by a biconical symmetry with respect to the underlying galactic disks (e.g., Veilleux et al. 1994; Shopbell & Bland-Hawthorn 1998). Hot tenuous winds show evidence for associated cold ($\sim 10^4$ K) and neutral clouds (e.g., Heckman 1994, p. 234; Stewart et al. 2000)—so they represent a multiphase interstellar medium (ISM). An important link with theoretical predictions has been made by recent observations revealing that outflow velocities appear to correlate with the galaxy's stellar mass or its SFR (e.g., Martin 2005).

A number of driving mechanisms have been proposed for galactic winds, such as radiation pressure in the UV lines and on dust, thermal pressure from the SN and OB stellar heating (e.g Larson 1974; Castor et al. 1975; Dekel & Silk 1986; Mac Low & McCray 1988; Ostriker & McKee 1988; Shull & Saken 1995; Costa et al. 2014; Vogelsberger et al. 2014), and pressure from cosmic rays (e.g., Uhlig et al. 2012). The contribution of AGNs is debatable at present. The driving of the outflows by the SN and stellar winds occurs when ejecta of individual sources form a bubble of hot gas, at $\sim 10^{7-8}$ K, which expands due to the strong overpressure down the steepest pressure gradient and enters the “blow-out” stage. Besides the mass, energy, and momentum injected by such winds, this is probably the main way for the highly enriched material to be placed into the halo and further out into the IGM (e.g., Cen et al. 2005).

Galactic bulges may provide a testing ground for our understanding of various mechanisms that regulate the star formation and angular momentum redistribution in a forming disk galaxy. An unexpectedly large fraction, ~76%, of massive ($\gtrsim 10^{10} M_\odot$) galactic disks can be fit with the Sérsic index $\lesssim 2$, and ~69% have a bulge-to-total disk mass (B/T) ratio of $\lesssim 0.2$, both in barred and unbarred galaxies (e.g., Weinzirl et al. 2009). Such low Sérsic indexes represent disky rather than spheroidal stellar distribution (e.g., Kormendy & Kennicutt 2004). On the other hand, bulges obtained in numerical simulations are dominated by massive spheroids. This contradiction becomes even stronger in the recent study of the bulge population within a sphere of 11 Mpc radius around the Milky Way using Spitzer 3.6 µm and Hubble Space Telescope (HST) data (e.g., Fisher & Drory 2011). The dominant galaxy type in the local universe has been found to possess pure disk properties, i.e., having a disky bulge or being bulgeless. These results reinforce the opinion that additional physical processes are needed to explain much less massive spheroidal components in disk galaxies.

Overall, strong arguments exist in favor of a process that lowers the efficiency of conversion of gas into stars (e.g., Fukugita et al. 1998). Such a process (or processes) may resolve the discrepancy between the observed galaxy luminosity function, for both high- and low-mass galaxies, and the computationally obtained DM halo mass function (HMF). Simulations of baryon evolution within the $\Lambda$ cold DM ($\Lambda$CDM) framework appear to overproduce low- and high-mass galaxies. Feedback from stellar evolution in low-mass galaxies and from AGNs in high-mass ones is expected to quench the star formation, depleting the host galaxy from its ISM, and resolving this discrepancy (e.g., Khochfar & Silk 2006; Naab et al. 2007; Somerville et al. 2008; McCarthy et al. 2012; Hilz et al. 2013); see also reviews by Veilleux et al. (2005) and Shlosman (2013).

In this work we focus on rare overdense regions (by construction) in the universe that host massive ($\sim 10^{12} h^{-1} M_\odot$, where $h$ is the Hubble constant in units of $100 \text{ km s}^{-1} \text{ Mpc}^{-1}$) DM halos by $z \sim 6$, possibly the host halos of QSOs, and quantify the evolution of DM and baryons in these regions, in comparison with the “normal” regions of average density. We compare numerical models of galactic winds by means of high-resolution cosmological simulations of galaxy formation and evolution at these high redshifts, and analyze the efficiency of these winds in carrying metals outside galaxies, into the halo and IGM environment. The effect of these winds on individual galaxies is discussed as well.

This paper is structured as follows. Section 2 describes the numerical methods, initial conditions used, and the details of the wind models. Resulting halo and galaxy properties are given in Section 3, and the effects of winds on the large-scale environment are analyzed in Section 4. Discussion and conclusions are given in the last section.

## 2. SIMULATIONS

### 2.1. The Code

We use a modified version of the tree-particle-mesh smoothed particle hydrodynamics (SPH) code GADGET-3, originally described in Springel (2005), in its conservative entropy formulation (Springel & Hernquist 2002). Our conventional code includes radiative cooling by H, He, and
metals (Choi & Nagamine 2009), a recipe for star formation and SN feedback, a phenomenological model for galactic winds, and a sub-resolution model for the multiphase ISM (Springel & Hernquist 2003, hereafter SH03). In the multiphase ISM model, star-forming SPH particles contain the cold phase that forms stars and the hot phase that results from SN heating. The cold phase contributes to the gas mass, and the hot phase contributes to gas pressure. Metal enrichment is also taken into account according to the recipe from SH03 (see also Choi & Nagamine 2009), in which the increase in metallicity of star-forming gas particles is related to the fraction of gas in the cold phase, the fraction of stars that turn into SNe, and the metal yield per SN explosion. Although metal diffusion is not explicitly implemented, metals can still be transported outside galaxies by wind particles and can enrich the halos and the IGM.

Since we focus only on $z \gtrsim 6$, before the full reionization, and do not implement on-the-fly radiative transfer of ionizing photons, the UV background is not included in our simulations. One expects this omission to have a larger effect on the overdense regions. To quantify this would require the introduction of additional free parameter(s), because the precise level of the UV background is unknown at present. We therefore refrain from doing so. For the same reason we have neglected the AGN feedback (e.g., Sijacki et al. 2009), which has been predicted to affect the more massive galaxies.

For star formation, we use the “Pressure model” (Choi & Nagamine 2010), which implements the relationship between the surface density and volume density of the gas using the prescription of Schaye & Dalla Vecchia (2008). The characteristic timescale for the SF becomes $\tau_{SF} \sim A^{-1}$ ($1 M_{\odot} \text{pc}^{-2})^n(\gamma P/G)^{(1-n)/2}$, where $A = 2.5 \pm 0.7 M_{\odot} \text{yr}^{-1}$, $n = 1.4 \pm 0.15$, $\gamma = 5/3$ (Kennicutt 1998), and P is the total effective gas thermal pressure (including the contribution from both cold and hot phases). This model reduces the high-$z$ SFR relative to the standard recipe described in SH03. Star formation is triggered when the gas density is above the threshold value $n_{crit}^{SF} = 0.6 \text{ cm}^{-3}$ (e.g., Springel 2005; Choi & Nagamine 2010; Romano-Díaz et al. 2011a). The value of $n_{crit}^{SF}$ is based on the translation of the threshold surface density in the Kennicutt-Schmidt law, SFR $\sim \Sigma_{\text{gas}}^{\alpha}$, where SFR is the disk surface density of star formation, $\Sigma_{\text{gas}}$ is the surface density of the neutral gas, and $\alpha \sim 1-2$, depending on the tracers used and on the relevant linear scales.

## 2.2. Galactic Wind Models

We have considered three different galactic outflow models: a constant-velocity wind (CW) model based on the method of SH03 and a variable-velocity wind (VW) model from Choi & Nagamine (2011, hereafter CN11), supplemented with a model without outflows (NW). All of our three wind models share the same thermal feedback as in SH03.

**No-Wind model (NW).** This model has no kinematic feedback, but maintains thermal feedback by the SNe.

**Constant Wind model (CW).** The galactic wind has been triggered by modifying some gas particles into “wind” particles. Those were not subject to hydrodynamical forces and have experienced the initial kick from the SNe. All CW particles had the same constant velocity, $v_w = 484 \text{ km s}^{-1}$, and the same mass loading factor, $\beta_w \equiv \dot{M}_w/M_{SF}$, where $\dot{M}_w$ is the mass loss in the wind and $M_{SF}$ is the SFR. The mass loading

factor has been fixed to a value $\beta_w = 2$ in agreement with the CW prescription used in CN11.

**Variable Wind model (VW).** The variable wind model from CN11 has introduced a more flexible subgrid physics compared to the simple recipe of SH03. In this paper, we adopt the 1.5ME wind model described in CN11. Briefly, the model assumes that all gas particles in a given galaxy have the same chance to become part of the wind. The probability that a gas particle becomes a wind particle is based on the values of $\beta_w$, galaxy mass, and galaxy SFR. The main parameters are the wind load $\beta_w$, defined above, and the wind velocity $v_w$—both have been constrained by observations that express these two parameters in terms of the host galaxy stellar mass, $M_*$, and the galaxy SFR.

The wind velocity is calculated as a fraction of the escape speed from the host galaxy, $v_w = \zeta v_{esc}$, where $\zeta = 1.5$ for momentum-driven winds and $\zeta = 1$ for energy-driven winds in the current setting. The empirical relation between galaxy SFR and $v_{esc}$ becomes

$$ \text{SFR} = 1.0 \left( \frac{v_{\text{esc}}}{130 \text{ km s}^{-1}} \right)^3 \left( \frac{1+z}{4} \right)^{-3/2} M_{\odot} \text{yr}^{-1}, \quad (1) $$

which is consistent with observations (e.g., Martin 2005; Weiner et al. 2009). The mass loading factor, $\beta_w$, is assumed to represent the energy-driven wind ($\beta_w \propto v_{esc}^{-1}$) in the low-density case, $n < n_{crit}^{SF}$, and the momentum-driven wind ($\beta_w \propto v_{esc}^{-2}$) for $n > n_{crit}^{SF}$. Hence, for the low-density gas, i.e., away from star-forming regions (but inside the galaxy!), we apply the parameters of the energy-driven wind for the gas particles. This procedure requires that simulations compute $\beta_w$ and $v_w$ using an on-the-fly group finder, which, in our runs, is a simplified version of the SUBFIND algorithm (Springel et al. 2001).

When an SPH particle is converted to a wind particle, it receives a $v_w$ kick and decouples from hydrodynamic forces. The direction of the kicks is chosen to be preferentially perpendicular to the angular momentum of the particle such that, overall, the winds are launched above or below the galactic plane. The wind particles are turned back to normal SPH particles when the ambient gas density becomes lower than $0.1 n_{crit}^{SF}$, or when they have traveled further than $20 \text{ kpc h}^{-1}$, whichever comes first.

## 2.3. Initial Conditions

The initial conditions have been generated using the constrained realization (CR) method (e.g., Bertschinger 1987; Hoffman & Ribak 1991; Romano-Díaz et al. 2007) and are those used by Romano-Díaz et al. (2011a), being downgraded from $2 \times 1024^3$ to $2 \times 512^3$. A CR of a Gaussian field is a random realization of such a field constructed to obey a set of linear constraints imposed on the field. The algorithm is exact, involves no iterations, and is based on the property that the residual of the field from its mean is statistically independent of the actual numerical value of the constraints (for more details see Romano-Díaz et al. 2011b). The main advantage of such a method is to bypass the sampling problem of highly overdense regions, which are rare and thus require simulation of a large volume ($\gtrsim \text{Gpc}^3$) of the universe.

The constraints were imposed on a grid of $1024^3$ within a cubic box of size $20 \text{ h}^{-1} \text{Mpc}$ to create a DM halo seed of
<table><caption>**Table 1**
Summary of the Different Simulation Runs and Their Properties</caption><thead><tr><th>Name</th><th>Initial Conditions</th><th>Wind Model</th><th>N<sub>eff</sub><sup>a</sup></th><th>R<sub>zoom</sub><sup>b</sup><br>(h<sup>-1</sup> Mpc)</th><th>R<sub>inner</sub><sup>c</sup><br>(h<sup>-1</sup> Mpc)</th><th>m<sub>DM</sub><sup>d</sup><br>(10<sup>5</sup> h<sup>-1</sup> M<sub>☉</sub>)</th><th>m<sub>gas</sub><sup>d</sup><br>(10<sup>5</sup> h<sup>-1</sup> M<sub>☉</sub>)</th><th>m<sub>*</sub><sup>d</sup><br>(10<sup>5</sup> h<sup>-1</sup> M<sub>☉</sub>)</th><th>ε<sub>grav</sub><sup>e</sup><br>(h<sup>-1</sup> kpc)</th><th>v<sub>w</sub><sup>f</sup><br>(km s<sup>-1</sup>)</th></tr></thead><tbody><tr><td>CW</td><td>CR</td><td>SH03</td><td>2 × 512<sup>3</sup></td><td>3.5</td><td>2.7</td><td>37.3</td><td>8.88</td><td>4.44</td><td>0.14</td><td>484</td></tr><tr><td>VW</td><td>CR</td><td>CN11</td><td>2 × 512<sup>3</sup></td><td>3.5</td><td>2.7</td><td>37.3</td><td>8.88</td><td>4.44</td><td>0.14</td><td>1–1.5v<sub>esc</sub></td></tr><tr><td>NW</td><td>CR</td><td>...</td><td>2 × 512<sup>3</sup></td><td>3.5</td><td>2.7</td><td>37.3</td><td>8.88</td><td>4.44</td><td>0.14</td><td>0</td></tr><tr><td>UCW</td><td>UCR</td><td>SH03</td><td>2 × 512<sup>3</sup></td><td>7.0</td><td>4.0</td><td>37.3</td><td>8.88</td><td>4.44</td><td>0.14</td><td>484</td></tr><tr><td>UVW</td><td>UCR</td><td>CN11</td><td>2 × 512<sup>3</sup></td><td>7.0</td><td>4.0</td><td>37.3</td><td>8.88</td><td>4.44</td><td>0.14</td><td>1–1.5v<sub>esc</sub></td></tr></tbody></table>

**Notes.**

a Effective resolution in the highest refinement (zoom-in) region.

b Comoving radius of the zoom-in region.

c Comoving radius of the central region used to identify halos and galaxies.

d Mass resolution in the DM, gas, and stellar component.

e Comoving gravitational softening length.

f Outflow velocity used in the galactic wind model.

~$10^{12} h^{-1} M_{\odot}$ collapsing by $z \sim 6$, according to the top-hat model. We assume the $\Lambda$CDM cosmology with WMAP5 parameters (Dunkley et al. 2009), $\Omega_{m} = 0.28$, $\Omega_{\Lambda} = 0.72$, $\Omega_{b} = 0.045$, and $h = 0.701$. The variance $\sigma_{8} = 0.817$ of the density field convolved with the top-hat window of radius $8h^{-1} Mpc^{-1}$ was used to normalize the power spectrum. The overdensity in the CR models corresponds to ~$5\sigma$ with respect to the average density of the universe.

We also evolved the same parent random realization of the density field used to construct the CR runs but without any imposed constraints to represent an average region of the universe. In order to compare the effect of winds in different environments, we have run the unconstrained simulation using the CW and VW models. In the following, we refer to individual runs by the name of the wind model used in the simulation (see Table 1): CW, VW, and NW for constrained runs (CR runs), and UCW and UVW for their unconstrained counterparts (UCR runs).

The evolution has been followed from $z = 199$ down to $z = 6$. Our simulations were run in comoving coordinates and with vacuum boundary conditions. We used the multimass approach and have downgraded the numerical resolution outside the central region using three different levels of resolution in order to reduce the computational time. The highest refinement region has a radius of $3.5 h^{-1} Mpc$ and an effective resolution of $2 \times 512^3$ in DM and SPH particles. The total mass within the computational box is $6.19 \times 10^{14} h^{-1} M_{\odot}$ and the mass of the high-resolution region is $1.45 \times 10^{13} h^{-1} M_{\odot}$. Within this region, we obtain particle masses of $3.73 \times 10^6 h^{-1} M_{\odot}$ (DM), $8.88 \times 10^5 h^{-1} M_{\odot}$ (gas), and $4.44 \times 10^5 h^{-1} M_{\odot}$ (stars). The gravitational softening is $\epsilon_{grav} = 140 h^{-1} pc$ (comoving), which is about $20 h^{-1} pc$ in physical units at $z = 6$.

In order to avoid contamination of heavy particles from the outer refinement levels, we have focused on the central inner spherical region within a radius of $2.7 h^{-1} Mpc$ for our CR simulations, and within $4 h^{-1} Mpc$ for the unconstrained (UCR) simulations. The size of the inner region in the UCR runs has been chosen by matching the cumulative mass profile with the CR models. The properties of the different simulation runs are summarized in Table 1.

Finally, the CW model is switched on at the initial redshift, while the VW model is run as NW (i.e., without outflows) until $z = 12$, at which point the winds are switched on. This is done

<table><caption>**Table 2**
Properties of the Halo and Galaxy Samples</caption><thead><tr><th>Name</th><th>M<sub>h, min</sub><sup>a</sup><br>(h<sup>-1</sup> M<sub>☉</sub>)</th><th>N<sub>h, z = 6</sub><sup>b</sup></th><th>M<sub>*</sub><sup>c</sup><br>(h<sup>-1</sup> M<sub>☉</sub>)</th><th>N<sub>g, z = 6</sub><sup>d</sup></th></tr></thead><tbody><tr><td>CW</td><td>10<sup>8</sup></td><td>6004</td><td>3 × 10<sup>7</sup></td><td>509</td></tr><tr><td>VW</td><td>10<sup>8</sup></td><td>5768</td><td>3 × 10<sup>7</sup></td><td>425</td></tr><tr><td>NW</td><td>10<sup>8</sup></td><td>5960</td><td>3 × 10<sup>7</sup></td><td>1969</td></tr><tr><td>UCW</td><td>10<sup>8</sup></td><td>9135</td><td>3 × 10<sup>7</sup></td><td>245</td></tr><tr><td>UVW</td><td>10<sup>8</sup></td><td>8470</td><td>3 × 10<sup>7</sup></td><td>92</td></tr></tbody></table>

**Notes.**

a Minimum total halo mass considered.

b Number of halos at $z = 6$.

c Minimum galaxy stellar mass considered.

d Number of galaxies at $z = 6$.

in order to speed up the simulation until the SF starts to rise and the SN feedback becomes important.

## 2.4. Group Finder

We use the group-finding algorithm HOP (Eisenstein & Hut 1998) to identify halos and galaxies. Halos are isolated according to a criterion based purely on particle density—the (local) *total* particle densities, i.e., DM + baryons, are calculated with an SPH kernel. A halo is defined as the region enclosed within a total iso-density contour of $\Delta_c \rho_{crit}(z)$, where $\Delta_c = 80$ and $\rho_{crit}(z)$ is the critical density at redshift $z$. Hence, no geometry of the particle distribution is assumed with this approach. This is a change compared to Romano-Díaz et al. (2011a), who used the standard pure DM definition for the halo, but is exactly as in Romano-Díaz et al. (2014). The benefit of including the baryonic contribution to the density field when identifying halos, while not significantly modifying the resulting halo catalog, is to automatically identify the baryonic component inside each halo. In order to avoid resolution effects, we adopt a minimum mass threshold of $M_h = 10^8 h^{-1} M_{\odot}$.

Galaxies are identified with respect to the *baryonic* density field with the outer boundary corresponding to an iso-density contour of $0.01 n_{crit}^{SF}$, which also differs from that of Romano-Díaz et al. (2011a) but is as in Romano-Díaz et al. (2014). This definition ensures the inclusion of regions that host star-forming gas ($n > n_{crit}^{SF}$), as well as lower density non-star-forming gas ($n < n_{crit}^{SF}$), which is roughly bound to the galaxy.
<figure><img src="image_2.png" /><figcaption>**Figure 1.** Left: evolution of the distribution of wind particles as a function of their outflow velocities in the VW run. The vertical dashed line corresponds to the fixed velocity of 484 km s⁻¹ used in the CW model. Right: the median wind velocities as a function of galaxy stellar mass in the VW model at z = 6. The horizontal dashed line shows the velocity used in the CW model (as above). Filled circles represent individual galaxies in the VW run, color-coded by their median SFR. The solid line shows the median relation between wind velocity and stellar mass for these galaxies, whereas the dashed and dotted-dashed lines represent the expected trends for energy-driven and momentum-driven winds, respectively.</figcaption></figure>

Following Romano-Díaz et al. (2014), we also exclude galaxies below a minimum stellar mass of $3 \times 10^7 h^{-1} M_{\odot}$, in agreement with recent observations (Ryan et al. 2014), corresponding to ~70 stellar particles. In order to verify that our simulations are not affected by numerical resolution, we have compared our CW run with the higher-resolution version presented in Romano-Díaz et al. (2014). We found that the two runs are in excellent agreement in terms of galaxy stellar mass functions, SFRs, and metallicities in galaxies, indicating that our simulation runs are numerically converged.

### 3. RESULTS: EFFECTS OF WINDS AND ENVIRONMENT ON HALO AND GALAXY POPULATION

In this section, we compare the wind models in the overdense and average-density environments and analyze their effects on the mass functions of halos and galaxies, as well as on the intrinsic properties of these objects, such as gas and stellar fractions, SFRs, and the mass-metallicity relation at $z \gtrsim 6$.

#### 3.1. Wind Velocities

Section 2.2 provided the general description of the wind models implemented in this work. We start by verifying the kinematic properties of these winds, and especially of the VW run. In Figure 1, we plot the resulting velocities of the wind particles in the VW and CW runs. The left panel exhibits the distribution of wind particles as a function of their outflow velocity, $v_w$, at different redshifts. The vertical dashed line is the value used in the CW model, i.e., 484 km s⁻¹. Even though, overall, the total number of particles in the wind is steadily increasing with time, we observe three interesting trends in the VW model: the number of wind particles corresponding to velocity at the peak of the distribution increases with redshift until $z \sim 8$ and falls off thereafter; the average wind velocity increases with time, and the tail of the distribution extends gradually toward higher velocities. Despite this latter increase, we note that, at all redshifts, the majority of wind particles in the VW run have a lower velocity than the constant velocity adopted in the CW model.

The right panel of Figure 1 displays the relation between the median $v_w$ and stellar mass in galaxies at $z=6$ in the VW run. The value of $v_w$ used in the CW model is shown as the horizontal dashed line. Individual galaxies are represented by filled circles with color scaling based on their SFR. The solid black line is the median trend calculated from these objects. The black dotted-dashed and dashed lines are the expected trends for the momentum-driven and energy-driven winds, respectively. They have been calculated from Equation (1), assuming a specific SFR (i.e., sSFR) of $2.5 \text{ Gyr}^{-1}$, in accordance with results shown in Figure 7.

We find that the scaling of wind velocities with stellar mass closely follows the energy-driven case for nearly the entire mass range considered, although the scatter increases toward low masses. This is not entirely surprising since, in the VW model, winds coming from non-star-forming regions ($n < n_{\text{crit}}^{\text{SF}}$) are assumed to be energy-driven. Given our definition of a galaxy, which includes gas well below the star formation threshold (Section 2.4), these regions usually dominate the total gas content in a given galaxy. Since all gas particles have an equal probability of turning into wind, the majority of wind particles emitted from a galaxy come from those low-density regions and thus produce energy-driven winds, as seen in Figure 1.

We also find that, except for the most massive systems, nearly all galaxies in the VW model have lower median wind velocities than the CW case, in agreement with the velocity distributions presented in the left panel. The most massive objects appear to have stronger winds in the VW case than in the CW case. However, since they constitute the deepest potential wells, these objects are also less affected by winds than the intermediate- and low-mass galaxies. Thus, overall, we expect mechanical feedback from winds in the VW run to increase with time, peak at $z \sim 7-8$, and decline afterwards, as shown by the redshift evolution in the left panel. At the same time, this feedback is less efficient than in the CW run at all times. Given the direct coupling between wind and galaxy properties in the VW model, this non-monotonic evolution is accompanied by a similar behavior in galaxy gas content and SFR as we show below. In terms of stellar feedback, the VW
model, therefore, represents an intermediate situation between
the other two cases, NW and CW.

## 3.2. Mass Distribution and Mass Functions

Figure 2 shows the mass distribution of DM, gas, and stars in all the simulation runs at $z=6$ on a scale of $1.3 h^{-1} \text{Mpc}$ (comoving), corresponding to a physical scale of $\sim 186 h^{-1} \text{ kpc}$. In the same environment (CR or UCR), no substantial differences are expected in the DM distribution on large scales since they differ only in the physics of galactic outflows implemented. On the other hand, we do observe that the distribution of baryons differs between the runs on progressively smaller spatial scales.

For a more quantitative analysis, we start by examining the differential HMFs, $\Phi = dn/d \log M_h$, i.e., the number of halos per unit volume per unit logarithmic total halo mass interval, as shown in Figure 3 (upper frame). As expected, the HMFs of the different wind models in a given environment (CR or UCR) are essentially identical over the entire mass range of $10^8 h^{-1} M_{\odot} - 10^{12} h^{-1} M_{\odot}$. Slight differences are observed between wind models because our definition of a halo is based on the *total* mass, which includes baryons (Section 2.4).

On the other hand, the effect of the overdensity is readily visible from this figure. First, the HMFs for the three CR models are shifted up with respect to the UCR HMFs because of the overdense region sampled by these models. Second, the slopes of the three models are also somewhat shallower than that of the theoretical Sheth-Tormen slope (Sheth & Tormen 1999), as expected around density peaks (e.g., Romano-Díaz et al. 2011b). Third, they extend toward higher masses than do the UCR halos. The HMFs in the UCR runs appear to be in reasonable agreement with the HMF of Sheth & Tormen (1999), except at the low-mass end, where we see some deviation from the theoretical prediction.

We have analyzed the possible causes of this excess of low-mass halos in the UCR simulations over the Sheth-Tormen HMF. First, it could follow from the spurious detection of gas-dominated blobs, which might be identified as low-mass halos by HOP. We have therefore computed the baryon fraction $f_b = M_{\text{baryons}}/M_{\text{DM}}$ for each halo in the UCR runs at $z=6$, but found that $f_b$ is always below 40%. Second, the difference is still present even when we use only the DM component to find and define our halos. Third, as a complementary test, we have also analyzed the HMFs in a *DM-only* version of our UCW model (run in a larger volume) presented in Romano-Díaz et al. (2011b) and found a good agreement with the theoretical HMFs from Jenkins et al. (2001) over the entire mass range considered here. For this reason, we are confident that the deviation between our HMFs and the theoretical prediction from Sheth & Tormen (1999) is due to the inclusion of baryonic physics, which is expected to have the most effect on low-mass halos (work in preparation).

The bottom panel of Figure 3 shows the galaxy mass functions (GMFs) for all models as a function of galaxy stellar mass, $M_*$, at different redshifts. Here we observe a clear impact of the outflow model on the shape and evolution of the GMF, reflecting the effect of winds on the growth of galaxies. At every redshift, the normalization (amplitude) of the GMF in the NW model is the highest overall among the CR runs, which means that galaxies of a given stellar mass are more abundant in the NW model than in the other two CR runs (CW and VW). However, since the CR runs have the same underlying HMFs,

an alternative way to characterize the differences between GMFs is to compare their relative shift at a given number density $\Phi$ (i.e., at a given halo mass). In this case, we find that the GMF in the NW model is shifted toward higher stellar masses, meaning that, in a given halo, galaxies have produced more stars in the NW case because of the absence of feedback from winds, unlike in the models including outflows (CW and VW).

In comparison, the effect of winds is prominent in the evolution of the VW and CW GMFs at different redshifts. At $z=11$, VW and NW GMFs closely follow each other as expected, since these two models are identical until winds are turned on at $z \sim 12$ in the VW run (Section 2.2). With decreasing redshift, the outflows become stronger, which causes both VW and CW GMFs to become shallower, especially at the low-mass end. At $z=6$, the VW GMF is more closely related to the CW one, but their slopes differ, which results in a smaller number of low-mass ($\sim 10^8 h^{-1} M_{\odot}$) galaxies and a larger population of objects at intermediate masses in the VW run. The differences mentioned above are no longer visible at the massive end ($M_* \gtrsim 10^{10} h^{-1} M_{\odot}$) of the GMFs, where all CR runs agree reasonably well. This means that winds are inefficient in removing baryons from the deepest potential wells and thus affect mostly the low-mass end of the GMFs, in agreement with previous models and observations (e.g., De Young & Heckman 1994; Mac Low & Ferrara 1999; Scannapieco et al. 2001; Benson et al. 2003; Choi & Nagamine 2011).

The evolution of the GMFs in the region of average density (UCW and UVW models) appears to resemble that of their wind-model counterparts (CW and VW) in the CR runs, but is shifted toward lower stellar masses. In particular, we observe a flattening with time of the low-mass end of the UVW GMF compared to the UCW case. However, the main difference between the two environments lies in the low-mass end of the GMFs, which appears to be steeper for UCR than for CR runs. The difference is especially pronounced when comparing CW and UCW runs, which have slopes at the low-mass end of $\sim -2$ and $\sim -1$ respectively. Because UCR and CR runs differ only in the presence of the imposed overdensity, it is clear that these trends reflect the difference in the underlying HMFs, as seen in the top panels of Figure 3. This is an important result concerning galaxy evolution in overdense environments at high redshifts, which should have implications for the expected contribution of low-mass galaxies to reionization at $z \gtrsim 6$, since measurements of the faint-end slope of the UV luminosity function at high $z$ indicate a steep trend with a slope of $\sim -2$ (Dressler et al. 2015; Song et al. 2015), in agreement with the trend we find in the regions of average density.

## 3.3. Gas Fractions: Effect of Winds

Next, we analyze the impact of the different prescriptions for galactic wind on the gas content of halos and galaxies. Figure 4 (top) displays the halo gas fraction, $f_{\text{gas,h}} = M_{\text{gas,h}}/M_h$ as a function of redshift in four mass bins, from $M_h = 10^8$ to $10^{12} h^{-1} M_{\odot}$ with a bin size of 1 dex. We keep these mass bins constant in time, which means that halos are assigned to a bin based on their total mass at a given redshift. In each panel, solid lines represent the median trend of each model, and shaded regions enclose the 20th and 80th percentiles of the population in each bin.
<figure><img src="image_3.png" /><figcaption>**Figure 2.** Snapshots of the mass distribution within the central 1.3 $h^{-1}$ Mpc (comoving) at $z = 6$ in the three different models considered here: constrained simulations (CR) with constant-velocity wind (CW, top row), variable-velocity wind (VW, second row), no wind (NW, third row), and unconstrained simulations with constant-velocity wind (UCW, fourth row) and variable-velocity wind (UVW, bottom row). Each column represents the view for a different component: dark matter (left column), gas (middle column), and stars (right column). In all panels, the same color scale represents column density in the corresponding component with an arbitrary normalization and is used only for visual representation of the mass distribution.</figcaption></figure>
<figure><img src="image_4.png" /><figcaption>**Figure 3.** Top: evolution of the halo mass functions for the three CR wind runs considered here (CW, VW, and NW) and for the UCR runs (UCW and UVW) at different redshifts. $M_h$ is the total halo mass including DM + baryons. The dashed line corresponds to the theoretical Sheth-Tormen HMF. HMFs are insensitive to the outflow model used but show an overproduction of structures in the high-density region with respect to the average. Bottom: evolution of the galaxy mass functions as a function of their stellar mass, for CR and UCR runs. The dashed lines show the trend expected from a power-law relation scaling respectively as $M^{-1}$ and $M^{-2}$. In both top and bottom panels, the mass functions are normalized by the bin size, which is set to a constant 0.3 dex in log M. Error bars correspond to the expected Poisson 1$\sigma$ deviation in each bin.</figcaption></figure>

For all models, we observe that the evolution of the gas content in halos shows a clear dependence on halo mass. Low-mass halos start with a gas fraction close to the universal baryon fraction, $f_{\text{bar,0}} = \Omega_b/\Omega_m \sim 16\%$ (dashed line), and exhibit a steeper decline with z than the more massive ones. As halo mass increases, the curve flattens, with the highest-mass bin, $\gtrsim 10^{11} h^{-1} M_{\odot}$, showing little variation in gas fraction, from z = 9 to z = 6, although trends are more noisy in this bin because of the low number of objects at each z. At later times, $f_{\text{gas,h}}$ is in the range of ~5%–15%, below the universal value $f_{\text{bar,0}}$. Note that the highest-mass bin does not contain any UCR halos since such massive objects have not yet formed in the region of average density.

The CR runs (CW, VW, and NW) show clear differences from each other, reflecting the impact of the outflow model on the gas content of halos. On average, at each z, gas fractions are highest in the NW model for lower-mass halos, $\lesssim 10^{10} h^{-1} M_{\odot}$, and become lower than in the CW and VW models for the most massive objects, $\gtrsim 10^{11} h^{-1} M_{\odot}$. Although it seems contradictory at first that NW halos have less gas, we find this trend to

be accompanied by a similar trend in SFRs. The SF is more efficient and starts earlier in the NW model because of the absence of substantial feedback from galactic outflows. The SFRs are also mass-dependent with very little SF at the low-mass end in all models. As a consequence, since halos in the NW run are better able to retain their gas than halos in the CW and VW runs, more of it gets consumed to form stars in higher-mass objects, which explains the above trend. Such a high efficiency in converting the gas into stars is one of the reasons why the NW model can be ruled out.

Figure 4 (bottom) displays the evolution of the gas fraction, $f_{\text{gas,gal}} = M_{\text{gas}}/(M_{\text{gas}} + M_{\star})$, in galaxies in four mass bins. The lines and shaded areas have the same meaning as in the case for halos (top) but now the galaxies are binned in terms of their stellar mass $M_{\star}$ at a given z.

The comparison between the CR models shows a diverse population of objects, whose gas fraction evolves in a much more complex way than in halos. The NW model shows a monotonic decline of $f_{\text{gas,gal}}$ with time, except in the highest-mass bin for $z \gtrsim 9$. By $z \sim 6$, the gas fraction is decreasing
<figure><img src="image_5.png" /><figcaption>**Figure 4.** Evolution of gas fraction in halos (top) and galaxies (bottom) as a function of redshift in different mass bins. Objects are assigned to one of the bins based on their *current* total mass (halos) or stellar mass (galaxies) at each redshift. Bin limits are indicated in each panel in units of $h^{-1}M_{\odot}$. Solid lines are the median trends in each bin. Shaded regions show the scatter around the median between the 20th and 80th percentiles of the population. In the top frame, the dashed horizontal line indicates the universal baryon fraction $\Omega_{b}/\Omega_{m}$.</figcaption></figure>

with increasing mass bin. At the same time, $f_{\text{gas,gal}}$ is independent of mass in the CW model. For the VW model, the final $f_{\text{gas,gal}}$ is also lower for higher galaxy masses. The overall evolution of CW and VW models, however, differs substantially. For almost the entire mass range, the gas fraction stays constant for $z \gtrsim 8-9$, and exhibits a mild decline thereafter. On the other hand, $f_{\text{gas,gal}}$ in the VW model shows a sharp decline for $z \gtrsim 9$, followed by a sharp increase, and saturation. In the highest-mass bin, another sharp decline can be clearly observed.

Overall, by $z \sim 6$, the NW galaxies appear to be most gas-poor, followed by VW galaxies and CW ones. Obviously, this trend correlates with the efficiency of feedback in the form of galactic outflows. The weakest feedback is that of the NW

model. The CW feedback is strong but fixed in time, while that of the VW is increasing with time, but not monotonically. So averaged over the whole galaxy population, the effect of the VW model falls roughly between that of the CW and NW models, as already mentioned in Section 3.1.

In comparison with the galaxies of Choi & Nagamine (2011), the spread between various models is smaller at $z \sim 5$, and $f_{\text{gas,gal}} \sim M_{*\text{,gal}}^{-1}$, with $f_{\text{gas,gal}} \sim 0.9$ at $\log M_{*\text{,gal}} \sim 7$. Our CR models show smaller gas fractions for all galaxies already at $z \sim 6$, which is the result of elevated SFR caused by the higher accretion rates and gas supply in the overdense region. This trend persists despite the overdense region, especially at the high-mass end. This is also true in comparison with the models of Thompson et al. (2014). The explanation for the evolution of gas fraction that we observed in the CR models lies in the interplay between the histories of SF and mass accretion onto galaxies (e.g., Romano-Díaz et al. 2014). As we discuss later on, the kinetic luminosity of these winds has a substantial effect on the temperature of the IGM and the halo gas—an effect that lowers accretion rates profoundly.

### 3.4. Gas Fractions: Effect of Environment

As is evident from Figure 4, the evolution of gas content in CW and UCW halos and galaxies is nearly identical in all mass bins, except in the mass bin $10^{10}-10^{11} h^{-1} M_{\odot}$, where the most massive objects are found in the UCR runs. Since both models use the same outflow prescription, this strongly suggests that, for a given wind model, the overdensity has little effect on the amount of gas inside these halos and galaxies.

When comparing VW and UVW halos, we find that their gas fractions tend to follow similar trends in each mass bin, but with larger relative differences than between CW and UCW runs. This indicates that the VW model couples galaxy evolution with the large-scale environment more strongly than the CW case.

However, addressing the effect of the environment on gas fraction is difficult from a direct comparison between CR and UCR evolution from Figure 4 because objects in the same mass bin can reside in different environments. For example, low-mass galaxies can be found in underdense regions in CR and UCR models, but also as satellites around the massive halos in the overdense region. This effectively causes mixing among objects of similar masses and can thus wash out the differences induced by the presence of the overdensity.

To investigate further the effects of the environment on the gas content in galaxies, it is thus easier to directly divide objects according to their large-scale environment. Figure 5 (top) shows the redshift evolution of gas fractions in galaxies in different density bins for the CR and UCR models that include outflows. The bins are defined based on the smoothed relative density contrast $\delta = \rho/\bar{\rho}_{\text{b}} - 1$ where $\bar{\rho}_{\text{b}} = \Omega_{b}\rho_{\text{crit}}$ is the universal mean baryonic density. The smoothed galaxy density field $\rho$ has been obtained by convolving the local galaxy densities calculated using an SPH method with a top-hat filter of scale 250 $h^{-1}$ kpc. Note that we use the universal mean baryonic density $\bar{\rho}_{\text{b}}$ in our definition of $\delta$ instead of the mean baryonic density inside the simulation box in order to have a common normalization factor for both CR and UCR models. Galaxies are divided between underdense ($\delta < 0$, top left panel) and overdense ($\delta > 0$, top middle panel) bins and we also show the evolution for the entire galaxy population for comparison (top right panel). Based on our definition of galaxy
<figure><img src="image_6.png" /><figcaption>**Figure 5.** Top: evolution of gas fraction in galaxies in the CR and UCR runs as a function of redshift, similar to Figure 4, but objects are now binned by density contrast $\delta = \rho/\bar{\rho}_b - 1$ where $\bar{\rho}_b$ is the universal mean baryonic density. The bins are defined as: $\delta < 0$ for underdense (left panel), $\delta > 0$ for overdense (middle panel), and all galaxies (right panel). Bottom: mass-weighted average gas fractions in galaxies as a function of stellar mass at $z = 10$ (left panel), $z = 8$ (middle panel), and $z = 6$ (right panel).</figcaption></figure>

density, we found no galaxies with $\delta > 0$ in the UCR models. The comparison between CR and UCR runs is therefore possible only for galaxies residing in underdense regions.

The top left panel of Figure 5 shows that the evolution of gas fraction in these galaxies residing in underdense regions follows a similar trend in both CW and UCW models with a modest decline from ~65%–70% at $z = 12$ to ~55% at $z = 6$. Differences between the two models are not significant given the noise introduced by the small number of galaxies formed in the UCW model. In comparison, the gas fraction in galaxies residing in overdense regions in the CW model exhibits a much steeper decline from $z = 10$ to $z = 6$ (top middle panel) due to the accelerated evolution in these regions. The UVW galaxies show a more complicated behavior, but decline as well after $z \sim 8$, in tandem with CW. The evolution for the entire galaxy population (top right panel) shows larger differences between CW and UCW due to the mixing of both underdense and overdense environments in the CW case but still shows a similar decreasing trend with redshift and comparable gas fraction levels at all $z$. The VW and UVW models in the underdense bin exhibit much less correlation than CW and UCW among themselves, as already noted in the previous section.

The similar evolution of gas fraction for CW and UCW galaxies residing in the same environment further confirms that the gas content is mainly affected by the different wind models and not by the presence of the overdensity. The gas fractions in VW and UVW models differ more among themselves, which can be explained by the weaker winds they produce. To more fully understand the reason for this effect, we need to

investigate how the galaxy gas fractions depend on the outflow prescription at a given redshift. Figure 5 (bottom) shows the mass-weighted average galaxy gas fractions as a function of galaxy stellar mass in all the simulation runs at $z = 10, 8,$ and 6. Since we are only interested in the average trends, the curves have been smoothed with a boxcar filter of size 0.5 dex in stellar mass to reduce Poisson noise caused by the small number of objects in each mass bin.

We find significant differences among the various wind models at all $z$, but the main discrepancy can be observed between CW and the other two runs, VW and NW. In particular, there is a clear trend of decreasing gas fraction with increasing stellar mass in both VW and NW, which is not seen in CW. The latter shows nearly constant gas fractions across the entire mass range $M_* = 10^{7.5}-10^{11} h^{-1} M_\odot$. The UCW model also exhibits nearly mass-independent gas fractions at all $z$, which follows the CW trend as expected. This is in agreement with results shown in the top panels of Figures 5 and 4, although there are clear signs of divergence between the two models for $z \le 8$. These results thus suggest that the likely explanation for the apparent similarities in galaxy gas content between CW and UCW can be attributed, for the most part, to the prescription of constant and strong outflow used in both runs. Since winds scale independently of galaxy properties (mass and SFR) in this outflow model, the gas content inside galaxies varies weakly with stellar mass. As a consequence, the CW and UCW runs produce galaxies with nearly constant gas fractions, which explains the similar behavior observed between the two models despite their different environments.

### 3.5. Star Formation Histories (SFHs)

To construct the SFH for a given galaxy selected at $z = 6$, we identify the progenitors at earlier redshifts by matching the unique identification numbers (IDs) associated with each particle. The ancestor galaxy is taken to be the progenitor containing the largest number of particles with IDs matching the ones in the descendant galaxy. We repeat this process iteratively to obtain the list of ancestors of the selected galaxy, from $z = 6$ to $z = 12$. Figure 6 shows the redshift evolution of galaxy SFRs in different mass bins. Galaxies are assigned to one of the bins based on their stellar mass $M_*^{z=6}$ at $z = 6$, and their SFH is then followed back to $z = 12$ using the above method. We stress that this mass binning is different from the one in which galaxies are selected by mass at each redshift, as used in Sections 3.3 and 3.6. For this reason, correlating the evolution of SFR and gas fraction cannot be attempted using a straightforward comparison between Figures 4 and 6, because they represent different galaxy populations at each redshift.

We find a strong dependence of SFHs on stellar mass in all models. Overall, massive galaxies exhibit a strongly rising SFR, while the lowest-mass galaxies maintain flat SFRs on average. These trends are well correlated with the rates of mass accretion onto the galaxies, analyzed in Romano-Díaz et al. (2014). More specifically, in the lowest mass range, $7.5 < \log M_*^{z=6} < 8.5$, SFRs are essentially constant with $z$ for all models, and lie in the range of ~0.1–1 $M_\odot$ yr$^{-1}$. The intermediate-mass bins, $8.5 < \log M_*^{z=6} < 9.5$ and $9.5 < \log M_*^{z=6} < 10.5$, exhibit an increase in SFR by a factor of ~7–8 between $z = 12$ and $z = 6$. Finally, the most massive galaxies display the strongest increase, roughly a factor of 30 over this redshift interval. Note, however, that small SF levels seen for low-mass galaxies likely result from poorly resolved
<figure><img src="image_7.png" /><figcaption>**Figure 6.** Star formation histories for galaxies in different stellar mass bins. Galaxies are selected based on their stellar mass $M_{\star}^{z=6}$ at $z=6$ (which is different than selecting galaxies based on their mass at every redshift) and their ancestors are identified at earlier redshifts to be the most massive object among all the progenitors. The bin limits are indicated in each panel in units of $h^{-1}M_{\odot}$. Solid lines show the median SFR evolution of all the ancestors at each $z$. The shaded regions are bounded by the 20th and 80th percentiles.</figcaption></figure>

outflows, with only a few particles for these objects. For this reason, we focus mostly on SFRs in intermediate- and high-mass galaxies in the following.

At high $z$, the SFRs for the most massive galaxies appear similar among all the CR runs (note that these galaxies are absent from the UCR runs). The rates start to diverge below $z \sim 9$, with the CW model becoming the dominant one, followed by the VW and the NW, each lower by a factor of 2–3 by $z=6$. For these galaxies, the SFR has reached $\sim 10^3 M_{\odot} \text{ yr}^{-1}$, bringing them into the observed regime of LBGs at high $z$ (see also Yajima et al. 2015). We also observe a similar trend in the intermediate-mass bins, where SFR in the NW model starts to decrease after $z \sim 8$, while models including galactic outflows continue to increase their SFRs. By the end of the simulations, the CW model has the highest median SFR with $\sim 30–40 M_{\odot} \text{ yr}^{-1}$ for $9.5 < \log M_{\star}^{z=6} < 10.5$, and $\sim 600–800 M_{\odot} \text{ yr}^{-1}$ for $10.5 < \log M_{\star}^{z=6} < 11.5$. The corresponding median SFRs in the NW run are $\sim 10 M_{\odot} \text{ yr}^{-1}$ and $\sim 100 M_{\odot} \text{ yr}^{-1}$, and the VW model ends up with intermediate values between the CW and NW runs.

These results indicate that, even at high $z$, galactic outflows do play an important role in controlling galaxy growth in overdense regions. In the NW run, the absence of these outflows means the absence of a strong feedback that can delay or quench star formation. As a consequence, the peak of star formation in the NW run happens at higher redshifts, $z > 6$, as noted above, and thus the available gas supply for SF in galaxies is already consumed by $z=6$. Subsequent gas accretion could in principle increase SFRs again at later times

in the NW run but, at $z=6$, we see that accretion has not yet been able to replenish the gas content in galaxies to maintain high SFRs. For this reason, SFRs are higher in CW and VW models at $z=6$ since the feedback from winds in those models acts to effectively delay the peak of star formation activity.

This is in agreement with Choi & Nagamine (2011) and Anglés-Alcázar et al. (2014), which have analyzed galaxy evolution in regions of average density, but lack massive galaxies at these redshifts. However, since galaxy evolution is accelerated in the overdense regions, the peak of SF activity happens at earlier times than in an environment of average density. When comparing CW and VW models, results published in the literature show that the CW SFRs are generally the lowest ones, as measured in simulations at low $z \sim 3–4$ (e.g., Choi & Nagamine 2011). However, the cosmic SFR evolution presented in Choi & Nagamine (2011, see their Figure 5) shows that the difference between CW and VW SFRs decreases with increasing redshift at $z > 6$, and that these two models eventually reach similar SF levels at high $z$.

### 3.6. Specific Star Formation Rates

For a more detailed picture of the SF efficiency in different runs, we investigate the specific star formation rate, which represents SFR per unit mass, and essentially indicates the production rate of stellar mass per stellar mass. Figure 7 displays the relation between the median sSFR and stellar mass in galaxies at $z=6$ (left panel) and the evolution of sSFR with redshift (right panel) in all runs. Observational data from Labbé et al. (2013) ($z=8$), Duncan et al. (2014) ($z=6$ and 7), Salmon
<figure><img src="image_8.png" /><figcaption>**Figure 7.** Median specific star formation rate as a function of stellar mass (left) and redshift (right). Observational data points are taken from Labbé et al. (2013) at $z=8$, Salmon et al. (2015) at $z=6$, and Duncan et al. (2014) at $z=6$ and 7. The yellow hexagon (left frame) is an observational constraint at $z \sim 5-5.5$ from the VUDS survey (Tasca et al. 2015, see their Figure 3). The relation between sSFR and stellar mass in the left panel is shown for galaxies at $z=6$. The dashed gray line is derived from the fitted SFR-$M_*$ relation at $z=6$ from Salmon et al. (2015). In the right panel, the solid lines correspond to the median sSFR of all the galaxies at each redshift. The dashed lines are the median sSFR for galaxies with stellar masses in the range $10^9-10^{9.5} h^{-1} M_{\odot}$, which approximately matches the mass range of observed galaxies corresponding to the data points.</figcaption></figure>

et al. (2015) ($z=6$), and Tasca et al. (2015) ($z=5.5$) have been added to both panels. In the left panel, galaxies are binned by stellar mass with intervals of 0.3 dex.

We observe a nearly constant sSFR with $M_*$ in the NW and VW models, with sSFR ~1.5 and 2.5 Gyr⁻¹, respectively. At the same time, the CW model shows an overall positive correlation between sSFR and $M_*$. This upward turn above the characteristic mass of ~$10^9 h^{-1} M_{\odot}$ is in contradiction with the observational trend shown in Figure 7.

We also find that the sSFR in the CW run falls within the observed range, but the existence of the peak in the mass range of $10^9-10^{10} h^{-1} M_{\odot}$ does not agree with the observations of Salmon et al. (2015). On the other hand, for the NW run, the overall trend roughly agrees with the observed one, but the sSFR lies outside the observed range. The best fit to observations is that of the VW model—both the trend with $M_*$ and the sSFR range agree, although the median observed fit is slightly shifted away from the model values.

The right-hand panel of Figure 7 shows the evolution of the sSFR with redshift. Solid lines correspond to the median trends for the whole galaxy population. The dashed lines show the evolution only for galaxies with stellar masses in the range $10^9-10^{9.5} h^{-1} M_{\odot}$, corresponding approximately to the estimated masses of objects targeted by observers. We find the median evolution of all galaxies to be dominated by the low- and intermediate-mass objects because of a higher number of these galaxies compared to the massive ones. The median sSFRs of all galaxies in the CW and UCW runs are in good agreement—again underlining that gas content and SF in low-mass galaxies are not affected significantly by the overdensity.

The sSFR in the NW run is the highest for $z \gtrsim 7$, where it falls below the VW model. In all models except VW, the sSFR drops by almost an order of magnitude, from ~10 Gyr⁻¹ at $z=12$ to ~1 Gyr⁻¹ at $z=6$. The evolution of sSFR in the VW case exhibits a trend similar to the one found for the gas fraction in low-mass galaxies (Figure 4). As we have argued in the previous section, This trend is caused by scaling of the feedback strength with the SFR in the VW run. At high z, outflows quench the SF, which in turn decreases the feedback

from the winds and allows an elevated SFR at later times. The evolution of the median sSFR for the whole galaxy samples does not reproduce correctly the observations for $z < 8$. The explanation of the mismatch lies in the dominant contribution from low-mass galaxies ($M_* \lesssim 10^9 h^{-1} M_{\odot}$) to the global sSFR. When comparing the evolution for galaxies in the observed stellar mass bin (dashed lines), both CW and VW runs, as well as the UCR run, appear to be in good agreement with the observed increase in sSFR with redshift.

### 3.7. Mass-Metallicity Relation

Next, we focus on the mass-metallicity relation in halos and galaxies. Figure 8 (left panels) shows the metallicity of gas and stars in halos as a function of the total halo mass at $z=6$. The highest metallicities are found in the most massive halos because their star formation starts earlier and they exhibit a higher SFR than lower-mass objects. We observe a rapid increase in the metallicity in the low-mass halos, and it levels off at the intermediate/high-mass end. The CW and VW runs reach close to solar metallicity for $M_h \gtrsim 10^{11} h^{-1} M_{\odot}$ by the end of the simulations. For smaller $M_h$, the CW and VW models display a slow decline in Z, which steepens dramatically below $M_h \sim 10^{10} h^{-1} M_{\odot}$ for the VW, and below ~$10^9 h^{-1} M_{\odot}$ for the CW. This “knee” is moving gradually toward higher masses with redshift.

While halos with $M_h \sim 10^{11}-10^{12} h^{-1} M_{\odot}$ have a median $Z/Z_{\odot} \sim 0.6$ (gas) and ~0.8 (stars), the lower-mass halos, with ~$10^8-10^9 h^{-1} M_{\odot}$, show $Z/Z_{\odot} \sim 0.01-0.03$ (gas) and ~0.03 (stars), respectively. Note that the NW model achieves the highest metallicity in all $M_h$, for both gas and stars. Furthermore, compared with the metallicities of halo gas, stellar metallicities display higher median Z—this means that the gas metallicity is diluted by the cold accretion from the IGM, in agreement with Romano-Díaz et al. (2014).

Galaxy metallicities in the gaseous and stellar components are shown in Figure 8 (right panels). As expected, the metallicity of the NW model exceeds substantially those of CW and VW in both components, due to the metals being locked in galaxies and their immediate vicinity. A possible
<figure><img src="image_9.png" /><figcaption>**Figure 8.** Left: gas (top) and stellar (bottom) metallicities of halos (including gas in galaxies) as functions of total halo mass at $z = 6$ in all simulation runs. Right: gas (top) and stellar (bottom) metallicities of galaxies as functions of stellar mass at the same redshift. The solid lines indicate the median values, while the colored regions correspond to the values between the 20th and 80th percentiles. Metallicity profiles of the IGM gas are shown in Figure 11.</figcaption></figure>

<figure><img src="image_10.png" /><figcaption>**Figure 9.** Maps of temperature (left) and metallicity (right) of the IGM gas (excluding the gas inside halos) within the central $2.7 h^{-3} Mpc^{-3}$ inner region, at $z = 10, 8, \text{ and } 6$ (top to bottom) for the three CR models (CW, VW, and NW).</figcaption></figure>

exception may be the highest-mass bin for the ISM, where CW model shows marginally higher metallicity. The bimodal behavior observed for the halos can be seen here as well. Stars also remain more metal-rich than the gas, again in tandem

with the halos. Future observations will constrain these models. We point out that, in Figure 8, we have computed only the average metallicity inside individual galaxies and ignored the scatter in metallicity in a given galaxy.
## 4. RESULTS: EFFECT OF WINDS ON ENVIRONMENT

In this section we focus on the effects of winds on the thermal evolution and metal enrichment in the circumgalactic and intergalactic medium. The effects of winds on the ISM have been discussed in Sections 3.3 and 3.4.

### 4.1. Temperature and Metallicity Maps

The temperature maps of the IGM gas, i.e., gas outside galaxies and their host halos, in the inner central region of the simulation box clearly display the growing difference between the wind models with redshift. Figure 9 (left panels) underlines profound differences existing between the CR runs already at $z \sim 10$. The NW run lacks gas with $T \gtrsim$ few $\times 10^4$ K, except for a few spots near the central peak in density corresponding to the growing massive halo. The warmest gas can be found in the large-scale filaments, and it roughly delineates them in all CR models. The gas at a temperature of a few $\times 10^4$ K is more closely associated with large-scale structures in the VW case, because the feedback from outflows is lower than in the CW run, where this gas “spills out” of the filaments. The rest of the volume hosts gas with $T \lesssim 10^4$ K.

The redshift evolution of these maps shows gradually heating and cooling regions, with gas at a temperature of a few $\times 10^5$ K at $z = 8$, and a few $\times 10^6$ K at $z = 6$. By the end of the run, the massive central halo is surrounded by a hot bubble of gas at $\sim 10^7$ K in CW and VW wind models. Away from the most massive halo, the filament temperature drops, most profoundly in the NW model. By $z \sim 6$, the CW gas in nearly the whole volume is heated up to $10^5-10^7$ K. In this model, the gas at $10^5-10^6$ K is present in the main filament, and gas at $\sim 10^5$ K in other filaments. In VW and NW models, the gas outside the filaments remains cold.

Hence, at higher redshifts, $z \sim 10$, we observe increasing temperatures of the IGM gas, along the sequence from NW to VW and to CW. With decreasing $z$, the central bubble heats up, and so do the filaments. The outflows in the CW model are capable of heating up the IGM, while the VW outflows have a dramatically lower impact on the temperature of the IGM at these redshifts, but nevertheless heat it up toward $z \sim 6$.

Figure 9 (right panels) shows maps of the metallicity distribution in the IGM gas at similar redshifts for all the CR models. At $z = 6$, the IGM peak metallicity has reached $\sim 0.1 Z_\odot$ inside the central hot bubble and along the main filament corresponding to the regions surrounding the massive central structure in the CW and VW models. The peak metallicity in the NW is a factor of a few lower even though the SFRs are higher in this model at high $z$. The extended region in the CW model shows a negative metallicity gradient away from the central bubble and the main filament. Still, most of the volume is polluted by metals and has $Z \gtrsim 10^{-2} Z_\odot$. The VW exhibits a much lower metallicity outside the filaments, $Z \lesssim 10^{-3} Z_\odot$, whereas the NW IGM remains mostly pristine. The reason for this is the absence of outflows, which causes metals to be locked in galaxies and their halos. In comparison, solar metallicities have been reached only within galaxies of the NW model (see Figure 8).

Figure 10 displays the temperature and metallicity maps of the IGM in the UCR runs at $z = 6$, and can be compared directly with the CW and VW maps in Figure 9, as both figures use the same color scales. Since the UCR runs represent a region of average density and lack massive halos (e.g.,

<figure><img src="image_11.png" /><figcaption>Figure 10. Maps of temperature (top) and metallicity (bottom) of the IGM in the UCR runs at $z = 6$ within the central $4.0 h^{-3}$ Mpc$^3$ inner region. The color scales are the same as in Figure 9.</figcaption></figure>

Figure 3), both the peak and average IGM temperatures are dramatically lower than in the CR runs. The peak metallicity in the IGM and its average are also lower than in the CR runs. Thus the impact of the overdensity is readily visible from these two figures—the universe of average density remains relatively cool at $z \sim 6$.

### 4.2. Temperature and Metallicity Profiles

A more quantitative comparison between the different runs on the state of the IGM at $z = 6$ is displayed in Figure 11, which shows the spherically averaged profiles of IGM temperature (top) and metallicity (bottom) from the center of the computational box to the radius $R_{\text{inner}}$ of the inner high-resolution region, where $R_{\text{inner}} = 2.7 h^{-1}$ Mpc and $R_{\text{inner}} = 4 h^{-1}$ Mpc for the CR and UCR models respectively. The profiles are computed by averaging the temperature and metallicity of the IGM gas on spherical shells of thickness $50 h^{-1}$ kpc.

Both the temperature and metallicity profiles agree with the qualitative results from Figures 9 and 10. The CR runs have similar IGM temperature within the central $1 h^{-1}$ Mpc, corresponding approximately to the size of the imposed constraint, which is visible as a hot bubble on Figure 9. In this central region, the IGM temperature reaches $T \sim 10^7$ K, which is consistent with the virial temperature for a halo of mass $\sim 10^{12} M_\odot$, corresponding to the mass of the central halo that forms in the overdense region. On the other hand, the IGM temperature in the UCR runs remains constant with radius at much lower values, around $\sim 10^{3.5}-10^4$ K and $\sim 10^{4.5}-10^5$ K for the UVW and UCW models respectively. Thus, it appears that shock heating caused by the presence of the massive structure in the overdense region is largely responsible for the elevated IGM temperature within $1 h^{-1}$ Mpc in the CR runs. This also
<figure><img src="image_12.png" /><figcaption>Figure 11. Radial profiles of temperature (top panel) and metallicity (bottom panel) of the IGM gas at z = 6 calculated from the center of the computational domain to the radius R<sub>inner</sub> of the inner high-resolution region (see Table 1).</figcaption></figure>

explains why the CR runs show similar IGM temperature in the central region, because the differences in temperature due to the different wind models (which are much smaller) get washed out.

The metallicity profiles (Figure 11, bottom) are in general much flatter than the temperature ones. The lowest metallicities are seen in the NW run, followed by the VW and CW models, in agreement with Figure 9. The UCW IGM metallicity is smaller than both CR models including winds (CW and VW) within ~2 h⁻¹ Mpc and then rises on larger scales to reach levels comparable with the CW run. This increase is associated with the presence of a massive filament in the UCR runs (Figure 10), because a similar rise can be seen in the UCW temperature profile at the same distance from the center.

### 4.3. Gas Phases at z = 6

A supplementary global view of the temperature and metallicity of the gas can be taken from Figure 12, where we

show the gas phase diagrams in the temperature–density plane at z = 6 for all gas particles found in the high-resolution region. Gas particles are also colored based on their metallicity. The densities are now normalized by the average baryon density in the high-resolution region, as opposed to using the universal value like we did in Section 3.4. This normalization is thus different between CR and UCR runs since the average density is higher in the overdense region. For this reason, the vertical dashed line, which separates the non-star-forming gas (to its left) from star-forming gas (to its right), is slightly offset in the bottom panels (UCR) from its position in the top panels (CR). We observe the highest metallicities, around solar values, in the star-forming gas that is confined deep within galaxies. One can easily identify the wind gas in the lower right part of the diagrams, because it is present in the VW, CW, UCW, and UVW cases, and is metal-rich, but not in the NW case.

The fraction of low-density gas polluted with metals (green and red) compared to pristine (blue), and the highest metallicities attained there, depend strongly on the wind model. The metal-rich, low-density gas originates in the galactic outflows, and heats up as a result of the shock-stopping of these high-speed outflows. This gas lies at temperatures of ~10⁵–10⁷ K and densities log ρ<sub>b</sub>/ρ̃<sub>b</sub> ≲ 4. Here the baryon density ρ<sub>b</sub> is normalized by the average density in the high-resolution box. We note that a large amount of the ~10⁴ K star-forming gas found in our models at high densities is due to the use of the delayed star formation algorithm, the Pressure model (Section 2.1), and because of the absence of H₂ cooling in our simulations. The amount of heated gas increases substantially from the NW to VW and CW runs. The hot and low-density pristine gas at T ≳ 3 × 10⁶ K in the CR runs consists mostly of virially shock-heated material infalling onto the most massive halo, and, as such, is not observed in the UCR models. The amount of intermediate-metallicity (Z/Z<sub>☉</sub> ~ 10⁻¹) gas at T ≳ 10⁶ K is steadily increasing in all CR models.

Figure 13 provides the metallicity distribution for different wind models in various gas phases—the IGM, halo gas, and the ISM, as well as the total distribution for gas particles in the high-resolution region. The ISM includes all the gas inside galaxies, the halo gas includes the gas within halos but outside galaxies, and the IGM includes the gas outside the halos. In the IGM panel, we observe that the CW models for CR and UCR exhibit the highest metallicities, basically indistinguishable from each other, with maxima located at log(Z/Z<sub>☉</sub>) ~ (-0.8) – (-0.7). The least efficient injection of metals happens in the NW model (the maximum at log(Z/Z<sub>☉</sub>) ~ -2.2), while VW occupies broadly the region between the CW and NW models. It displays a nearly flat top at log(Z/Z<sub>☉</sub>) ~ (-2) – (-1). At the other extreme, the metal distributions in the ISM are basically identical, and the maxima are found between log(Z/Z<sub>☉</sub>) ~ (-1) and 0. The metal distributions in the halo gas are flat-topped for NW and VW models, and show maxima at log(Z/Z<sub>☉</sub>) ~ -1 for the UCW and ~0.7 for CW. From Figure 13, we conclude that the ISM is polluted by metals in all models. Moreover, there is a large difference in expelling the metals to the halo and the IGM. The CW and UCW models appear to be the most efficient in this process in a given environment.
<figure><img src="image_13.png" /><figcaption>**Figure 12.** Distributions of all gas particles found in the high-resolution region in the temperature–density plane at $z = 6$ in all the simulation runs: (a) CW, (b) VW, (c) NW, (d) UCW, and (e) UVW. Gas particles are colored based on their metallicity as indicated by the color scale on the right. The density is normalized by the average baryon density in the high-resolution computational box, which is different between CR and UCR runs. For this reason, the vertical dashed line corresponding to the threshold density for star formation $n_{\text{crit}}^{\text{SF}}$ is offset in the bottom panels (UCR runs) from its position in the top panels (CR runs).</figcaption></figure>

<figure><img src="image_14.png" /><figcaption>**Figure 13.** Metallicity distributions for different gas phases at $z = 6$ in the high-resolution region: all gas particles (top left), IGM (top right), halo gas (bottom left) and ISM (bottom right) for all the simulation runs: CW (red), VW (green), NW (blue), UCW (black), and UVW (gray).</figcaption></figure>

4.4. Evolution of Gas Phases

To take an alternative view of the evolution of gas, we subdivide the gas phases in a different way. In this section, we alter definitions for the IGM and other gas phases used in the previous section in order to facilitate the comparison with previous works. Figure 14 compares the evolution of various gas phases, namely, diffuse, warm-hot + hot, IGM, and condensed gas. These phases are defined based on the temperature of the gas and the density fluctuation $\delta = \rho_b/\bar{\rho}_b - 1$. For this purpose, we employ the definitions from Davé et al. (1999) and Choi & Nagamine (2011). The “diffuse” gas is defined by $T < 10^5$ K and $\delta < 1000$, the “hot” gas by $T > 10^7$ K, and the “warm-hot” gas by $10^5 < T < 10^7$ K. The “condensed” phase consists of non-star-forming gas with $T < 10^5$ K and $\delta > 1000$, and of a star-forming gas as well as stars.

We observe that the mass fraction of the “condensed” phase is increasing with time for all models, especially for the NW, while the IGM (non-star-forming) fraction decreases. The fraction of “diffuse” filamentary gas is decreasing faster for the CW model than for the VW one. At the same time, the fractions of both the UCW and UVW diffuse gas remain very high at all $z$, while the condensed fraction stays very low. While in the warm-hot + hot phase the NW and VW evolution are similar, in the condensed and IGM phases it is the VW and CW runs that evolve in tandem.
<figure><img src="image_15.png" /><figcaption>**Figure 14.** Evolution of the mass fractions of the four different baryon phases (hot, warm-hot, diffuse, and condensed)—see Section 4.4 for definitions) for the CR and UCR simulations (see the text for definitions). Panel (b) shows the addition of the two phases, “warm-hot + hot,” and panel (c) shows the total mass in the IGM, i.e., “diffuse + warm-hot + hot.” Most of the UCR baryons are in the “diffuse” phase, while in the CW the majority are found in the “warm-hot” phase. The “condensed” phase is also more populated in the CR runs, especially for the NW model.</figcaption></figure>

While the mass fractions in our UCRs and the models of Choi & Nagamine (2011) appear to be reasonably similar, the CR models differ substantially, both from the UCR cases and among themselves. This difference between the CR and the UCR models is expected, because of the overdensity present in the CR runs. But we also observe that the NW model exhibits a much higher fraction of cold gas, as it evolves toward $z \sim 6$, and consequently it exhibits smaller fraction of IGM gas. The mass fraction of the condensed phase, which also includes stars, in the NW run reaches ~50%, while for the VW and CW models it is only ~32%. The diffuse component shows that the fraction drops to ~43% in NW, to ~51% in CW, and to only ~57% in VW.

## 5. DISCUSSION AND CONCLUSIONS

We have used high-resolution cosmological zoom-in simulations of average and overdense regions in the universe

in order to follow the evolution of their halo and galaxy populations at $z \gtrsim 6$. Specifically, we compared this evolution by applying different galactic wind models in order to analyze the effects of outflows on:

* galaxy evolution, including their environment;

* the thermodynamic state of the IGM;

* the spread of metals in the ISM, DM halos, and the IGM.

We have tested three different prescriptions for galactic winds: a constant-velocity wind (CW) model from Springel & Hernquist (2003), a variable-velocity wind (VW) model based on the method of Choi & Nagamine (2011), and a model without outflows (NW). Each model has been followed down to $z = 6$. In order to get an insight into the environmental dependence of galaxy evolution, the CW prescription for galactic outflows has been applied to both constrained (CR) and unconstrained (UCR) simulations, representing overdense and average regions in the universe.
The most general conclusions that emerge from this work are that, regardless of the applied wind prescription, (1) galaxies in the overdense region evolve much faster than their counterparts in the average regions (e.g., Figures 2 and 3, see also Romano-Díaz et al. 2014), and (2) the low-mass end of the GMF is shallower in the overdense regions—an important issue to be addressed elsewhere. This is a direct corollary of the imposed constraint that has been designed to collapse by $z \sim 6$, based on the top-hat model. As a consequence, we observe an accelerated evolution of its environment within the correlation length of the constrained halo (e.g., van de Weygaert & Bertschinger 1996). DM halos in the vicinity of the constraint experience a faster evolution than their average-density counterparts, inducing earlier baryon collapse and setting the stage for star formation. We find that the galaxy population in such regions is much more evolved by the end of simulation. This is indicated not only by the presence of more massive galaxies (up to $1-2 \times 10^{11} h^{-1} M_{\odot}$ in stellar mass) that are absent from the region of average density (where $M_* \lesssim 10^{9.5} h^{-1} M_{\odot}$), but also by a substantially larger number of galaxies over the mass range considered here.

Finally, we observe that galactic winds have a substantial effect on the galactic environment, which includes modifying its thermodynamic state and metallicity. While CWs show little dependence on the overdensity, the corollaries of VWs are more diverse. This is clearly related to the amount of energy deposited by the VWs in DM halos and the IGM, which varies with galaxy properties, and in general it is smaller than that deposited by the CWs. In other words, the VWs appear less “destructive” than their constant-velocity counterparts, and more fine-tuned to the mode of galaxy growth in the form of cold accretion flows—their overall effect differs between overdense and average regions. On the other hand, the CWs weaken the flows of cosmological filaments around galaxies, irrespective of surrounding densities.

Galactic winds, as a form of feedback, serve as main agents in distributing recycled, metal-enriched gas, which is used to form new stars. The effect of winds on the evolution of our modeled galaxies can be most clearly observed in the GMF (Figure 3, bottom panel), gas fractions (Figure 4), and SFRs (Figure 6). The IGM properties can be affected by galactic winds, so comparing the observed IGM properties would provide the clearest justification for the wind models. However, the IGM at $z \gtrsim 6$ is not fully reionized, and the conventional technique of IGM observation, i.e., measurement of quasar absorption lines, cannot be used in this epoch. This limits the resources for observing the IGM that can constrain our wind models. However, the pre-reionized H I status at $z \gtrsim 6$ provides an interesting implication that can discriminate between our CW, VW, and NW models, as discussed below.

The model with no wind feedback produces the largest number of galaxies because the gas is hardly expelled from their DM halos. In this case, the feedback is limited to the thermal feedback from the SNe, and therefore contributes to the survival of gas and galaxies within the DM substructure. Most of this gas remains locked up within their DM halos, cools down, and is accreted by the galaxies (e.g., Kereš et al. 2005; Dekel & Birnboim 2006; Romano-Díaz et al. 2014). As a consequence, these galaxies become very metal-rich early-on (Figure 8), enhance their SFRs (Figure 4), consume their gas very quickly (Figure 4), and become gas-poor and very compact already at $z \sim 6$. An overall steady decay in the gas

<figure><img src="image_16.png" /><figcaption>**Figure 15.** Evolution of star formation density with respect to the total computational box for all galaxies that have $M_* \gtrsim 10^9 h^{-1} M_{\odot}$ at each redshift. This is different from Figure 6, which provides star formation histories for galaxies chosen at $z=6$. The dashed line is the best-fitting line from Madau & Dickinson (2014), and is based on UV and IR data. The observational points are from Bouwens et al. (2015)—blue points are uncorrected for dust, while red points have been corrected for the dust extinction.</figcaption></figure>

fraction within the NW galaxies is correlated with their SFRs (Figure 6). The smallest effect is observed in the lowest-mass bin, i.e., $7.5 < \log M_*^{z=6} < 8.5$. This happens because such objects are less prone to the feedback effects because of their already very low SFRs.

Our conclusions regarding the NW model pertain mostly to the decline in the gas fraction in intermediate-mass and massive galaxies, and consequently the decline in their respective SFRs. This implies reddening of such galaxies due to aging of their stellar populations at relatively high $z$ (e.g., $z \gtrsim 3$). However, such a trend seems to be in contradiction with the observed rise in the cosmic SFR, which reaches its peak at $z \sim 2-3$ (see Madau & Dickinson 2014 for the latest review), although the latter SFR is averaged over all possible environments. Furthermore, this model is in disagreement with the measured SFR at $z \sim 6-8$ (Figure 7) due to an early gas depletion, and deviates most from the empirical sSFR-$M_*$ median shown in Figure 7.

For a more direct comparison with observations of high-$z$ galaxies, we have calculated the global star formation density in our computational box for galaxies with $M_* \gtrsim 10^9 h^{-1} M_{\odot}$ (Figure 15). We include all the wind models in the overdense and normal regions. Furthermore, we have added the best-fitting curve from Madau & Dickinson (2014), given by their Equation (15) (see also their Figure 9), as well as available observational points from Bouwens et al. (2015).

A number of conclusions can be drawn based on Figure 15. First, the overdense models lie well above the normal-density models. However, their slope is in good agreement with the observational data points, while their amplitude is above that of the data, as expected (Romano-Díaz et al. 2014; Yajima et al. 2015). Second, an extrapolation of the trend from the UCW and UVW models to higher $z$ shows that their respective curves have a slope that agrees well with the observational
points at $z \sim 7-10$, although it is somewhat steeper for $6 \lesssim z < 7$. Lastly, this agreement is better than that provided by the best-fitting curve (dashed line) at $z \gtrsim 8$. Hence, our UCR curves, although lying above the dashed curve at $z \sim 6-7$, provide a better match to observations at higher redshifts. Note that the latest corrections to the observational points (both dust-corrected and uncorrected) move them up steadily, e.g., when comparing them in Bouwens et al. (2014) with Bouwens et al. (2015). Given the uncertainties of the current observations, our small computational box, and the absence of any fine-tuning or calibration in our simulations, the match between our UCR models and observations is acceptable, which is important when one discusses the overdense regions.

An immediate consequence of the NW model and the locked-up gas content within DM halos at $z \gtrsim 6$ is related to the reionization process of the universe and the metal pollution/enrichment of the IGM. If the metal-rich gas is locked within the DM halos and is highly concentrated, the dust content within those halos and galaxies will be relatively high, making it difficult for the UV photons to escape and ionize their surroundings. The escape fraction of photons will be lower than the expected fraction of ~0.3 required to keep the universe fully ionized by $z \sim 6$ (Kashikawa et al. 2011; Finkelstein et al. 2012), resulting in a delay in the process of reionization.

In comparison, the CW model points to a relatively high-temperature IGM gas (Figure 9) that extends well beyond the halo boundaries, deep into low-density regions, where most of the IGM has reached $T \gtrsim 2-3 \times 10^4$ K by $z \sim 8$. This temperature can collisionally ionize the H I, even in the absence of UV photons. This is in contradiction with the currently accepted scenario for reionization, in which the neutral gas in the IGM is photoionized by the UV photons coming from galaxies and quasars. By $z \sim 6$, the CW model shows an overheated IGM gas, and we observe the formation of a hot bubble consisting of a shock-heated gas with a temperature up to ~$10^7$ K in the central region. This CW model gas, potentially detectable as an X-ray bubble in this overdense environment, is absent from its UCW counterpart (Figure 14). Importantly, both the CW overdense model and the UCW model exhibit large volume-filling ionized gas with $T \gtrsim 10^{4.5}$ K, which contradicts the reionization constraint.

Both NW and CW models exhibit a minimal dependence on the properties of the host galaxy. However, galactic outflows are expected to be driven by various mechanisms, including radiation pressure, SN and stellar winds, and photoheating of H II regions (e.g., Choi & Nagamine 2011; Hopkins et al. 2012; Agertz & Kravtsov 2015). The significance of each process depends on the galaxy's properties, and one should consider all mechanisms in the simulations of galaxy formation to properly account for the cumulative effects. Within this framework, the galactic outflows are difficult to describe by a simple scaling relationship. Consequently, we turn to the effects of a more complex type of model, the VW (Section 2.2).

In the VW model, the wind's properties depend on the SFR of the host galaxy, calibrated at low $z$. Although this is a rather simplistic parameterization of a very complicated process, it introduces self-regulation in the growth of the galaxy. This growth is based on the local, instantaneous SFR rather than being independent of host galaxy properties, as in the CW model.

In our simulations, the VW model exhibits the most complicated (e.g., non-monotonic) behavior in terms of interplay between galactic outflows, gas consumption, and SF among all the wind models. This behavior is especially pronounced for the gas fraction, $f_{gas,gal}$. On average, the gas fraction in VW galaxies lies between that in the gas-poor NW galaxies and that in the gas-rich CW ones, at the end of the simulations. This is also true for the median SFRs. The VW objects also show a strong decline of $f_{gas,gal}$ with the galaxy mass, $M_*$, in a sharp disagreement with the CW models.

When comparing with observational data (e.g., Figure 7), a case can be made in favor of the VW prescription. Due to the direct coupling between galaxy and wind properties imposed in the VW prescription, the SFR is correlated with the mass loss rate. This self-regulation is very important for the intermediate-mass galaxies, $10^{8.5} < M_* < 10^{10}$, which are the most abundant population at these high $z$, and probably the main contributors to reionization, but see Madau & Haardt (2015). We find that VW galaxies exhibit a tight correlation between sSFR and $M_*$, which lies within $1\sigma$ of the observationally determined median shown in Figure 7.

Observational results concerning the state of the IGM metallicity at high $z$ make it difficult to draw conclusions in favor of specific wind models in identical environments, based on the differences between their efficiency in spreading metals on large scales. Recent studies focusing on the statistics of C IV absorbers in spectra of distant quasars have shown a decline in their occurrence for $z > 5.3$ (Becker et al. 2009; Ryan-Weber et al. 2009). Ryan-Weber et al. (2009) found a cosmic mass density of C IV ions at $z \sim 5.7$ of $\Omega_{C\text{IV}} = 4.4 \times 10^{-9}$, a factor of 3.5 lower than the corresponding value for lower redshifts, and derived a corresponding lower limit for the IGM metallicity of $Z_{\text{IGM}} \gtrsim 10^{-4} Z_{\odot}$.

The enrichment levels in the IGM at $z \sim 6$ in all of our simulation runs seem to satisfy this criterion based on the maps of IGM metallicity shown in Figure 9. Even the NW model is able to reach peak metallicities of $Z \sim 10^{-3} Z_{\odot}$ in the IGM (Figure 13) despite the clear lack of metal-spread on large scales and a very small volume filling factor (Figure 9). The relatively efficient IGM pollution even in the NW is likely caused by the accelerated evolution due to the presence of the overdensity, which boosts SFR and, consequently, the metal production in this region.

It is worth noting that other simulations used different outflow prescriptions to model the observed IGM metallicity at low $z$ (e.g., Oppenheimer et al. 2011, p. 78). These of course do not include the no-wind model. Oppenheimer et al. (2011, p. 78) found that C IV absorbers are associated with $T \sim 10^4$ K gas in halos, for which, in our models, this material is found at similar metallicity levels (Figure 12).

The metal content of the IGM is affected as well by the presence of the overdensity because of the development of massive galaxies. This is in contrast with the models of average density (UCRs) where no metal-rich outflows are observed. Nevertheless, we do find that the peak metallicity attained in a few regions in the IGM of the UCW run (visible as high-metallicity spots in Figure 10) is at a similar level to the CW model. So, even though the galaxy population is less developed, the metallicity content is of the same order, albeit in limited volumes—a direct consequence of galaxy evolution. When averaged over the computational volume, the IGM
metallicity in UCR models is lower than in the CR models, with the exception of NW.

Heating the IGM to $T \gtrsim 10^5$ K will destroy the dust and facilitate reionization at much earlier times, which is in contradiction with the latest estimates of reionization time (Becker et al. 2015; McGreer et al. 2015). Furthermore, the HI opacity would decrease rapidly from even earlier redshifts, in contradiction with the results from Fan et al. (2006). At the same time, the IGM metallicity is observed to be very high in our CW models, nearly sub-solar, which should leave an imprint in the spectra of the high-z objects, such as quasars and Ly$\alpha$-emitting galaxies. However, this has not been detected yet. Furthermore, it will also help the SF process at later stages, favoring the formation of protogalaxies in the surroundings in large numbers. Therefore, galaxies, such as the so-called CR7 galaxy (e.g., Sobral et al. 2015), could be found in the vicinity of similar objects, forming protocluster regions (e.g., Trenti et al. 2012).

Although one could state that these are merely results from a high-density environment, this situation is not very different when comparing the metallicity levels between average and high-density environments. Their respective overall metallicity profiles lie within 0.5 dex, but their distributions for the different media (i.e., the IGM, halo gas, and the ISM) are very similar (Figure 13). Their overall temperature profiles are similar as well (Figure 11).

Overall, the temperature in the VW model lies within around $10^4$ K and metallicity at $log(Z/Z_\odot) \approx -2$, a much more acceptable range value, which is also in agreement with reionization. Figure 7 shows that the VW model is able to reproduce both the observed sSFR evolution between $z \sim 8$ and $z \sim 6$ and the sSFR-stellar mass relation at $z \sim 6$. The CW and NW models fail to satisfy these two constraints simultaneously. Moreover, the NW run shows a clear disagreement with observed levels of star formation at $z \sim 8-6$. At $z \sim 6$, the scaling of sSFR with stellar mass is underestimated over the entire mass range $10^{8.8}-10^{10} h^{-1} M_\odot$ spanned by the observations of Salmon et al. (2015) (Figure 7, left). The sSFR values are also off at $z \sim 8$ and $z \sim 7$, not just for the entire galaxy population, but even when considering galaxies only within the observed stellar mass range.

Unfortunately, computational expense prohibits implementing high-resolution studies in large-scale cosmological simulations. The full cosmological simulation still needs to rely on a rather simple scaling. As we discussed in Section 2.2, the VW model possesses one important feedback parameter—the dependence of wind velocity on the host galaxy SFR, which constrains the wind properties better than the CW model. Although its value is tuned to the low-z observation (Martin 2005), the median VW speed appears lower than that of the CW. As a result, the VW does not overheat the IGM, unlike the CW, while it regulates the SF and galaxy growth, and shows agreement with the sSFR observations. This implies that the scaling relation that the VW adopts still holds in galaxy evolution at $z \gtrsim 6$, and the predictions from the VW simulation, such as the GMF, SFR, sSFR, and $f_{gas,gal}$, would be more reliable than those in CW.

The low SF levels at $z < 7$ in the NW case are caused by the absence of winds, which allows efficient consumption of gas by SF at high $z$ and gradually depletes galaxies of their gas. This also causes the metallicity in galaxies to reach very high values by $z \sim 6$ for a broad range of stellar masses (Figure 8).

In particular, for low-mass galaxies, we find unrealistic metallicity levels in the stellar and gaseous components, causing the resulting mass-metallicity relation to be nearly flat in the NW case, as opposed to increasing with mass as shown by the other two wind models. Furthermore, the CW model is able to reproduce the observed sSFR but fails to capture correctly the decreasing trend of sSFR with stellar mass, because it shows a positive correlation between those two quantities. In addition, the CW model appears in contradiction with observations at lower redshifts, such as an excessive metal enrichment of low-density gas (e.g., Oppenheimer & Davé 2008; Choi & Nagamine 2011). We therefore conclude that the VW wind model more closely reproduces the properties of galactic wind, IGM, and halo gas at high $z$, although a thorough analysis using radiative transfer is needed to confirm this conclusion.

In summary, we find that galaxy evolution at $z \sim 6-12$ depends strongly on the type of feedback mechanisms applied. While specific prescriptions for galactic winds allow us to reproduce some of the observed properties of galaxies at these redshifts, the best fit is attained by using the VW model, which couples wind properties to the ongoing star formation and to the depth of the galactic potential well. No perfect wind model can be determined at this moment; perhaps a new wind model is necessary that is more tightly related to the underlying galaxy's properties. The VW model, although quite simplistic, appears best equipped to reproduce the galaxy's properties, in both overdense and average regions in the universe. Future improvements of galactic wind models must include a dependence on morphological and radiative properties of host galaxies.

We thank Volker Springel for providing us with the original version of GADGET-3. We are grateful to Ken Nagamine for fruitful discussions, to Steven Finkelstein, Richard Bouwens, and Hidenobu Yajima for helpful comments, and to Yehuda Hoffman for producing the initial conditions by means of Constrained Realizations. This work has been partially supported by the NSF grant AST-080776, the HST/STScI grant AR-12639.01-A, and JSPS KAKENHI grant \#16H02163 (to I.S.). R.S. is partially supported by NSF grant AST-1208891. I.S. is grateful for support from International Joint Research Promotion Program at Osaka University. J.H.C. acknowledges support from NASA ATP NNX11AE09G, NSF AST-1009799, and Caltech/JPL SURP Project No. 1515294 through the UT Austin (P.I. Paul Shapiro). E.R.D. thanks DFG for support under SFB 956. Support for HST/STScI grant was provided by NASA through a grant from the STScI, which is operated by the AURA, Inc., under NASA contract NAS5-26555. Simulations have been performed on the University of Kentucky DLX Cluster and using a generous allocation on the XSEDE machines to I.S.

## REFERENCES

* Adelberger, K. L., Steidel, C. C., Shapley, A. E., & Pettini, M. 2003, ApJ, 584, 45

* Agertz, O., & Kravtsov, A. V. 2015, ApJ, 804, 18

* Anglés-Alcázar, D., Davé, R., Özel, F., & Oppenheimer, B. D. 2014, ApJ, 782, 84

* Arav, N., Shlosman, I., & Weymann, R. J. 1997, in ASP Conf. Ser. Vol. 128, Mass Ejection from Active Galactic Nuclei (San Francisco, CA: ASP)

* Barkana, R., & Loeb, A. 2004, ApJ, 609, 474

* Becker, G. D., Bolton, J. S., Madau, P., et al. 2015, MNRAS, 447, 3402

* Becker, G. D., Rauch, M., & Sargent, W. L. W. 2009, ApJ, 698, 1010

* Benson, A. J., Bower, R. G., Frenk, C. S., et al. 2003, ApJ, 599, 38
Bertschinger, E. 1987, ApJL, 323, L103

Blandford, R. D., & Payne, D. G. 1982, MNRAS, 199, 883

Bouwens, R. J., Bradley, L., Zitrin, A., et al. 2014, ApJ, 795, 126

Bouwens, R. J., Illingworth, G. D., Oesch, P. A., et al. 2015, ApJ, 803, 34

Cano-Díaz, M., Maiolino, R., Marconi, A., et al. 2012, A&A, 537, L8

Castor, J., McCray, R., & Weaver, R. 1975, ApJL, 200, L107

Cecil, G., Bland-Hawthorn, J., Veilleux, S., & Filippenko, A. V. 2001, ApJ, 555, 338

Cen, R., Nagamine, K., & Ostriker, J. P. 2005, ApJ, 635, 86

Choi, J.-H., & Nagamine, K. 2009, MNRAS, 393, 1595

Choi, J.-H., & Nagamine, K. 2010, MNRAS, 407, 1464

Choi, J.-H., & Nagamine, K. 2011, MNRAS, 410, 2579

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

Combes, F., García-Burillo, S., Casasola, V., et al. 2013, A&A, 558, A124

Costa, T., Sijacki, D., Trenti, M., & Haehnelt, M. G. 2014, MNRAS, 439, 2146

Croston, J. H., Hardcastle, M. J., Kharb, P., Kraft, R. P., & Hota, A. 2008, ApJ, 688, 190

Davé, R., Hernquist, L., Katz, N., & Weinberg, D. H. 1999, ApJ, 511, 521

Davé, R., Oppenheimer, B. D., & Sivanandam, S. 2008, MNRAS, 391, 110

De Young, D. S., & Heckman, T. M. 1994, ApJ, 431, 598

Dekel, A., & Birnboim, Y. 2006, MNRAS, 368, 2

Dekel, A., & Silk, J. 1986, ApJ, 303, 39

Dressler, A., Henry, A., Martin, C. L., et al. 2015, ApJ, 806, 19

Duncan, K., Conselice, C. J., Mortlock, A., et al. 2014, MNRAS, 444, 2960

Dunkley, J., Komatsu, E., Nolta, M. R., et al. 2009, ApJS, 180, 306

Eisenstein, D. J., & Hut, P. 1998, ApJ, 498, 137

El-Zant, A., Hoffman, Y., Primack, J., Combes, F., & Shlosman, I. 2004, ApJL, 607, L75

El-Zant, A., Shlosman, I., & Hoffman, Y. 2001, ApJ, 560, 636

Emmering, R. T., Blandford, R. D., & Shlosman, I. 1992, ApJ, 385, 460

Emonts, B. H. C., Morganti, R., Tadhunter, C. N., et al. 2005, MNRAS, 362, 931

Fan, X., Carilli, C. L., & Keating, B. 2006, ARA&A, 44, 415

Fan, X., Narayanan, V. K., Lupton, R. H., et al. 2001, AJ, 122, 2833

Finkelstein, S. L., Papovich, C., Ryan, R. E., et al. 2012, ApJ, 758, 93

Fisher, D. B., & Drory, N. 2011, ApJL, 733, L47

Franx, M., Illingworth, G. D., Kelson, D. D., van Dokkum, P. G., & Tran, K.-V. 1997, ApJL, 486, L75

Fukugita, M., Hogan, C. J., & Peebles, P. J. E. 1998, ApJ, 503, 518

Heckman, T. M. 1994, in Mass-Transfer Induced Activity in Galaxies, ed. I. Shlosman, (Cambridge: Cambridge Univ. Press), 234

Heckman, T. M., Armus, L., & Miley, G. K. 1990, ApJS, 74, 833

Hilz, M., Naab, T., & Ostriker, J. P. 2013, MNRAS, 429, 2924

Hoffman, Y., & Ribak, E. 1991, ApJL, 380, L5

Hopkins, P. F., Quataert, E., & Murray, N. 2012, MNRAS, 421, 3522

Hummel, E., Beck, R., & Dettmar, R.-J. 1991, A&AS, 87, 309

Jenkins, A., Frenk, C. S., White, S. D. M., et al. 2001, MNRAS, 321, 372

Kashikawa, N., Shimasaku, K., Matsuda, Y., et al. 2011, ApJ, 734, 119

Kennicutt, R. C., Jr. 1998, ARA&A, 36, 189

Kereš, D., Katz, N., Weinberg, D. H., & Davé, R. 2005, MNRAS, 363, 2

Khochfar, S., & Silk, J. 2006, ApJ, 648, L21

Kirby, E. N., Simon, J. D., Geha, M., Guhathakurta, P., & Frebel, A. 2008, ApJL, 685, L43

Konigl, A., & Kartje, J. F. 1994, ApJ, 434, 446

Kormendy, J., & Kennicutt, R. C., Jr. 2004, ARA&A, 42, 603

Kunth, D., Mas-Hesse, J. M., Terlevich, E., et al. 1998, A&A, 334, 11

Labbé, I., Oesch, P. A., Bouwens, R. J., et al. 2013, ApJL, 777, L19

Larson, R. B. 1974, MNRAS, 169, 229

Lynds, C. R., & Sandage, A. R. 1963, ApJ, 137, 1005

Mac Low, M.-. M., & Ferrara, A. 1999, ApJ, 513, 142

Mac Low, M.-. M., & McCray, R. 1988, ApJ, 324, 776

Madau, P., & Dickinson, M. 2014, ARA&A, 52, 415

Madau, P., & Haardt, F. 2015, ApJL, 813, L8

Maller, A. H., & Dekel, A. 2002, MNRAS, 335, 487

Martin, C. L. 2005, ApJ, 621, 227

Mathews, W. G., & Baker, J. C. 1971, ApJ, 170, 241

McCarthy, I. G., Schaye, J., Font, A. S., et al. 2012, MNRAS, 427, 379

McGreer, I. D., Mesinger, A., & D'Odorico, V. 2015, MNRAS, 447, 499

Morselli, L., Mignoli, M., Gilli, R., et al. 2014, A&A, 568, A1

Murray, N., Chiang, J., Grossman, S. A., & Voit, G. M. 1995, ApJ, 451, 498
Naab, T., Johansson, P. H., Ostriker, J. P., & Efstathiou, G. 2007, ApJ, 658, 710

Oppenheimer, B. D., Davé, E., Kereš, D., et al. 2011, GAFO, 78

Oppenheimer, B. D., & Davé, R. 2006, MNRAS, 373, 1265

Oppenheimer, B. D., & Davé, R. 2008, MNRAS, 387, 577

Osterbrock, D. E. 1960, ApJ, 132, 325

Ostriker, J. P., & McKee, C. F. 1988, RvMP, 60, 1

Overzier, R. A., Guo, Q., Kauffmann, G., et al. 2009, MNRAS, 394, 577

Pettini, M., Steidel, C. C., Adelberger, K. L., Dickinson, M., & Giavalisco, M. 2000, ApJ, 528, 96

Proga, D. 2003, ApJ, 585, 406

Robertson, B., Yoshida, N., Springel, V., & Hernquist, L. 2004, ApJ, 606, 32

Romano-Díaz, E., Choi, J.-. H., Shlosman, I., & Trenti, M. 2011a, ApJL, 738, L19

Romano-Díaz, E., Hoffman, Y., Heller, C., et al. 2007, ApJ, 657, 56

Romano-Díaz, E., Shlosman, I., Choi, J.-. H., & Sadoun, R. 2014, ApJL, 790, L32

Romano-Díaz, E., Shlosman, I., Hoffman, Y., & Heller, C. 2008, ApJL, 685, L105

Romano-Díaz, E., Shlosman, I., Trenti, M., & Hoffman, Y. 2011b, ApJ, 736, 66

Ryan, R. E., Jr., Gonzalez, A. H., Lemaux, B. C., et al. 2014, ApJL, 786, L4

Ryan-Weber, E. V., Pettini, M., Madau, P., & Zych, B. J. 2009, MNRAS, 395, 1476

Salmon, B., Papovich, C., Finkelstein, S. L., et al. 2015, ApJ, 799, 183

Scannapieco, C., Tissera, P. B., White, S. D. M., & Springel, V. 2006, MNRAS, 371, 1125

Scannapieco, E., Thacker, R. J., & Davis, M. 2001, ApJ, 557, 605

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

Schaye, J., & Dalla Vecchia, C. 2008, MNRAS, 383, 1210

Schaye, J., Dalla Vecchia, C., Booth, C. M., et al. 2010, MNRAS, 402, 1536

Shapley, A. E., Steidel, C. C., Pettini, M., & Adelberger, K. L. 2003, ApJ, 588, 65

Sheth, R. K., & Tormen, G. 1999, MNRAS, 308, 119

Shlosman, I. 2013, in Secular Evolution of Galaxies, ed. J. Falcon-Barroso & J. H. Knapen, (Cambridge: Cambridge Univ. Press), 555

Shlosman, I., Vitello, P. A., & Shaviv, G. 1985, ApJ, 294, 96

Shopbell, P. L., & Bland-Hawthorn, J. 1998, ApJ, 493, 129

Shull, J. M., & Saken, J. M. 1995, ApJ, 444, 663

Sijacki, D., Springel, V., & Hernquist, L. 2009, MNRAS, 400, 100

Smail, I., Chapman, S. C., Ivison, R. J., et al. 2003, MNRAS, 342, 1185

Sobral, D., Matthee, J., Darvish, B., et al. 2015, ApJ, 808, 139

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

Song, M., Finkelstein, S. L., Ashby, M. L. N., et al. 2016, ApJ, 825, 5

Springel, V. 2005, MNRAS, 364, 1105

Springel, V., & Hernquist, L. 2002, MNRAS, 333, 649

Springel, V., & Hernquist, L. 2003, MNRAS, 339, 289

Springel, V., White, S. D. M., Tormen, G., & Kauffmann, G. 2001, MNRAS, 328, 726

Stewart, S. G., Fanelli, M. N., Byrd, G. G., et al. 2000, ApJ, 529, 201

Tasca, L. A. M., Le Fèvre, O., Hathi, N. P., et al. 2015, A&A, 581, A54

Thompson, R., Nagamine, K., Jaacks, J., & Choi, J.-H. 2014, ApJ, 780, 145

Tremonti, C. A., Heckman, T. M., Kauffmann, G., et al. 2004, ApJ, 613, 898

Trenti, M., Bradley, L. D., Stiavelli, M., et al. 2012, ApJ, 746, 55

Trenti, M., & Stiavelli, M. 2008, ApJ, 676, 767

Uhlig, M., Pfrommer, C., Sharma, M., et al. 2012, MNRAS, 423, 2374

Utsumi, Y., Goto, T., Kashikawa, N., et al. 2010, ApJ, 721, 1680

van de Weygaert, R., & Bertschinger, E. 1996, MNRAS, 281, 84

Veilleux, S., Cecil, G., & Bland-Hawthorn, J. 2005, ARA&A, 43, 769

Veilleux, S., Cecil, G., Bland-Hawthorn, J., et al. 1994, ApJ, 433, 48

Vogelsberger, M., et al. 2014, MNRAS, 444, 1518

Weiner, B. J., Coil, A. L., Prochaska, J. X., et al. 2009, ApJ, 692, 187

Weinzirl, T., Jogee, S., Khochfar, S., Burkert, A., & Kormendy, J. 2009, ApJ, 696, 411

Yajima, H., Shlosman, I., Romano-Díaz, E., & Nagamine, K. 2015, MNRAS, 451, 418