- Split View
-
Views
-
Cite
Cite
Paramita Barai, Matteo Viel, Stefano Borgani, Edoardo Tescari, Luca Tornatore, Klaus Dolag, Madhura Killedar, Pierluigi Monaco, Valentina D’Odorico, Stefano Cristiani, Galactic winds in cosmological simulations of the circumgalactic medium, Monthly Notices of the Royal Astronomical Society, Volume 430, Issue 4, 21 April 2013, Pages 3213–3234, https://doi.org/10.1093/mnras/stt125
- Share Icon Share
Abstract
We explore new observationally constrained subresolution models of galactic outflows and investigate their impact on the circumgalactic medium (CGM) in the redshift range z = 2–4. We perform cosmological hydrodynamic simulations, including star formation, chemical enrichment and four cases of supernovae-driven outflows: no wind (NW), an energy-driven constant velocity wind (CW), a radially varying wind (RVWa) where the outflow velocity has a positive correlation with galactocentric distance (r) and a RVW with additional dependence on halo mass (RVWb). Overall, we find that the outflows expel metal-enriched gas away from galaxies, significantly quench the star formation, reduce the central galactic metallicity and enrich the CGM. At z = 2, the radial profiles of gas properties around galaxy centres are most sensitive to the choice of the wind model for halo masses in the range (109–1011) M⊙.
We infer that outflows in the RVWb model are least effective, with results similar to the NW case, except that the CGM is enriched more. Moreover, we find that the models CW and RVWa are similar, both showing the impact of effective winds, with the following notable differences. RVWa causes a greater suppression of star formation rate at z ≤ 5, and has a higher fraction of low-density (δ < 10), warm-hot (104–106 K) gas than in CW. Outflows in CW produce a higher and earlier enrichment of some intergalactic medium phases than in RVWa. By visual inspection, we note that the RVWa model shows galactic discs more pronounced than all the other wind models. We predict that some observational diagnostics are more promising to distinguish between different outflow driving mechanisms in galaxies: ZC of the CGM gas at r ∼ (30–300) h−1 kpc comoving, and C iv fraction of the inner gas at r < (4–5) h−1 kpc comoving.
1 INTRODUCTION
Energy feedback from star formation (SF) and powerful supernovae (SNe) explosions are believed to eject gas from galaxies driving outflows called galactic winds1 (e.g. Johnson & Axford 1971; Mathews & Baker 1971; Larson 1974; Veilleux, Cecil & Bland-Hawthorn 2005). SNe-driven outflows are detected in many observations at low redshifts (e.g. Lynds & Sandage 1963; McCarthy, van Breugel & Heckman 1987; Bland & Tully 1988; Heckman, Armus & Miley 1990; Lehnert & Heckman 1996; Strickland et al. 2000, 2004; Rupke, Veilleux & Sanders 2005), reaching velocity as high as 1000 km s−1 (Diamond-Stanic et al. 2012), and at high-z (e.g. Franx et al. 1997; Pettini et al. 2001; Adelberger et al. 2003; Shapley et al. 2003; Erb et al. 2012; Newman et al. 2012b), extending up to distances of 60–130 kpc physical (e.g. Lundgren et al. 2012).
We use the term circumgalactic medium (CGM) to denote the gas-phase structures, excluding star-forming gas, within the virial
radius of galaxies. Galactic winds are considered to be the primary mechanism by which metals are ejected out of star-forming regions in galaxies and deposited into the CGM and the intergalactic medium (IGM; e.g. Larson & Dinerstein 1975; Aguirre et al. 2001; Scannapieco, Ferrara & Madau 2002; Aracil et al. 2004; Bouche et al. 2007; Fox et al. 2007; Kawata & Rauch 2007; Pieri, Martel & Grenon 2007; Pinsonneault, Martel & Pieri 2010; Gauthier & Chen 2012; Hummels et al. 2012). The IGM can also be metal enriched by active galactic nuclei (AGN) feedback (e.g. Khalatyan et al. 2008; Germain, Barai & Martel 2009) and heated by blazars (e.g. Puchwein et al. 2012).
SNe- and starburst-driven outflow is an important source of feedback in galaxy evolution, and constitute a key ingredient of current galaxy formation studies, both in hydrodynamic simulations (e.g. Springel & Hernquist 2003, hereafter SH03; Schaye et al. 2010) and in semi-analytical models (e.g. Benson et al. 2003; Bertone, Stoehr & White 2005). The winds act by removing gas available to make stars, hence quench SF, and are argued to suppress the formation of low-mass galaxies, flattening the low-mass end of the luminosity function in simulations (e.g. Theuns et al. 2002; Rasera & Teyssier 2006; Stinson et al. 2007), which is closer to observations. Feedback is invoked to reproduce the realistic disc galaxies (e.g. Weil, Eke & Efstathiou 1998; Sommer-Larsen, Gotz & Portinari 2003; Governato et al. 2004; Robertson et al. 2004; Okamoto et al. 2005). These outflows are also argued to affect the large column density parts of the Lyman α absorption line forest [Lyman limit systems (LLSs), sub-damped Lyman α systems (sub-DLAs) and DLAs] seen in the spectra of distant quasars, which traces the IGM matter distribution (e.g. McDonald et al. 2005; Kollmeier et al. 2006; Tescari et al. 2011; Viel, Schaye & Booth 2013).
The detailed internal physics underlying the origin and driving of galactic outflows is complex, occurring on scales proper of the multiphase structure of the interstellar medium (ISM; e.g. Heckman 2003; Stringer et al. 2012). The gas is likely accelerated either by thermal pressure (e.g. Chevalier & Clegg 1985), radiation pressure (e.g. Murray, Quataert & Thompson 2005; Sharma, Nath & Shchekinov 2011; Chattopadhyay et al. 2012; Zhang & Thompson 2012), ram pressure, cosmic rays (e.g. Samui, Subramanian & Srianand 2010; Uhlig et al. 2012) or a combination of them (e.g. Nath & Silk 2009; Everett, Schiller & Zweibel 2010; Sharma & Nath 2012). Cosmological simulations do not resolve the scales at which these physical processes actually happen. Outflows are hence implemented in a phenomenological way in the simulations using subresolution prescriptions and their interplay with the larger-scale environments is investigated.
Energy ejection by SNe and starburst in galaxies, often in the form of powerful outflows, can be incorporated in a number of ways as a source of feedback in the numerical scheme of cosmological simulations. A thermal feedback, where SNe energy is distributed as heating energy of neighbouring gas, is well known to be ineffective (e.g. Katz 1992; Friedli & Benz 1995; Steinmetz & Muller 1995; Katz, Weinberg & Hernquist 1996), because the dense star-forming gas has high cooling rate, and the injected thermal energy is radiated away quickly before it can significantly impact the gas. Therefore depositing the SNe energy in the kinetic form is a more popular implementation in the literature, which has been shown to have significant feedback effects (e.g. Navarro & White 1993; Mihos & Hernquist 1994; Cen & Ostriker 2000; Kawata 2001; SH03; Dalla Vecchia & Schaye 2008; Dubois & Teyssier 2008; Oppenheimer et al. 2012). In our current work, we adopt this kinetic feedback by giving a velocity kick to the affected gas.
Alternative forms of SNe feedback have also been implemented for example by considering that the affected gas undergoes adiabatic evolution (e.g. Mori et al. 1997; Thacker & Couchman 2000; Brook et al. 2005); by turning off radiative cooling temporarily for part of the neighbouring gas (e.g. Governato et al. 2007; Piontek & Steinmetz 2011) or by distributing SNe energy to hot and cold gas phases separately (e.g. Marri & White 2003; Scannapieco et al. 2006; Murante et al. 2010). Note that few studies have proposed prescriptions for efficient thermal feedback in smoothed particle hydrodynamics (SPH) simulations (e.g. Kay, Thomas & Theuns 2003; Dalla Vecchia & Schaye 2012).
Our kinetic feedback models (Sections 2.2 and 2.3) are based on the so-called energy-driven wind scenario originally proposed by SH03, and subsequently used by others (e.g. Tornatore et al. 2004; Nagamine et al. 2007; Dalla Vecchia & Schaye 2008; Tescari et al. 2009; Fabjan et al. 2010; Barnes et al. 2011). Here, a fraction of SNe energy provides the outflow kinetic energy, and the wind speed is independent of galaxy mass or SF rate. If we consider the underlying physics, which occurs on scales orders of magnitude below the scales resolved in cosmological simulations, such an outflow is likely to be driven by the thermal pressure of SNe, and might be (more physically) called a thermally driven wind.
Other wind models in cosmological simulations generally involve a different formulation of outflow velocity and wind mass loading factor (described more in Section 2.2) in terms of galaxy properties [mass, velocity dispersion, star formation rate (SFR)], sometimes suggested by observations. Oppenheimer & Davé (2006) implemented momentum-driven wind, (also e.g. Tescari et al. 2009, 2011), using analytical prescriptions from Murray et al. (2005) for momentum injection provided by radiation pressure of photons and SNe (or, radiatively driven wind). Okamoto et al. (2010) investigated models in which the outflow properties are dependent on the local velocity dispersion of the dark matter (DM). Choi & Nagamine (2011) developed a multicomponent and variable velocity galactic outflow model. Puchwein & Springel (2013) implemented a halo-mass-dependent energy-driven outflow model. Schaye et al. (2010) compared a number of models of both energy-driven and momentum-driven outflows: few constant-mass loading and constant-velocity cases; mass loading and velocity dependent on gas density, gravitational potential and halo circular velocity; wind particles temporarily decoupled hydrodynamically versus not decoupled and thermal injection of SNe energy.
In this paper, we explore new wind models where the galactic outflow speed is a function of radius (distance from the galaxy centre), as motivated by recent observations (Steidel et al. 2010). We investigate both the cases of a mass-independent dependence on radius, and a halo mass-dependent parametrization of the radially varying wind. To our knowledge, this is the first numerical implementation of both of these outflow models.
This paper is organized as follows: we describe our numerical code and simulation set-up in Section 2, in Section 3 we present and discuss our results, while in Section 4 we give a summary of the main findings and discuss possible future applications.
2 NUMERICAL METHOD
We use a modified version of the TreePM (particle mesh)–SPH code gadget-3, whose first version was described in Springel, Yoshida & White (2001). The public release of the code (gadget-2; Springel 2005) contains a time integration using energy and entropy conserving formulation of SPH (Springel & Hernquist 2002), uses fully adaptive smoothing lengths, and a standard SPH artificial viscosity prescription (Monaghan 1997). gadget-3 includes a more efficient domain decomposition to improve the work-load balance over gadget-2.
Some of the additional subgrid physics2 included in the semi-public version of gadget-3 code we use are outlined in Section 2.1. The general galactic wind feedback model is described in Section 2.2, our new outflow models in Section 2.3 and the code implementations in Section 2.4. The simulations we perform are presented in Section 2.5.
2.1 Subresolution physics: radiative processes, SF, chemical evolution
Radiative cooling and heating is based on the original implementation of Katz et al. (1996) and improved by adding metal-line cooling which is implemented by adopting the cooling rates from the tables of Wiersma, Schaye & Smith (2009a). Net cooling rates are computed element-by-element tracking 11 species: H, He, C, Ca, O, N, Ne, Mg, S, Si and Fe. A spatially uniform time-dependent photoionizing background radiation is considered from the cosmic microwave background (CMB) and the Haardt & Madau (2001) model for the ultraviolet/X-ray background produced by quasars and galaxies. Contributions from the 11 elements are interpolated as a function of density, temperature and redshift from tables that have been pre-computed using the public photoionization code cloudy (last described by Ferland et al. 1998), assuming the gas to be dust free, optically thin and in (photo-)ionization equilibrium.
SF and SNe feedback are implemented following the effective subresolution model by SH03. In this model, the physics of multiphase structure of the ISM, on scales unresolved in cosmological simulations, is modelled using spatially averaged properties describing the medium on scales that are resolved. Gas particles with density above a limiting threshold, nSF = 0.13 cm−3 (in units of number density of hydrogen atoms), are considered to contain two phases: cold condensed clouds that are in pressure equilibrium with ambient hot gas. Each gas particle represents a region of the ISM, where the cold clouds supply the material available for SF. Star particles are collisionless, and are spawned from gas particles undergoing SF, according to the stochastic scheme introduced by Katz et al. (1996). We allow a gas particle to spawn up to four generations of star particles; therefore, a typical star particle mass is about one-fourth of the initial mass of gas particles.
Stellar evolution and chemical enrichment feedback are implemented following the chemical evolution model of Tornatore et al. (2007). Production of nine different metal species (C, Ca, O, N, Ne, Mg, S, Si and Fe) are accounted for using detailed yields from Type Ia SN (SN-Ia), Type II SN (SN-II), along with low- and intermediate-mass stars (LIMS) in the thermally pulsating asymptotic giant branch (TP-AGB) phase. Contributions from both SN-Ia and SN-II to thermal feedback are considered. Mass-dependent time delays with which different stellar populations release metals are included, adopting the lifetime function by Padovani & Matteucci (1993). Different stellar yields are used: for SN-Ia taken from Thielemann et al. (2003), SN-II from Woosley & Weaver (1995) and LIMS from van den Hoek & Groenewegen (1997). The mass range for SN-II is considered to be M/M⊙ > 8, while that for SN-Ia originating from binary systems is 0.8 < M/M⊙ < 8 with a binary fraction of 10 per cent.
We include a fixed stellar initial mass function (IMF) according to the formalism given by Chabrier (2003), which is a power law at M/M⊙ > 1 and has a log-normal form at masses below. However, we use power-law IMFs with different slopes over the whole mass range of 0.1–100 M⊙, which mimics the log-normal form of Chabrier (2003) at lower masses, as tests indicate. In our model the functional form, ϕ(M) = KM−y, is composed of three slopes and normalizations: y = 0.2 and K = 0.497 for stellar masses 0.1 ≤ M/M⊙ < 0.3; y = 0.8 and K = 0.241 for 0.3 ≤ M/M⊙ < 1 and y = 1.3 and K = 0.241 for 1 ≤ M/M⊙ < 100. Stars within a mass interval [8–40] M⊙ become SNe first before turning into black holes (BHs) at the end of their lives, while stars of mass >40 M⊙ are allowed to directly end in BHs.
The chemical evolution model also incorporates mass loss through stellar winds and SNe explosions, which are self-consistently computed for a given IMF and lifetime function. A fraction of a star particle's mass is restored as diffuse gas during its evolution, and distributed to the surrounding gas particles. There is no AGN feedback in our simulations.
2.2 Galactic wind feedback
For our adopted Chabrier power-law IMF (Section 2.1), ϵSN = 1.1 × 1049 erg M|${^{- 1}_{\odot}}$|. We initialize the wind energy fraction using vw = 400 km s−1, which is the constant velocity of our simulation run CW (Section 2.5). It corresponds to χ = 0.29 of SNe energy being carried away by the wind.
Outflow models implemented in the gadget code generally involve different scaling relations of vw and η in terms of galaxy velocity dispersion, mass and/or SFR. For the original energy-driven outflows, vw and η are constant. The momentum-driven wind prescription by Oppenheimer & Davé (2006, 2008) has |$v_{\rm w} \propto \sigma \sqrt{(L/L_{\rm crit}) - 1}$| and η∝1/σ, where σ is the galaxy velocity dispersion, and L/Lcrit is its luminosity in units of a critical value. In Choi & Nagamine (2011)'s model, vw∝SFR1/3 and η is a function of galaxy stellar mass. Puchwein & Springel (2013) assumed vw and η to be proportional to halo mass.
However, any dependence of vw on galactocentric radius (or distance from the galaxy centre) has not been explored before and this is what we are going to do in this work.
2.3 Wind velocity dependent on galactocentric distance, |$\boldsymbol {v_{\rm w}(r)}$|
Such a trend of vw(r) correlating positively with r has also been found in other observations: for example Veilleux et al. (1994) parametrized outflow velocity as vw∝rn with 2 < n < 3 for the spiral galaxy NGC 3079.
We note that a scenario in which the wind velocity is not unique but has instead any value (fast and/or slow outflow) could possibly result in faster ejecta at large distances from the galaxy centre and slower ones closer to the centre, since the former would have had more time to travel out to large distances. The model presented here could mimic such a behaviour.
2.4 Implementation in the gadget-3 code
Our new wind model described in Section 2.3 requires an extra parameter: the distance of gas particles from their host galaxy centre. We identify galaxies by running a Friends-of-Friends (FOF) group finder on-the-fly within our gadget-3 simulations. Given a linking length, Lgr, the FOF algorithm finds a group by selecting all particles that lie within a distance Lgr from any other particle in the group.
First, an average interparticle distance, Lp, is calculated using the total mass, number of particles and the mean density. The linking length for the default on-the-fly FOF in gadget-3 is Lgr/Lp = 0.16, which is commonly used to find DM haloes. However, sometimes the FOF can link together DM substructures which actually belong to multiple galaxies with distinct stellar structures, into a single halo. Therefore, the stellar components of the FOF groups are considered to be more stable and thereby identified as galaxies.
In order to select galaxies, we obtain stellar groups by modifying the FOF group finder as follows: we allow the FOF to first link over stars as the primary particle type, and then link over (gas + DM) particles as the secondary type. To facilitate selection of stellar groups, which are generally smaller in size, we use a linking length three times smaller than the default one: Lgr/Lp = 0.16/3. The FOF is executed on-the-fly within a simulation at time intervals of a multiplicative factor 1.001 of the scale factor a, or, anext/aprev = 1.001. All the groups above a minimum length of ≥32 particles are finally selected, which provide us with a set of galaxies at a given simulation time. The code execution comprises a FOF overhead of (6–9) per cent of the total processing time.
Our preliminary tests indicate that the centre-of-mass position of the FOF groups is offset from the gas density peak location. Hence we decide not to consider the centre-of-mass positions but to record the position of the member gas particle with maximum SPH density in each group: this maximum gas density location is considered as the group centre. We tag in the code the group (representing a galaxy) and the group centre each relevant gas particle belongs to. We define the galactocentric radius, r, as the distance of each gas particle from its group centre.
In order to enable the outflow to escape from dense SF regions without affecting the SF, a new wind particle (gas particle just receiving a wind kick following equation 11) is decoupled from hydrodynamic interactions, for a maximum duration tdec = 0.025tH(z), where tH(z) is the Hubble time at a relevant simulation redshift. Then the wind particle does not enter hydrodynamic force computation in the code, but is included in gravity and SPH density calculation. If the particle's density has fallen below ndec = 0.25nSF, then the decoupling is stopped before tdec, and full hydrodynamics are enabled again. In our simulations, the majority of wind particles end their decoupling phase following the density condition.
Summarizing, the main free parameters in the wind model are vw and η, along with tdec and ρdec for the decoupling prescription. The wind velocity is either constant (equation 3, Section 2.2), or dependent on radius (equation 5) and halo mass (equation 8).
In this section, we have shown the actual implementation of different forms of vw, including radial and mass dependence of velocity which are observationally motivated. We later study the implications such models have on the CGM and IGM properties in the simulated volume.
2.5 Simulations
The series of simulations we perform are listed in Table 1. A cosmological volume, with periodic boundary conditions, is evolved starting with an equal number of DM and gas particles from z = 99, where the initial conditions have been generated using the camb3 software (Lewis, Challinor & Lasenby 2000), up to z = 2. A concordance flat Λ cold dark matter (ΛCDM) model is used with the following parameters: ΩM, 0 = 0.2711, ΩΛ, 0 = 0.7289, ΩB, 0 = 0.0463, nS = 0.96, σ8 = 0.809, H0 = 70.3 km s−1 Mpc−1, in agreement with recent observations of the CMB radiation, weak gravitational lensing, Lyman α forest and mass function evolution of galaxy clusters (e.g. Lesgourgues et al. 2007; Vikhlinin et al. 2009; Komatsu et al. 2011).
Run . | Lbox . | Galactic wind feedback . |
---|---|---|
name . | (h−1 Mpc) . | . |
Smaller box runs: SB | ||
NWt | 5 | No wind |
CWt | 5 | Energy-driven constant velocity v = 400 km s−1 |
RVWat | 5 | Radially varying with fixed parameters |
RVWbt | 5 | RVW with halo-mass-dependent parameters |
Larger box runs: LB | ||
NW | 25 | No wind |
CW | 25 | Energy-driven constant velocity v = 400 km s−1 |
RVWa | 25 | Radially varying with fixed parameters |
RVWb | 25 | RVW with halo-mass-dependent parameters |
Run . | Lbox . | Galactic wind feedback . |
---|---|---|
name . | (h−1 Mpc) . | . |
Smaller box runs: SB | ||
NWt | 5 | No wind |
CWt | 5 | Energy-driven constant velocity v = 400 km s−1 |
RVWat | 5 | Radially varying with fixed parameters |
RVWbt | 5 | RVW with halo-mass-dependent parameters |
Larger box runs: LB | ||
NW | 25 | No wind |
CW | 25 | Energy-driven constant velocity v = 400 km s−1 |
RVWa | 25 | Radially varying with fixed parameters |
RVWb | 25 | RVW with halo-mass-dependent parameters |
Run . | Lbox . | Galactic wind feedback . |
---|---|---|
name . | (h−1 Mpc) . | . |
Smaller box runs: SB | ||
NWt | 5 | No wind |
CWt | 5 | Energy-driven constant velocity v = 400 km s−1 |
RVWat | 5 | Radially varying with fixed parameters |
RVWbt | 5 | RVW with halo-mass-dependent parameters |
Larger box runs: LB | ||
NW | 25 | No wind |
CW | 25 | Energy-driven constant velocity v = 400 km s−1 |
RVWa | 25 | Radially varying with fixed parameters |
RVWb | 25 | RVW with halo-mass-dependent parameters |
Run . | Lbox . | Galactic wind feedback . |
---|---|---|
name . | (h−1 Mpc) . | . |
Smaller box runs: SB | ||
NWt | 5 | No wind |
CWt | 5 | Energy-driven constant velocity v = 400 km s−1 |
RVWat | 5 | Radially varying with fixed parameters |
RVWbt | 5 | RVW with halo-mass-dependent parameters |
Larger box runs: LB | ||
NW | 25 | No wind |
CW | 25 | Energy-driven constant velocity v = 400 km s−1 |
RVWa | 25 | Radially varying with fixed parameters |
RVWb | 25 | RVW with halo-mass-dependent parameters |
The smaller box size Lbox = 5 h−1 Mpc comoving series (or, SB box) has Npart = 2 × 1283 DM and gas particles in the initial condition. This is the higher resolution series with gas particle mass mgas = 7.66 × 105 h−1 M⊙, and DM particle mass mDM = 4.49 × 106 h−1 M⊙. The Plummer-equivalent softening length for gravitational forces is set to Lsoft = 0.98 h−1 kpc comoving. In the larger box size Lbox = 25 h−1 Mpc comoving series (or, box LB) the numbers are Npart = 3203, mgas = 6.13 × 106 h−1 M⊙, mDM = 3.59 × 107 h−1 M⊙, Lsoft = 1.95 h−1 kpc. The minimum gas smoothing length is set to a fraction 0.001 of Lsoft in all the simulations. However the minimum smoothing which is actually achieved in the simulations depends on the resolution, and in our runs the gas smoothing lengths went down to ∼0.2Lsoft. The second set of simulations, at slightly worse resolution compared to the first set, has been chosen in order to increase the statistics, since more number of haloes and more massive haloes are expected to form in the larger box.
In each series four runs are performed incorporating the same non-wind subgrid physics described in Section 2.2, and investigating different galactic wind models:
NW: no wind;
CW: energy-driven wind with constant velocity (Section 2.2, equation 3) vw = 400 km s−1;
RVWa: radially varying wind with fixed parameters (Section 2.3) rmin = 1 h−1 kpc, Reff = 100 h−1 kpc, vmax = 800 km s−1, α = 1.15 and
RVWb: radially varying wind with parameters dependent on halo mass.
The smaller box size runs have an extra ‘t’ appended to the names. Few results are presented for both the box sizes, namely SB and LB runs. While the remaining results are shown only for one box size, with the reason given in the relevant parts of Section 3.
The new velocity prescriptions are plotted in Fig. 1. The black curve in the left-hand panel shows the outflow velocity as a function of galactocentric radius for the parameters in run RVWa. For comparison, the blue horizontal line denotes the constant value of velocity, v = 400 km s−1, used in run CW. The right-hand panel depicts the parametrization of vmax, or the velocity at Reff in run RVWb, as a function of halo mass for different redshifts of collapse.
We have analysed the carbon content of the gas in the box, since carbon is one of the most abundant heavy elements in the Universe, and the spectral lines produced by ionized carbon are relatively easy to observe. The metallicity of carbon, ZC, is computed as the ratio of carbon mass to the total particle mass for each gas particle. Abundance ratios are expressed in terms of the solar metallicity, which is ZC, ⊙ = 0.002177 (mass fraction of carbon in the Sun) derived from the compilation by Asplund, Grevesse & Sauval (2005).
3 RESULTS AND DISCUSSION
3.1 Outflow velocity
The different implementations of the outflow velocity in the SB runs are shown in the panels of Fig. 2: NWt (top-left), CWt (top-right), RVWat (bottom-left) and RVWbt (bottom-right). Here the smaller box size is chosen for visualization purposes, to unambiguously distinguish the wind velocity from the remaining gas in the most massive galaxy within the box volume which has the largest SFR. The trends for the LB runs are the same as the SB results shown here; however, for galaxies of this mass range there are fewer wind particles because of lower SFR, since the LB runs produce more massive galaxies with higher SFR in a larger cosmological volume than SB. The velocity magnitude of the gas particles is plotted in Fig. 2 as a function of distance from the centre of galaxy at z = 2.44. As described in Section 2.4, galaxies are identified on-the-fly in the simulations using a FOF algorithm linking over star particles as the primary type using a linking length three times smaller than the default one, and the galaxy centre corresponds to the location of the maximum gas density.
All the gas particles lying within a radius of 100 h−1 kpc from the group centre are shown as black points in Fig. 2. The particles which have recently received a wind kick velocity can be tracked, for the duration of their hydrodynamic decoupling, through their positive decoupling time (Section 2.4). We call these the wind-phase particles and we mark them with the red plus symbols. The median wind-phase particle's velocity in radial bins is plotted as the green asterisks joined by the green dashed curve. The input subgrid wind velocity (vw) implemented in the code is instead represented by the blue solid curve. In order to mark the differences clearly, in each run results are stacked over the 20 most massive galaxies having total halo (DM + gas + stars) mass in the range Mhalo = (0.9–15) × 1010 M⊙, the exact range for each case is written in the respective panel. The escape velocity, |$v_{\rm esc} = \sqrt{2\,G M_{\rm halo} / r}$|, for a representative halo mass Mhalo = 1011 M⊙, is shown by the cyan dashed curve. At a given radius, particles moving faster than vesc are unbound to the halo.
Run NWt gives the radial structure of gas velocity without any wind, reaching a maximum of v ∼ 400–500 km s−1 around galaxies. The other three runs have winds implemented and present a distribution of red points. Wind particles are kicked over a radial range from the centre up to r ∼ (30–50) h−1 kpc, marking the size of star-forming galaxy disc. These particles gain an additional velocity of vw during their wind kick (Section 2.4). Runs CWt and RVWat demonstrate this as some red points clustered around or somewhat above the blue curve. Afterwards they move away from the galaxy centre and get decelerated, causing a fraction of the red points to fall well below vw. Therefore, when considering the wind particles at given r, the median velocity (green dashed curve) is somewhat lower than vw (blue). Run CWt has particles kicked by vw = 400 km s−1 at all r, while in our new implementation RVWat the upper envelope of the red points follow the vw(r) radial profile of equation (5). In run RVWbt (bottom right-hand panel), the blue solid curve shows vw for the halo mass Mhalo = 1011 M⊙. This turns out to be a weak wind, because most of the gas is moving at velocity above the blue curve at r ≤ 10 h−1 kpc. Hence the wind particles do not move outward appreciably, and end up following a distribution very close to the remaining gas (black points).
Fig. 2 also demonstrates that even though wind particles are kicked only one time according to the original prescription of SH03, since this procedure is done at every time-step (stochastically selecting from all the star-forming particles), the desired vw is reproduced by multiple gas particles at any given time. The wind particles slow down, but in the next time-steps more particles are kicked out with vw, as long as there is SF, thus roughly mimicking a continuous outflow of gas moving with a given subgrid input speed around a galaxy.
Velocity magnitude histograms of gas particles at z = 2.91 are plotted in Fig. 3, showing mass fraction in velocity bins; wind (gas) particles are represented by solid (dashed) curves. In each case the histogram is normalized to the total gas mass in the respective simulation. The SB runs are shown in the left-hand panel where the maximum particle velocity is ∼700 km s−1, and LB runs in the right-hand panel where particles reach ∼1100 km s−1 because of higher mass haloes forming in a larger box.
Runs CWt and CW have peaks in their wind particle velocity histogram at v ∼ 400 km s−1, since in this case vw is kept constant at this value. The wind velocity histogram is much flatter in runs RVWat and RVWa with no well-defined peak, because vw has a range of values depending on galactocentric distance. Runs RVWbt and RVWb also have a peak wind velocity, which is at a lower value v ∼ 100–200 km s−1, corresponding to the outflow of lower mass haloes which are more numerous in hierarchical structure formation. Considering all the gas particles, the NW run has slightly higher mass fraction at v ≤ 300 km s−1 (not distinguishable in the plot), while runs with wind (CW, RVWa, RVWb) have higher gas mass with v ≥ 300 km s−1.
3.2 Star formation rate
3.2.1 Global star formation rate density
The global star formation rate density (SFRD) as a function of redshift is plotted in Fig. 4, for the SB runs in the left-hand panel, and for the LB runs in the right one; in each panel the four wind models are labelled by the different colours and plotting symbols. The SFRD is computed by summing over all the SF occurring in the whole simulation box and dividing it by the time-step interval and box volume to obtain the rate density in M⊙ yr−1 Mpc−3. Galactic wind feedback clearly has significant impact on SFRD, reducing SF several times depending on the outflow model.
Between z ∼ 8 and 10, there is already a trend of suppression of SFRD with winds as the CW run has 1.2 times lower SFRD than NW, while at z < 8 the SFRDs of the runs diverge more. The differences are small between runs RVWb and NW; in RVWb the SFRD values are 1.5–2 times smaller than NW between z ≈ 3 and 8, and agree with NW at z ≈ 2–3, making RVWb the wind model least effective in suppressing SF. There is however a substantial reduction of SFRD in runs CW and RVWa. The SB box has a stronger wind-driven suppression; SFRD is 12 times smaller in RVWat (and five times smaller in CWt) than NWt at z < 3. While in the LB box, the reduction is four times in RVWa and two times in CW at z < 4. Analysing the effective wind models, RVWa has 1.2 times higher SFRD than CW from high z up to z = 5, but later RVWa suppresses it more and produces two to four times lower SFRD at z = 2.
Effects of box size and resolution are visible in our results. The z = 12 SFRD is ∼two higher in the SB box than LB, because the SB series has higher resolution, and can track denser gas, forming more stars. However at z < 3–4, the RVWat run of SB series produces two to three times lower SFRD than RVWa of LB, one reason for which is the relevant box size; the haloes forming in SB box are less massive and fewer in number than in LB. Studying the shape of the SFRD evolution, the SB runs tend to show a peak in SFRD at a certain redshift: z = 2.5 in NWt and RVWbt, z ≤ 2 in CWt and z = 4 in RVWat, whereas all the LB runs have a plateau of SFRD between z = 2 and 3.5.
Observational data are overplotted in Fig. 4 with the cyan symbols and error bars. Each data set is shown with a different plotting symbol as listed next. These data are taken mainly from Cucciati et al. (2012) – filled circles, and the compilations therein originally from Steidel et al. (1999) – asterisks, Ouchi et al. (2004) – plus signs, Perez-Gonzalez et al. (2005) – inverted triangles, Schiminovich et al. (2005) – diamonds, Bouwens et al. (2009) – open squares, Reddy & Steidel (2009) – crosses, Rodighiero et al. (2010) – open circles, van der Burg et al. (2010) – upright triangles and Bouwens et al. (2012) – filled squares.
Comparing with observations, we see that no single simulation model can fit the data from a single observation. Taken collectively, there is a better match of simulations with observations at low z. At earlier cosmic epochs between z ∼ 4.5 and 10, SFRD in the simulations is systematically higher, reaching 2–10 times the observed values. At later times between z ∼ 2 and 4.5, most of the observations lie within the ranges of SFRD produced by the different wind models. Run CW of the LB series is the model providing the best fit to the observations at lower z.
3.2.2 Specific star formation rate
The specific star formation rate (sSFR) as a function of galaxy stellar mass at z = 1.98 is plotted in Fig. 5. The SB runs are in the left-hand panel and the LB runs in the right one; in each panel the four wind models are labelled by the different colours and plotting symbols. The solid curves denote the median value within a mass bin for each run, and the grey shaded area enclose the 70 percentiles above and below the median in runs RVWat and RVWa (red curves) showing the typical scatter. The black dashed line is observational data at z = 2 from Daddi et al. (2007), indicating the main sequence (MS) for star-forming galaxies, which is parametrized by the form: SFR |$= 200 ( \frac{M_{\rm stellar}}{10^{11} \,\mathrm{M}_{{\odot }}} )^{0.9} \,\mathrm{M}_{{\odot }}$| yr−1. The black dotted line marks the loci four times above the MS. According to Rodighiero et al. (2011), the region between MS and four MS encompass a majority fraction of observed galaxies.
Comparing simulation results with observations, we see reasonably good agreement. The median sSFR of all the models in the LB volume (right-hand panel in Fig. 5) lie between the observed MS and four MS. In the SB volume (left-hand panel), runs NWt and RVWbt lie within MS and four MS; however, CWt and RVWat are below MS, for a given stellar mass.
3.3 Galaxy mass function and mass fraction
The gas and stellar mass function (dN, in units of (dlog M)−1 Mpc−3) of galaxies in the LB runs at z = 2.23 is plotted as the solid curves in Fig. 6, left-hand panels. The top left-hand panel shows the mass function of the gas component in galaxies. All the four runs have the same slope in the range Mgas ∼ (109–1011) M⊙, but shifted sideways as described next. NW and RVWb produce almost similar trends: dN ≈ 0.3 at (2–3) × 109 M⊙, with the RVWb mass function shifted rightwards by ∼1.2 mass units ( M⊙). Results in CW and RVWa are quite similar, shifted leftwards with respect to NW by (2–3) M⊙. The wind in these runs expels gas efficiently out of galaxies with Mgas ≥ 109 M⊙, producing a smaller number of objects with such high gas masses. Only at ≥1010 M⊙, CW has slightly more galaxies at the same mass than RVWa.
The bottom left-hand panel of Fig. 6 shows the mass function of the galactic stellar component. Runs CW and RVWa having dN ≈ 0.2–0.3 at Mstellar = (1.5–2) × 108 M⊙ produce a steeper stellar mass function than NW and RVWb where dN ≈ 0.1 at 109 M⊙. There is a small excess of galaxies in RVWa compared to CW between Mstellar = (2 × 108–109) M⊙, and a slight excess in NW than RVWb in the range (109–1010) M⊙. The cyan filled circles and error bars represent observational data of the stellar mass function from Marchesini et al. (2009), within 2 ≤ z < 3 (the results which are tabulated in their table 1). Our model RVWa (red curve) provides a reasonably good match to these high-z observations over the mass range Mstellar = (8 × 109–1011) M⊙. However our stellar mass functions are steeper than that of semi-analytic models by Bower, Benson & Crain (2012) at z = 0, and the observational data therein, along with other low-z observations (e.g. Papastergis et al. 2012).
Fig. 6 also shows results of the wind model RVWat from the SB box, plotted as the red dashed curve in each panel. Comparing RVWa and RVWat, the smaller higher resolution box has lower mass haloes forming, extending the gas mass function to Mgas ∼ 108 M⊙, and the stellar mass function down to Mstellar ∼ 4 × 107 M⊙.
The mass fractions of gas and stars in galaxies with respect to total mass of haloes are plotted in the right-hand panels of Fig. 6: gas mass fraction in the top-right, and stellar mass fraction in the bottom-right. Solid curves denote the median value within a mass bin for each run labelled by the colour and plotting symbol. The grey shaded area encloses the 70th percentiles above and below the median in run RVWa (red curve), and shows the typical scatter at given halo mass. With NW the gas fraction in galaxies decrease from Mgas/Mhalo ≈ 0.13 at Mhalo = 1010 M⊙ to 0.04 at 4 × 1012 M⊙. The wind (CW, RVWa, RVWb) flattens the Mgas/Mhalo trend making it oscillatory around (0.08–0.1) over the same Mhalo range. The gas fractions in galaxies are lower than the cosmic baryon to DM fraction used in our simulations (ΩB, 0/ΩM, 0 = 0.17). The gas mass fractions we obtain are roughly consistent with those found in massive cluster haloes (M500 > 1013 M⊙) at z = 0, in simulations including baryonic physics of cooling, SF and AGN feedback (Planelles et al. 2012).
The stellar fraction is largest in NW, because the SF rate is higher in the absence of wind (Section 3.2.1) producing more stars for the same halo mass. The wind in the other runs expel some star-forming gas out of galaxies, reducing the SF and consequently the stellar mass. Mstellar/Mhalo increases from ≈(0.01–0.025) at Mhalo = 1010 M⊙ to ≈(0.08–0.16) at 6 × 1011 M⊙, at a varying rate for the different runs, and remain almost flat at higher halo masses up to 4 × 1012 M⊙. NW and RVWb have the steepest increase, followed by CW, and RVWa is the flattest. Observations also find the stellar mass fractions of low-z galaxies increasing with halo mass in the range Mhalo ∼ (1010–1012) M⊙ (e.g. Papastergis et al. 2012), with a comparable slope, but the absolute values in observations are ∼10 times smaller than our simulations.
The gas and stellar fractions reveal the mass dependence of the feedback acting on the galaxies in the different wind models. In runs NW and RVWb, more massive haloes have lower gas fraction than stars because of efficient conversion of gas to stars, and less massive haloes have more gas than stars; causing Mgas/Mhalo decrease with increasing Mhalo. In the effective wind models, CW and RVWa baryons are expelled efficiently and can escape from low-mass haloes, causing the gas fraction (and star fraction) to be lower than NW and RVWb. At the high-mass end, outflows are not efficient to expel gas from haloes, causing a higher gas fraction compared to cases NW and RVWb. The result is a flatter Mgas/Mhalo fraction versus Mhalo in the effective outflow models. However the wind feedback affects the galaxy SFR; it is efficient enough and increases the gas cooling time, preventing the gas from turning into stars, causing the stellar fraction to be lower than NW and RVWb in massive haloes.
We thus conclude that the outflow models CW and RVWa have significant impact on the gas and stellar mass functions, as well as on the mass fractions of gas and stars in haloes.
3.4 Temperature–density phase diagram
The temperature versus density phase diagram of all the gas at z = 1.98 is plotted in Fig. 7, top row, for the LB runs in four horizontal panels, colour coded by the gas mass fraction. The larger box size is chosen for the two [T–ρ] diagrams we present in this section, because it gives a statistics over a larger range of haloes. The x-axis denotes the density contrast, δ = ρgas/〈ρB〉, which is the comoving density of gas ρgas as a fraction of the comoving mean baryon density, 〈ρB〉 = ΩB, 03H|${^{2}_{0}}$|/(8πG). The vertical dashed line is the SF threshold density (nSF, Section 2.1) above which a fixed effective equation of state is imposed on the gas in the multiphase SF model (SH03; Wiersma et al. 2009b). The hot-phase temperature is shown for the multiphase gas particles, which forms the constant-slope straight line above nSF in the top-right portion of each run. In the NW run (top left-hand panel), all the gas denser than nSF follows the equation of state, and there is no cold gas at those high densities.
In the runs with wind, a fraction of the gas denser than nSF is cooler (<104 K) than the NW case, forming a cold, dense tail in the bottom right-hand portion of the phase diagram. This is a consequence of the decoupling formalism of the galactic outflow implementation (Section 2.4). The wind particles are allowed to move away from SF regions and no longer follow the SF effective equation of state. They are however still very dense, where the cooling rate is high and cooling time short. Therefore they cool very fast, soon after being kicked into wind, and form the cold, dense extension above nSF in the phase diagram. Subsequently, if the kick velocity is higher than the neighbours, the wind particles exit the dense SF region soon, moving to lower density (<nSF), hotter (>104–105 K) regions. Thus in our galactic outflow models, the wind is launched from dense SF phase, and soon after that the wind goes through a cold phase.
The remaining of the phase diagram looks qualitatively similar in the four runs. The majority of the gas (red in Fig. 7) lies in underdense to moderately overdense (δ ∼ 10−1–102) and cool-warm (T ∼ 103–105 K) phase, which is the cosmic baryons in balance between adiabatic cooling and photoionization background heating. A small fraction of gas at a higher range of overdensities (δ ∼ 10−1–104) is hot (T ∼ 106–108 K), being heated by shocks during gravitational collapse. At δ ∼ 102–104, there is another cold (T < 105 K) phase which is the dense gas cooling into DM haloes to form galaxies.
The gas mass fraction in the phases varies between the four runs, forming different patterns in the phase diagram. Runs NW and RVWb have a higher fraction of gas in the multiphase SF branch of effective equation of state. Among the wind runs, RVWb has the highest mass fraction in the cold wind phase, followed by CW and finally RVWa. This is because of the different wind kick velocity in each case; RVWb has the lowest vw (seen in Fig. 2, Section 3.1) where gas is not able to escape the dense phase, while in RVWa the wind particles are able to escape the dense phase quickly. Runs CW and RVWa have a higher fraction of underdense (δ < 1), warm-hot (T ∼ 104–106 K) gas, denoted by the red and yellow contours in Fig. 7, composed of gas previously ejected from dense SF regions of galaxies as strong winds.
The temperature–density phase diagram of only the CGM gas at z = 1.98 is plotted in Fig. 7, bottom row, for the LB runs in four horizontal panels. Gas particles lying inside the virial radius R200 (computed analytically using equation 7) of all the FOF haloes are selected, excluding the star-forming ones (i.e. particles with SFR = 0 are only counted). The selected CGM gas is composed of a number fraction (0.07–0.1) of all the gas particles in a run, and is shown in Fig. 7, colour coded by the gas mass in [T–δ] bins as a fraction of the total mass plotted. The densities defining the boundaries of the CGM are denoted by the vertical dashed lines. The SF threshold density is the upper limit, above which gas forms stars and becomes part of galaxy. From Fig. 7, bottom row, we infer a lower limiting density δCGM = 5 of the CGM.
The NW run (bottom left-hand panel) has no CGM gas denser than nSF, because of imposing the SFR = 0 criterion for selecting CGM particles, as here all the gas denser than nSF is star forming. The runs with wind have a fraction of gas denser than nSF which is not star forming, and they are the wind-phase particles: mostly cold undergoing hydrodynamic decoupling, and a small fraction which recently received a wind kick and still in the hot multiphase branch. Noting the effect of feedback, we see that the CGM is separated into a cold T < 105 K phase and hotter gas in a more well-defined way in models CW and RVWa impacted by their effective outflows, than in NW and RVWb.
Overall, analysing the [T–δ] plane, we can conclude that the thermal properties of the IGM remain roughly the same between the different outflow models; with the winds producing a higher fraction of underdense, warm-hot gas. Furthermore, the CGM is well separated into cold and warm-hot phases by the winds.
3.5 Single galaxy gas kinematics
The most massive galaxies at z = 2.12 in the SB simulations are plotted in Fig. 8. The smaller box size is used here again for ease of visualization (like in Section 3.1), to clearly distinguish the outflowing gas in the plot of a massive galaxy within the box volume. Each of the four rows shows a unique gas property, the different wind model runs are labelled in the top left-hand corner of the first, third and fourth rows. The maximum gas density location is at the centre, and the panels show the projection of gas kinematics in the [x–y] plane of a surrounding (200 h−1 kpc)3 volume. The total halo mass (DM + gas + stars) of the galaxies are in the range Mhalo = (1.4–1.9) × 1011 M⊙. One can see the formation of a galaxy disc in runs NWt, RVWat and RVWbt of varying sizes, whereas in CWt the galaxy looks irregularly shaped in the projected plane. The galaxy disc in model RVWat is the largest in size, followed by NWt and then RVWbt. The discs in NWt and RVWbt are quite similar, with RVWbt producing a thicker disc.
In the top two rows the outflowing (radial velocity with respect to the galaxy centre, vr > 0) particles are denoted as red, and the inflowing (vr < 0) as black. The first row shows the positions of all the gas particles within the projected volume, and the second row depicts the velocity vectors of 10 per cent outflowing gas particles and 10 per cent inflowing gas. The number fraction of outflowing (fout) and inflowing (fin) gas are written in the top left-hand corner of the second row panels. Within this (200 h−1 kpc)3 projected volume, in runs NWt, RVWat and RVWbt, about half of the gas is inflowing and the other half outflowing, while in run CWt more gas (64 per cent) is undergoing infall than outflow (36 per cent). This is likely because of fall-back of gas which was kicked outward at earlier epochs, but did not escape the halo potential, and is infalling back. The in/out flow fractions in run CWt are similar to the result by Shen et al. (2012) who found that in the Eris2 zoom-in simulation, the CGM of a Milky-Way-type galaxy at z ∼ 3 has about one-third of all the gas within the virial radius as outflowing. In Fig. 8, the inflow is overplotted on the outflow in the first row depicting a prominence of black points over red at the centre, and vice versa in the second row showing a prominence of red points over black; this demonstrates the mixture of gas dynamics (i.e. contains both inflow and outflow) in the central 20–30 kpc regions.
The outflow is not well structured in the no-wind (NWt) and the low-velocity wind model (RVWbt) cases, some gas escaping along the direction of least resistance through a low-density void in the bottom-right half. In run CWt most of the outflowing gas lies inside r < 50 kpc, because here the wind kick velocity vw = 400 km s−1 is not enough to drive the outflow to larger distances; the strong gravitational potential of the halo causes the gas to fall back inward also resulting in a higher inflow fraction in this run. The radially varying subgrid wind velocity vw(r) of equation (5) in run RVWat produces a well-developed gas outflow propagating perpendicular to the galaxy disc, escaping to r > 100 kpc from the centre, seen as the red arrows along top-right and bottom-left in the third panel of the second row.
The third row of Fig. 8 shows the gas density contrast colour coded on a log-scale from red as the highest and black as the lowest values. All the runs have a central overdense region, where SF occurs forming a galaxy disc in some cases. Wind particles are kicked from these SF regions, carrying a fraction of the central gas to larger distances. In run CWt the constant-velocity wind forms an extended halo of gas at r ≤ 50 kpc, seen as yellow–green in the second panel, which is most likely bound to the galaxy potential and not able to escape. The strong radially varying wind in RVWat produces a more extended (r ∼ 100 kpc) but lower density diffuse outflow, which is likely escaping the galaxy in the third panel. The low-velocity wind run RVWbt has most of the gas concentrated in the central 10–20 kpc.
Carbon metallicity (ZC, defined in Section 2.5) is plotted in the fourth row of Fig. 8. Runs NWt and RVWbt have a larger central concentration of metals (red in the figure) originating from SF, because there is either no-wind or the wind is not effective to spread the metals around. The central metal content is smaller in runs CWt and RVWat, since winds carry the metals out from the SF regions and enrich the CGM and IGM. Consequently, in the second and third panels, the distribution of metals outside galaxy disc is aligned with the outflowing velocity (red arrows in the second row), and two extended gas outflow regions (density contrast in the third row) are present.
The large-scale structure of all the runs are similar as expected from the same initial condition, though there are small-scale differences as described above arising from different wind models. There is gas inflow from bottom/bottom-right regions of the panels into the centre and from top/top-right into the centre, seen as black points in the first row, black inflowing arrows in the second, overdense filaments in the third, which have lower metallicity as seen in the bottom row. These likely denote pristine gas infall along cosmological filaments that ‘feed’ the galaxy at the centre. Another common feature is the overdense gas clump 70–80 kpc to the top of the centre, most likely an infalling substructure which will eventually merge with the central galaxy.
3.6 Radial profiles around galaxy centres at |$\boldsymbol {z \sim 2}$|
The radial profiles of gas properties around centres of galaxies (found by FOF group finder) at z = 1.98 for the LB runs are plotted in Fig. 9. The larger box size is chosen for all the radial profiles we present in this section, since more number of haloes and more massive haloes form in the larger box; it hence gives a wider statistics.
Each row of Fig. 9 shows a property as a function of comoving radius (or distance from the maximum density position considered as the galaxy centre), by counting all the gas particles lying inside a distance Rlim from the centre, in the following format. The four horizontal panels denote total (DM + gas + star) halo mass (Mhalo/M⊙) ranges: 109–1010 (left) with number of haloes in the four plotted runs within the range Nhalo = 739–997; 1010–1011 (second from left) having Nhalo = 7453–7906; 1011–1012 (third from left) where Nhalo = 270–393; 1012–1013 (right) with Nhalo = 10–12. All the haloes within each mass range are stacked, and the plotted solid curve denotes the median value in radial bins for each run labelled by the colour and plotting symbol. The grey shaded area encloses the 70th percentiles above and below the median in run RVWa (red curve), showing the typical scatter at a given radius, since galaxies in general do not have spherically symmetric properties. The vertical dashed line is the halo virial radius R200 in comoving coordinates (analytical expression from equation 7) for the following masses in the horizontal panels from left: Mhalo/M⊙ = 3 × 109, 3 × 1010, 3 × 1011 and 3 × 1012, where the exact values are R200 = 32.2, 69.4, 149.6 and 322.2 h−1 kpc, respectively. The outer plotting radius Rlim is chosen to be 300 h−1 kpc for Mhalo/M⊙ = 109–1010. In the other halo mass ranges Rlim is scaled up according to the virial radius, making Rlim = 300, 646.6, 1393.8 and 3001.9 h−1 kpc.
3.6.1 Density
The gas density contrast (ratio of comoving density to the comoving mean baryon density, Section 3.4) radial profiles are plotted in the top-most first row of Fig. 9. Within the approximate virial radius r < R200, all the profiles are composed of two negative-sloped functions separated by a threshold radius dependent on halo mass and wind model. At r > R200 the density profiles tend to rise again and fall, forming a local peak. This occurs because of the presence of other smaller haloes and substructures, giving rise to a local density peak. For each halo mass range a sufficiently large volume within r ∼ 10R200 is plotted to reach surrounding structures. This trend is most prominent in the lowest mass haloes, and decreases at higher masses.
The inner parts r < 10 h−1 kpc of the two lower halo mass ranges (109–1010 and 1010–1011, left two panels) present significant differences: CW and RVWa produce a lower density than the NW and RVWb cases, because the strong wind is able to expel gas from the star-forming regions, consequently reducing the central gas density by 10–30 times. The trend is almost reversed in the outer parts (r > R200): CW, RVWa and RVWb have a somewhat higher density than NW, because of accumulation of gas expelled by wind.
The differences are smaller in the two higher halo mass ranges (1011–1012 and 1012–1013, right two panels), because the wind, not being dependent on halo mass in runs CW and RVWa, is less effective in ejecting gas out of massive galaxies. CW produces slightly higher density (by two to three times) than RVWa in the inner r < 10 h−1 kpc, while outside this distance they are quite similar. Comparing our density radial profiles in the halo mass range 1012–1013 (top right-hand panel in Fig. 9), with that of Hummels et al. (2012) Fig. 5, we find qualitative agreements; however, our density profiles are steeper.
3.6.2 Temperature
The temperature radial profiles are presented in the second row of Fig. 9, where the hot-phase temperature has been used for those gas particles which are multiphase (star forming). In the inner parts r ≤ 10 h−1 kpc of the galaxies, the T profiles follow the negative-sloped density profiles (first row of Fig. 9), for all the runs in the left three panels, and run NW in the right-hand panel. This represents the dense gas near the galaxy centre undergoing SF, and having a high temperature (∼105–107 K) as a result of following the SF effective equation of state (Section 3.4). In the wind runs there is also the presence of cooler (<104 K) gas, which has recently received a wind kick and undergoes hydrodynamic decoupling (Section 2.4). This gas was seen as the cold, dense tail in the bottom right-hand portion of the [T–ρ] phase diagram (Fig. 7, Section 3.4). Consequently, runs CW, RVWa and RVWb present a bimodal T distribution at r ≤ 10 h−1 kpc in some panels, composed of hot multiphase SF and cold wind.
There is a change in T slope in the outer parts, at r ≥ 6 h−1 kpc in the left two panels and at r ≥ (20–30) h−1 kpc in the right two panels, when the gas T increases with radius, likely because of shock heating at galaxy outskirts. Between (10 and 300) h−1 kpc, the NW model results in the highest T (∼105–106 K), followed by RVWb, then RVWa and CW. Our temperature profiles show a bump at 200–300 kpc, which is at a larger r than the peak of Hummels et al. (2012).
3.6.3 Carbon metallicity
Radial profiles of carbon metallicity are plotted in the third row of Fig. 9, showing the ratio of carbon mass fraction in the gas to that of the Sun. The solid curve medians and percentile values are computed considering all (both enriched and non-enriched) gas particles in radial bins. The dashed curves represent median ZC for the enriched particles only, i.e. those having ZC > 0, without counting particles with ZC = 0.
We also checked the gas fraction (by particle number) which has ZC > 0. All the gas is enriched inside a limiting radius, which is ∼7, 10, 40 and 120 h−1 kpc in the panels from left, defining the size of the central SF region or the sites of metal generation. The enriched fraction decreases outside the limiting radius, with a trend depending on halo mass and wind model: it reduces to (0.5–0.8) at R200 and (0.45–0.7) at 300 h−1 kpc. Runs CW and RVWa enrich a higher (up to 1.4 times) fraction of gas than NW and RVWb.
The dashed median ZC in Fig. 9, third row is indistinguishable from the solid median ZC in the inner parts, since most of the gas is enriched. While at large r the dashed median ZC is higher, because the contribution of the non-enriched (ZC = 0) particles reduces the solid median ZC, as can be seen in all the four panels. Especially in the second and third panels the enriched median ZCs are 100–5000 times higher than the total median ZCs.
Some features of ZC are similar to the gas density profiles (Section 3.6.1), because metals are produced during SF which occurs in dense regions. At r < (0.4–0.6)R200 all the profiles show decreasing ZC going outward from centre, with varying r-dependent negative slopes, which we later fit by a second-order polynomial. The ZC profiles start to rise again from r ≥ (0.4–0.6)R200 and fall at larger r, forming a local peak at r > R200. Such a trend is visible in the left two and rightmost panels for both solid and dashed curves, and in the other panels for the dashed curves only. It occurs because of a combination of reasons: the presence of surrounding substructures where more metals are produced in situ by ongoing SF, and the spreading of metals by wind from the central SF regions into the CGM. These simulation results are consistent with observations which show breaks (changes of slope) in the radial metallicity profiles, and/or rising metallicity gradients in the outer regions of galaxies (e.g. Scarano & Lepine 2013).
In the inner r < (5–7) h−1 kpc CW and RVWa profiles have a lower ZC than NW and RVWb, because wind suppresses central SF and transports some metal out. The trend reverses in the outer r > (5–7) h−1 kpc: CW, RVWa and RVWb produce a metallicity about 20–30 times higher than NW, because of accumulation of metal-enriched gas expelled by wind. The differences are most prominent in the lower halo mass ranges (109–1010 and 1010–1011, left two panels), and decreases at higher masses. Massive haloes of 1011–1012 (third from left) still show significantly different profiles at large r, while the profiles of the four wind models in the most massive haloes (right) are very similar.
Model RVWb produces noteworthy changes in the ZC profile compared to NW: metallicity in the central part is slightly lower in RVWb reaching about half of NW, but starts to become larger at r > 5 h−1 kpc up to 10 times higher than NW. This shows that the wind in RVWb, though least effective in other aspects, is substantially effective in transporting metals away from SF regions into lower density surroundings.
We infer from Fig. 9, third row that the CGM gas at galactocentric distances close to and beyond R200, within r ∼ (30–300) h−1 kpc comoving, around galaxies of masses Mhalo/M⊙ = 109–1011, can give the best ZC observational diagnostic to distinguish between different galactic outflow models. Our metallicity dashed curves (enriched-particle only median, which is analogous to mass-weighted metallicity) are comparable to that of Hummels et al. (2012).
3.6.4 C iv fraction
We compute the fraction of triply ionized carbon, C iv, of the gas particles in post-processing using photoionization tables derived from cloudy (last described by Ferland et al. 1998). The relevant quantities used are redshift, gas density and temperature from simulation snapshot, and an ionizing background from the HM05 tables (Haardt & Madau 2001), as included in version 07.02.00 of cloudy.
Fourth (bottom) row of Fig. 9 shows the radial profiles of the gas C iv fraction, |$f_{\rm C\,\small IV}$|. There are varying positive slopes at small r where |$f_{\rm C\,\small IV}$| increases with radius. |$f_{\rm C\,\small IV}$| attains a peak at an intermediate r, which occurs at (20–60) h−1 kpc depending on halo mass and wind model. Further out |$f_{\rm C\,\small IV}$| decreases from its peak at larger r, and have a negative slope up to 300 h−1 kpc.
In the left two panels, runs CW and RVWa produce higher |$f_{\rm C\,\small IV}$| than NW and RVWb; the differences are significant at r ≤ 7 h−1 kpc reaching 102–104 times, very small at larger r up to the peak of |$f_{\rm C\,\small IV}$|, and increases again to 10–70 times at r > 20 h−1 kpc. The small r positive-sloped regions in the right two panels show small differences: in the third RVWa has higher |$f_{\rm C\,\small IV}$| than NW and RVWb, which are in turn higher than CW, whereas in the fourth (right) NW has higher |$f_{\rm C\,\small IV}$| than CW, RVWa and RVWb. In the third panel |$f_{\rm C\,\small IV}$| in runs CW and RVWa reaches a ∼50 times higher peak at a larger r than runs NW and RVWb. The large r negative-sloped regions shows the same trends as the left two panels, runs CW and RVWa producing higher |$f_{\rm C\,\small IV}$| than NW and RVWb.
The overall behaviour of |$f_{\rm C\,\small IV}$| arises from the combined effects of density and temperature radial dependence of the gas (from first and second rows in Fig. 9) in the ionization tables, with the major role played by the temperature.
As a prediction for C iv fraction observations, we infer from Fig. 9, bottom row that the inner gas at galactocentric distances r < (4–5) h−1 kpc comoving, around galaxies of masses Mhalo/M⊙ = 1010–1011, can most effectively distinguish between strong-wind and no-wind cases.
3.7 Metallicity–density relation
The carbon metallicity as a function of gas density contrast is plotted in Fig. 10, at z = 1.98 in the top row, for the LB runs labelled by the colour and plotting symbol. The larger box size is selected again for increased statistics, since more number of haloes and more massive haloes form in it. Gas particles in three temperature bins are shown: 104–105 K in the top-left panel, 105–107 K in the top-middle and all the gas in the top-right. The solid curves denote the median ZC in δ-bins for each run. The grey shaded area enclose the 70th percentiles above and below the median in run RVWa (red curve), showing the scatter at a given density. Following a format similar to Fig. 9, the solid curve medians and percentiles are computed considering all (both enriched and non-enriched) gas particles, while the dashed curves represent median ZC for the enriched particles only, i.e. those having ZC > 0. Note that if ≥50 per cent of the particles are not enriched, then the solid median would be zero.
The dashed median ZC in Fig. 10, top row, coincide with the solid median ZC at high densities, δ ≥ 103–104, because most of the dense gas is star forming, generating metals and hence enriched. While at δ < 103 the dashed median ZCs are several orders of magnitude higher than the solid ones (most clear in the top left- and top right-hand panels), implying that most or all of the low-density and underdense gas is not enriched.
As an exception, the underdense gas having 0.1 ≤ δ ≤ 2 in the warm-hot 105–107 K phase (top-middle panel) is significantly enriched by the winds in runs CW and RVWa, so that all the particles (solid curves) have a median of ZC/ZC, ⊙ ∼ 0.05–0.08, just few times below the enriched-only dashed medians. Further trends discussed below analyses the dashed medians wherever the solid medians are below the plotting range.
Overall the metallicity–density relation shows a negative correlation at small δ and positive correlation at large δ. Similar [ZC–δ] relations were also found by Tornatore et al. (2010) for the average metallicity of the warm-hot IGM at z = 0. We find that the threshold density for the slope turnover depends on temperature bin and wind model, lying between δ = 30 and 1000. As a special feature, run CW has a slight increase in ZC from the smallest δ to δ ∼ 10 in all the temperature bins, while run RVWa first decreases and then remains constant over the same range of densities.
The wind models produce significantly different ZC in the low-density IGM at δ < 102–103, with runs CW and RVWa enriching up to few 100 times higher than NW and RVWb. The metallicity induced by the four runs becomes similar at δ > 104 in regions composed of SF gas. The warm phase of 104–105 K (top-left) presents the largest differences between the wind models: CW enriches the δ ∼ 30 gas about 10 times more than RVWa, which is 10 times more than RVWb, which in turn is 50 times more than NW.
RVWa produces higher ZC than CW at δ < 5 (top left-hand panel), but for larger δ values CW enriches more. The same trend (RVWa enriching more than CW at small δ) is visible in the other two panels (top-middle and top-right), but the resulting metallicity becomes almost comparable in RVWa and CW at larger densities.
Applying the pixel optical depth technique to quasar spectra, Schaye et al. (2003) found a positive gradient of observed CIV metallicity with density contrast, measuring a median [C/H]∝0.65log δ, in the range log δ = [−0.5, 1.8] and z = 1.8–4.1. Comparing with our Fig. 10, dashed curves in the top row, some of the wind runs show a weak similar trend in this δ range. In all the top panels of run CW and top-middle panel of RVWa, ZC show a weak positive correlation with δ, but with a slope ∼four times smaller than observed. On the other hand, runs NW and RVWb, top left- and top right-hand panels of RVWa, along with the all-particle solid medians of runs CW and RVWa in the top-middle panel, show a negative correlation of ZC with δ, or almost constant ZC. At the same time, this negative correlation is consistent with the results of Barai, Martel & Germain (2011). Studying IGM enrichment from anisotropic AGN outflows, they found that at z ≥ 2 the underdense regions (δ < 1) are enriched to higher metallicities, and the resulting [O/H] decreases with increasing IGM density, a trend more prominent with increasing anisotropy of the outflows.
3.8 Redshift evolution: metallicity at |$\boldsymbol {z \sim 4}$|
We investigate the redshift evolution of two metallicity correlations: ZC versus r, and ZC versus δ, of the LB runs.
The gas carbon abundance radial profiles of galaxies at z = 3.98 are plotted in Fig. 11, which is an earlier epoch than Fig. 9. The three panels in Fig. 11 denote three total halo mass ranges, Mhalo/M⊙ = 109–1010 (left-hand panel) with number of haloes in the four plotted runs between Nhalo = 800–1034, 1010–1011 (middle) having Nhalo = 4921–5336, and 1011–1012 (right) where Nhalo = 62–100. Because of the earlier time, here no halo has grown to the highest bin mass of the previous radial profile plots at z = 1.98 (Section 3.6).
We see that several of the median ZC versus r correlations are similar between z = 3.98 and 1.98, only the absolute metallicity values are lower at the earlier epoch. All the runs at z = 3.98 have ZC decreasing with r at r ≤ (0.5–2)R200, with r-dependent negative slopes. The trend of ZC profiles rising again at larger r values (because of in situ metal generation in surrounding structures and spreading of metals by wind from the central SF regions) is weaker at z = 3.98; it is only visible for the dashed curves in all the panels and few solid curves in the left-hand panel.
The differences between the wind models are smaller at z = 3.98 than z = 1.98; however, there are signatures of suppression of central SF by the impact of winds, which transport metals out and accumulate them in the CGM. Analysing the solid curves in the middle and right-hand panels, runs CW and RVWa have few times lower ZC than NW and RVWb in the inner r < (5–7) h−1 kpc. The trend reverses in the outer r > (5–7) h−1 kpc; CW, RVWa and RVWb have few times higher metallicity than NW. The enriched-only dashed curves show a feature different from their solid counterparts in the left and middle panels, for large r values, where they detach from the solid curves and increases with r; ZC in runs NW and RVWb are ∼two to three times higher than CW and RVWa.
The bottom row of Fig. 10 shows the gas carbon metallicity versus density contrast at z = 3.98, for three temperature ranges, following the top row. Analogous to [ZC–r], several of the median ZC versus δ correlations at z = 3.98 are similar to those at z = 1.98, with the absolute metallicity values being lower at earlier times. The differences between the wind models are smaller at z = 3.98 than z = 1.98.
The underdense gas with 0.1 ≤ δ ≤ 0.8 in the warm-hot 105–107 K phase (bottom-middle panel) is already enriched by the winds in run CW, making the all-particle (solid curve) median ZC/ZC, ⊙ ∼ 0.002 at z = 3.98, which is 40 times lower than at z = 1.98. However, the enrichment at similar δ values caused by the winds of model RVWa at z = 3.98 is substantially lower than that of CW by 1000 times or more.
Further trends at z = 3.98 discussed below analyse the dashed medians wherever the solid medians are below the plotting range. Starting from the lowest δ, the metallicity–density in all the outflow runs (except RVWb in the bottom-middle panel) shows a shallow positive correlation at δ ≤ 30. At z = 1.98, such a feature is only present in run CW. This positive-sloped [ZC–δ] relation that we obtain at z = 3.98 in the wind models is consistent with observations by Schaye et al. (2003) who found a positive gradient of C iv metallicity at similar overdensities.
4 SUMMARY AND CONCLUSION
We explore new models of galactic winds performing hydrodynamic simulations using the TreePM–SPH code gadget-3, and analyse their impact on the properties of the CGM at z = 2–4. Our outflow implementation imparts kinetic feedback in the energy-driven formalism where the wind velocity has a positive correlation with galactocentric radius vw(r) as seen in observations by Steidel et al. (2010). We further investigate a halo-mass-dependent parametrization of the radially varying wind, following observations by Martin (2005). The simulations include additional subgrid physics: metal-dependent radiative cooling and heating in the presence of photoionizing background radiation; star formation; stellar evolution; and self-consistent chemical enrichment using a fixed stellar IMF.
The implementation of the new wind models in the code involves finding the distance of gas particles from their host galaxy centre. We identify galaxies by running a FOF group finder on-the-fly within a simulation at intervals of 1.001 times the scale factor, and find stellar groups of at least 32 particles, by linking over stars as the primary particle type, using a linking length three times smaller than that for obtaining DM haloes. The location of the member gas particle with maximum density is considered as the galaxy centre. Multiphase gas particles undergoing SF are stochastically selected and kicked into wind by giving their speed a one-time vw boost in a direction perpendicular to the galaxy disc. Wind particles are also temporarily decoupled from hydrodynamical interactions.
We simulate two different cosmological volumes (in order to increase the statistics) a smaller box of (5 h−1 Mpc)3 comoving with 2 × 1283 DM and gas particles, and a larger box of (25 h−1 Mpc)3 with 2 × 3203 particles using the flat ΛCDM concordance model. For each volume, we perform four runs investigating different galactic wind models:
NW: no wind;
CW: energy-driven wind with constant vw = 400 km s−1;
RVWa: radially varying wind with fixed parameters;
RVWb: RVW with parameters dependent on halo mass.
The main analyses of our simulations reflect the following processes: the outflows expel gas away from the star-forming galaxies and suppress the SF; the gas is also metal enriched; thus winds carry metals out and accumulate them in the CGM and IGM, enriching these lower density regions with metals. Our results are summarized below.
Outflow speed
The outflow gas velocity magnitude as a function of galactocentric distance of multiple wind-phase gas particles, obtained in our simulations, follows the given input subgrid wind speed, constant or varying with radius, in agreement with observations.
Star formation rate
Galactic wind feedback quenches SF at z < 8. In the (25 h−1 Mpc)3 runs at z ∼ 2, the global SFRD is four times smaller in RVWa and two times smaller in CW than NW. At z ≤ 5 RVWa causes a greater suppression than CW, and produces two to four times lower SFRD than CW at z = 2. The SFRD at z ∼ 4.5–10 is systematically (up to 2–10 times) higher in the simulations, compared to observations. At later epochs z ∼ 2–4.5, most of the observations (e.g. Cucciati et al. 2012) lie within the ranges of SFRD produced by the different wind models.
The sSFR versus galaxy stellar mass at z ∼ 2 presents reasonably good agreement with observations, the models reproducing the observed main sequence (e.g. Daddi et al. 2007) for star-forming galaxies.
Gas and stellar mass functions and mass fractions
The gas mass function of galaxies in the (25 h−1 Mpc)3 runs at z = 2.23 has the same slope between Mgas ∼ (109–1011) M⊙, but in CW and RVWa is shifted leftwards with respect to NW and RVWb by (2–3) M⊙; because winds expel gas, causing a smaller number of objects with high gas masses. Runs CW and RVWa produces a steeper stellar mass function than NW and RVWb. There is a small excess of galaxies in RVWa compared to CW between Mstellar = (2 × 108–109) M⊙. Our model RVWa provides a reasonably good match to the observational data of the stellar mass function at 2 ≤ z < 3 (e.g. Marchesini et al. 2009), over Mstellar = (8 × 109–1011) M⊙.
In run NW the gas mass fraction decreases monotonically with halo mass; outflows (runs CW, RVWa, RVWb) flatten the Mgas/Mhalo trend making it oscillatory. The stellar mass fraction in galaxies is the largest in run NW, and increases with Mhalo, such that NW and RVWb have the steepest increase, followed by CW, and RVWa is the flattest.
Thermal state of the gas
The [T–ρ] phase diagram shows that the outflow, soon after leaving the dense SF phase, goes through a cold (<104 K) phase, due to the hydrodynamic decoupling. The thermal properties of the IGM look qualitatively similar in the four runs. Model RVWa has a higher fraction of underdense to slightly overdense (δ ∼ 0.1–10), warm-hot (T ∼ 104–106 K) gas, than CW; both of them are significantly higher than NW and RVWb.
Gas kinematics
Projection of gas kinematics around the centre of a galaxy of total halo mass Mhalo ∼ 2 × 1011 M⊙ at z = 2.12 shows that for NW and RVWb, the gas outflows along the direction of least resistance of a low-density void. Run CW has most of the outflow inside r ≤ 50 kpc, forming a metal-enriched gas sphere, vw = 400 km s−1, not being enough to drive the wind to larger distances. The vw(r) in run RVWa produces an extended diffuse enriched gas outflow propagating perpendicular to the galaxy disc (as expected from the model input), escaping to r > 100 kpc. The formation of a galaxy disc is visible in runs NW, RVWa and RVWb, where inside the (200 h−1 kpc)3 projection volume about half of the gas is inflowing and the other half outflowing. The CW galaxy looks irregularly shaped, and more gas (64 per cent) is undergoing infall than outflow (36 per cent). The galaxy disc in run RVWa is the largest in size and further investigations are needed to check if vw(r) can produce realistic disc galaxies.
Radial profiles
The radial profiles of gas properties around galaxy centres at z = 1.98 show most prominent differences between the models for lower halo masses (Mhalo/M⊙ = 109–1010 and 1010–1011), and almost uniform results at higher masses. The density and carbon metallicity profiles often form a local peak at r > R200 of galaxies, because of the presence of smaller haloes and surrounding substructures where more metals are produced by in situ SF, and the spreading of metals by wind.
Gas density profile. Runs CW and RVWa have a lower density in the inner r ≤ 10 h−1 kpc, by 10–30 times, than the NW and RVWb cases. In the outer r > R200, CW, RVWa and RVWb show a higher density than NW.
Temperature profile. The wind runs CW, RVWa and RVWb present a bimodal temperature distribution at r ≤ 10 h−1 kpc, composed of hot multiphase star-forming gas and cold winds. Between (10 and 300) h−1 kpc, the NW model causes the highest T (∼105–106 K), followed by RVWb, then RVWa and CW.
Carbon metallicity profile. Runs CW and RVWa have a lower ZC than NW and RVWb in the inner r < (5–7) h−1 kpc, while in the outer parts runs CW, RVWa and RVWb produce higher ZC by 20–30 times than NW. Metallicity in RVWb becomes larger than in NW at r > 5 h−1 kpc and up to 10 times higher at r ≥ R200. We perform second-order polynomial fits to ZC versus r, whose coefficients are listed in Table A1.
The carbon enriched (ZC > 0) gas fraction is 1 in the inner regions of galaxies; it reduces to (0.5–0.8) at R200 and to (0.45–0.7) at 300 h−1 kpc. Runs CW and RVWa enrich a higher (up to 1.4 times) fraction of gas than NW and RVWb.
C ivfraction around galaxies. C iv fraction profiles have varying positive and negative slopes at different r, attaining a peak between (20 and 60) h−1 kpc. Runs CW and RVWa produce higher |$f_{\rm C\,\small IV}$| than NW and RVWb. We fit the |$f_{\rm C\,\small IV}$| radial profiles with two first-order polynomials, and list the coefficients in Table A2. The results for Mhalo = (1010–1011) M⊙, fitted between r = (1–6) h−1 kpc, show the largest difference: |$f_{\rm C\,\small IV} \propto r^{6}$| for NW, and |$f_{\rm C\,\small IV} \propto r^{1}$| for RVWa.
Observational predictions. Inferred from our simulations, we predict that ZC observations of the CGM gas at galactocentric distances in the range r ∼ (30–300) h−1 kpc comoving, around galaxies of Mhalo = (109–1011) M⊙, can best distinguish between different galactic outflow scenarios. And C iv fraction observations of the inner gas in the range r < (4–5) h−1 kpc comoving, around galaxies of Mhalo = (1010–1011) M⊙, can most effectively distinguish between strong-wind and no-wind cases.
Metallicity as a function of density
The underdense δ = 0.1–2 warm-hot 105–107 K phase is significantly enriched by the winds in runs CW and RVWa, so that the all-particles median is ZC/ZC, ⊙ ∼ 0.02–0.08. ZC versus δ at z = 1.98 shows a negative correlation at small δ, and positive correlation at large δ. The only exception is run CW where ZC has a small increase between δ ∼ 0.1 and 10, while ZC in RVWa first decreases and then remains flat. The low-density IGM with δ < 102–103 is significantly enriched in runs CW and RVWa, up to few 100 times more than NW and RVWb. In the warm 104–105 K phase, CW enriches the δ ∼ 30 IGM 10 times more than RVWa, which enriches 10 times more than RVWb, which in turn enriches 50 times more than NW.
We perform two first-order polynomial fits to the [ZC–δ] correlation, whose coefficients are listed in Table A3. As an example, for the warm-hot (105–107) K gas, the fits within 102 ≤ δ < 104 are ZC∝δ0.656 for NW, and ZC∝δ0.546 for RVWa.
Redshift evolution
Several of the [ZC–r] and [ZC–δ] correlations at z = 3.98 are similar to those at z = 1.98, with lower metallicity values and smaller differences between the wind models at the earlier epoch. The trend of ZC radial profiles rising again at large r is almost absent at z = 3.98. For Mhalo/M⊙ = 109–1010 and 1010–1011, the enriched-only median ZC in runs NW and RVWb at large r are two to three times higher than in CW and RVWa.
The z = 3.98 underdense δ = 0.1–0.8 warm-hot 105–107 K gas is enriched by the winds in run CW (but not in RVWa) making the all-particle median ZC/ZC, ⊙ ∼ 0.002, which is 40 times lower than at z = 1.98. The [ZC–δ] relation at z = 3.98 in all the wind runs (except RVWb in the warm-hot phase) shows a shallow positive correlation over δ ∼ 0.1–30, consistent with observations (Schaye et al. 2003).
In summary, we have found that the wind model with the radially varying velocity dependent on halo mass (RVWb) is the least effective in modifying IGM related properties, with results similar to the NW case, except that it substantially enriches the low-density CGM.
The impact of the model RVWa, for which the velocity is increasing as a function of galactocentric distance, is instead similar to the energy-driven constant-velocity implementation CW. However, it shows interesting differences that deserve to be further investigated. RVWa causes a greater suppression of SFR than CW at z ≤ 5, this could have implications for the galaxy downsizing scenario. RVWa also produces galactic discs larger than all the other wind models, and one can study if the radially varying outflow formalism can generate more realistic disc galaxies. Run RVWa has a higher gas fraction than run CW in the low-density (δ ∼ 0.1–10) warm-hot (104–106 K) phase of the IGM, which could shed light on the missing baryon problem.
We see different trends of ZC versus δ between the CW and RVWa models at δ ≤ 10, and CW outflows generally produce a higher and earlier enrichment of some IGM phases than RVWa. To explore such IGM metal-enrichment differences, future progress in this field should include computing more observable statistics from the simulations, e.g. Lyman α flux, simulated quasar spectra, and compare them with observations of CGM and IGM at different impact parameters from galaxies.
We are grateful to Volker Springel for allowing us to use the gadget-3 code. Calculations for this paper were partly performed on the COSMOS Consortium supercomputer within the Dirac Facility jointly funded by STFC, the Large Facilities Capital Fund of BIS and the University of Cambridge, as well as the Darwin Supercomputer of the University of Cambridge High Performance Computing Service (http://www.hpc.cam.ac.uk/), provided by Dell Inc. using Strategic Research Infrastructure Funding from the Higher Education Funding Council for England. Simulations were also run at the CINECA Super Computer Center (CPU time assigned through an INAF-CINECA grant). We thank Olga Cucciati for sending us observational data for SFRD, and Martin Haehnelt, Giuseppe Murante for useful discussions. This work is supported by PRIN-MIUR, PRIN-INAF 2009, INFN/PD51 grant. PB and MV are supported by the ERC Starting Grant ‘cosmoIGM’. ET acknowledges the Australian Research Council Centre of Excellence for All-sky Astrophysics (CAASTRO), funded by grant CE110001020. MK acknowledges a fellowship from the European Commission's Framework Programme 7, through the Marie Curie Initial Training Network CosmoComp (PITN-GA-2009-238356). PM has been partially supported by a FRA2009 grant from the University of Trieste.
Here, we use the terms wind and outflow synonymously, meaning continuous outward flow of gas from a galaxy, which might or might not escape depending on its velocity and galactic potential.
By subgrid we mean subresolution, referring to physical processes occurring at length scales smaller than the resolved scales in our simulations.
REFERENCES
APPENDIX A:
Mhalo, min . | Mhalo, max . | rmin . | rmax . | Run . | Fit coefficients . | ||
---|---|---|---|---|---|---|---|
(M⊙) . | (M⊙) . | (h−1 kpc) . | (h−1 kpc) . | name . | A . | B . | C . |
109 | 1010 | 1 | 30 | NW | 0.578 | −0.685 | −1.60 |
109 | 1010 | 1 | 30 | RVWa | −0.519 | 1.08 | −1.85 |
1010 | 1011 | 1 | 40 | NW | −0.209 | 2.23 | −3.32 |
1010 | 1011 | 1 | 40 | RVWa | −0.180 | 0.487 | −1.54 |
1011 | 1012 | 1 | 300 | NW | 0.227 | 0.583 | −1.03 |
1011 | 1012 | 1 | 300 | RVWa | −0.265 | 1.02 | −1.03 |
1012 | 1013 | 1 | 300 | NW | 0.373 | −0.0354 | −0.384 |
1012 | 1013 | 1 | 300 | RVWa | 0.498 | −0.241 | −0.255 |
Mhalo, min . | Mhalo, max . | rmin . | rmax . | Run . | Fit coefficients . | ||
---|---|---|---|---|---|---|---|
(M⊙) . | (M⊙) . | (h−1 kpc) . | (h−1 kpc) . | name . | A . | B . | C . |
109 | 1010 | 1 | 30 | NW | 0.578 | −0.685 | −1.60 |
109 | 1010 | 1 | 30 | RVWa | −0.519 | 1.08 | −1.85 |
1010 | 1011 | 1 | 40 | NW | −0.209 | 2.23 | −3.32 |
1010 | 1011 | 1 | 40 | RVWa | −0.180 | 0.487 | −1.54 |
1011 | 1012 | 1 | 300 | NW | 0.227 | 0.583 | −1.03 |
1011 | 1012 | 1 | 300 | RVWa | −0.265 | 1.02 | −1.03 |
1012 | 1013 | 1 | 300 | NW | 0.373 | −0.0354 | −0.384 |
1012 | 1013 | 1 | 300 | RVWa | 0.498 | −0.241 | −0.255 |
Mhalo, min . | Mhalo, max . | rmin . | rmax . | Run . | Fit coefficients . | ||
---|---|---|---|---|---|---|---|
(M⊙) . | (M⊙) . | (h−1 kpc) . | (h−1 kpc) . | name . | A . | B . | C . |
109 | 1010 | 1 | 30 | NW | 0.578 | −0.685 | −1.60 |
109 | 1010 | 1 | 30 | RVWa | −0.519 | 1.08 | −1.85 |
1010 | 1011 | 1 | 40 | NW | −0.209 | 2.23 | −3.32 |
1010 | 1011 | 1 | 40 | RVWa | −0.180 | 0.487 | −1.54 |
1011 | 1012 | 1 | 300 | NW | 0.227 | 0.583 | −1.03 |
1011 | 1012 | 1 | 300 | RVWa | −0.265 | 1.02 | −1.03 |
1012 | 1013 | 1 | 300 | NW | 0.373 | −0.0354 | −0.384 |
1012 | 1013 | 1 | 300 | RVWa | 0.498 | −0.241 | −0.255 |
Mhalo, min . | Mhalo, max . | rmin . | rmax . | Run . | Fit coefficients . | ||
---|---|---|---|---|---|---|---|
(M⊙) . | (M⊙) . | (h−1 kpc) . | (h−1 kpc) . | name . | A . | B . | C . |
109 | 1010 | 1 | 30 | NW | 0.578 | −0.685 | −1.60 |
109 | 1010 | 1 | 30 | RVWa | −0.519 | 1.08 | −1.85 |
1010 | 1011 | 1 | 40 | NW | −0.209 | 2.23 | −3.32 |
1010 | 1011 | 1 | 40 | RVWa | −0.180 | 0.487 | −1.54 |
1011 | 1012 | 1 | 300 | NW | 0.227 | 0.583 | −1.03 |
1011 | 1012 | 1 | 300 | RVWa | −0.265 | 1.02 | −1.03 |
1012 | 1013 | 1 | 300 | NW | 0.373 | −0.0354 | −0.384 |
1012 | 1013 | 1 | 300 | RVWa | 0.498 | −0.241 | −0.255 |
Mhalo, min . | Mhalo, max . | Run . | r1, min . | r1, max . | Fit coeffs-1 . | r2, min . | r2, max . | Fit coeffs-2 . | ||
---|---|---|---|---|---|---|---|---|---|---|
(M⊙) . | (M⊙) . | name . | (h−1 kpc) . | (h−1 kpc) . | A1 . | B1 . | (h−1 kpc) . | (h−1 kpc) . | A2 . | B2 . |
109 | 1010 | NW | 1 | 6 | −3.92 | 2.18 | 6 | 20 | −5.27 | 3.69 |
109 | 1010 | RVWa | 1 | 6 | −2.68 | 0.613 | 6 | 20 | −4.54 | 3.13 |
1010 | 1011 | NW | 1 | 6 | −6.58 | 5.18 | 6 | 20 | −5.95 | 3.90 |
1010 | 1011 | RVWa | 1 | 6 | −3.11 | 1.01 | 6 | 20 | −5.43 | 3.67 |
1011 | 1012 | NW | 2 | 7 | −12.19 | 9.52 | 7 | 20 | −4.87 | 1.96 |
1011 | 1012 | RVWa | 2 | 7 | −10.63 | 8.05 | 7 | 20 | −5.94 | 2.90 |
1012 | 1013 | NW | 4 | 10 | −16.50 | 13.29 | 10 | 30 | −3.78 | 0.624 |
1012 | 1013 | RVWa | 4 | 10 | −14.13 | 9.66 | 10 | 30 | −7.17 | 3.10 |
Mhalo, min . | Mhalo, max . | Run . | r1, min . | r1, max . | Fit coeffs-1 . | r2, min . | r2, max . | Fit coeffs-2 . | ||
---|---|---|---|---|---|---|---|---|---|---|
(M⊙) . | (M⊙) . | name . | (h−1 kpc) . | (h−1 kpc) . | A1 . | B1 . | (h−1 kpc) . | (h−1 kpc) . | A2 . | B2 . |
109 | 1010 | NW | 1 | 6 | −3.92 | 2.18 | 6 | 20 | −5.27 | 3.69 |
109 | 1010 | RVWa | 1 | 6 | −2.68 | 0.613 | 6 | 20 | −4.54 | 3.13 |
1010 | 1011 | NW | 1 | 6 | −6.58 | 5.18 | 6 | 20 | −5.95 | 3.90 |
1010 | 1011 | RVWa | 1 | 6 | −3.11 | 1.01 | 6 | 20 | −5.43 | 3.67 |
1011 | 1012 | NW | 2 | 7 | −12.19 | 9.52 | 7 | 20 | −4.87 | 1.96 |
1011 | 1012 | RVWa | 2 | 7 | −10.63 | 8.05 | 7 | 20 | −5.94 | 2.90 |
1012 | 1013 | NW | 4 | 10 | −16.50 | 13.29 | 10 | 30 | −3.78 | 0.624 |
1012 | 1013 | RVWa | 4 | 10 | −14.13 | 9.66 | 10 | 30 | −7.17 | 3.10 |
Mhalo, min . | Mhalo, max . | Run . | r1, min . | r1, max . | Fit coeffs-1 . | r2, min . | r2, max . | Fit coeffs-2 . | ||
---|---|---|---|---|---|---|---|---|---|---|
(M⊙) . | (M⊙) . | name . | (h−1 kpc) . | (h−1 kpc) . | A1 . | B1 . | (h−1 kpc) . | (h−1 kpc) . | A2 . | B2 . |
109 | 1010 | NW | 1 | 6 | −3.92 | 2.18 | 6 | 20 | −5.27 | 3.69 |
109 | 1010 | RVWa | 1 | 6 | −2.68 | 0.613 | 6 | 20 | −4.54 | 3.13 |
1010 | 1011 | NW | 1 | 6 | −6.58 | 5.18 | 6 | 20 | −5.95 | 3.90 |
1010 | 1011 | RVWa | 1 | 6 | −3.11 | 1.01 | 6 | 20 | −5.43 | 3.67 |
1011 | 1012 | NW | 2 | 7 | −12.19 | 9.52 | 7 | 20 | −4.87 | 1.96 |
1011 | 1012 | RVWa | 2 | 7 | −10.63 | 8.05 | 7 | 20 | −5.94 | 2.90 |
1012 | 1013 | NW | 4 | 10 | −16.50 | 13.29 | 10 | 30 | −3.78 | 0.624 |
1012 | 1013 | RVWa | 4 | 10 | −14.13 | 9.66 | 10 | 30 | −7.17 | 3.10 |
Mhalo, min . | Mhalo, max . | Run . | r1, min . | r1, max . | Fit coeffs-1 . | r2, min . | r2, max . | Fit coeffs-2 . | ||
---|---|---|---|---|---|---|---|---|---|---|
(M⊙) . | (M⊙) . | name . | (h−1 kpc) . | (h−1 kpc) . | A1 . | B1 . | (h−1 kpc) . | (h−1 kpc) . | A2 . | B2 . |
109 | 1010 | NW | 1 | 6 | −3.92 | 2.18 | 6 | 20 | −5.27 | 3.69 |
109 | 1010 | RVWa | 1 | 6 | −2.68 | 0.613 | 6 | 20 | −4.54 | 3.13 |
1010 | 1011 | NW | 1 | 6 | −6.58 | 5.18 | 6 | 20 | −5.95 | 3.90 |
1010 | 1011 | RVWa | 1 | 6 | −3.11 | 1.01 | 6 | 20 | −5.43 | 3.67 |
1011 | 1012 | NW | 2 | 7 | −12.19 | 9.52 | 7 | 20 | −4.87 | 1.96 |
1011 | 1012 | RVWa | 2 | 7 | −10.63 | 8.05 | 7 | 20 | −5.94 | 2.90 |
1012 | 1013 | NW | 4 | 10 | −16.50 | 13.29 | 10 | 30 | −3.78 | 0.624 |
1012 | 1013 | RVWa | 4 | 10 | −14.13 | 9.66 | 10 | 30 | −7.17 | 3.10 |
Tmin . | Tmax . | Run . | δ1, min . | δ1, max . | Fit coeffs-1 . | δ2, min . | δ2, max . | Fit coeffs-2 . | ||
---|---|---|---|---|---|---|---|---|---|---|
(K) . | (K) . | name . | . | . | A1 . | B1 . | . | . | A2 . | B2 . |
104 | 105 | NW | 103 | 105 | −7.54 | 1.43 | – | – | – | – |
104 | 105 | RVWa | 103 | 105 | −5.59 | 0.937 | – | – | – | – |
105 | 107 | NW | 102 | 104 | −3.81 | 0.656 | 104 | 107 | −3.55 | 0.643 |
105 | 107 | RVWa | 102 | 104 | −3.22 | 0.546 | 104 | 107 | −2.70 | 0.456 |
102 | 108 | NW | 103 | 105 | −7.62 | 1.54 | 105 | 108 | −1.20 | 0.247 |
102 | 108 | RVWa | 103 | 105 | −5.98 | 1.19 | 105 | 108 | −1.74 | 0.300 |
Tmin . | Tmax . | Run . | δ1, min . | δ1, max . | Fit coeffs-1 . | δ2, min . | δ2, max . | Fit coeffs-2 . | ||
---|---|---|---|---|---|---|---|---|---|---|
(K) . | (K) . | name . | . | . | A1 . | B1 . | . | . | A2 . | B2 . |
104 | 105 | NW | 103 | 105 | −7.54 | 1.43 | – | – | – | – |
104 | 105 | RVWa | 103 | 105 | −5.59 | 0.937 | – | – | – | – |
105 | 107 | NW | 102 | 104 | −3.81 | 0.656 | 104 | 107 | −3.55 | 0.643 |
105 | 107 | RVWa | 102 | 104 | −3.22 | 0.546 | 104 | 107 | −2.70 | 0.456 |
102 | 108 | NW | 103 | 105 | −7.62 | 1.54 | 105 | 108 | −1.20 | 0.247 |
102 | 108 | RVWa | 103 | 105 | −5.98 | 1.19 | 105 | 108 | −1.74 | 0.300 |
Tmin . | Tmax . | Run . | δ1, min . | δ1, max . | Fit coeffs-1 . | δ2, min . | δ2, max . | Fit coeffs-2 . | ||
---|---|---|---|---|---|---|---|---|---|---|
(K) . | (K) . | name . | . | . | A1 . | B1 . | . | . | A2 . | B2 . |
104 | 105 | NW | 103 | 105 | −7.54 | 1.43 | – | – | – | – |
104 | 105 | RVWa | 103 | 105 | −5.59 | 0.937 | – | – | – | – |
105 | 107 | NW | 102 | 104 | −3.81 | 0.656 | 104 | 107 | −3.55 | 0.643 |
105 | 107 | RVWa | 102 | 104 | −3.22 | 0.546 | 104 | 107 | −2.70 | 0.456 |
102 | 108 | NW | 103 | 105 | −7.62 | 1.54 | 105 | 108 | −1.20 | 0.247 |
102 | 108 | RVWa | 103 | 105 | −5.98 | 1.19 | 105 | 108 | −1.74 | 0.300 |
Tmin . | Tmax . | Run . | δ1, min . | δ1, max . | Fit coeffs-1 . | δ2, min . | δ2, max . | Fit coeffs-2 . | ||
---|---|---|---|---|---|---|---|---|---|---|
(K) . | (K) . | name . | . | . | A1 . | B1 . | . | . | A2 . | B2 . |
104 | 105 | NW | 103 | 105 | −7.54 | 1.43 | – | – | – | – |
104 | 105 | RVWa | 103 | 105 | −5.59 | 0.937 | – | – | – | – |
105 | 107 | NW | 102 | 104 | −3.81 | 0.656 | 104 | 107 | −3.55 | 0.643 |
105 | 107 | RVWa | 102 | 104 | −3.22 | 0.546 | 104 | 107 | −2.70 | 0.456 |
102 | 108 | NW | 103 | 105 | −7.62 | 1.54 | 105 | 108 | −1.20 | 0.247 |
102 | 108 | RVWa | 103 | 105 | −5.98 | 1.19 | 105 | 108 | −1.74 | 0.300 |