The formation of hot gaseous haloes around galaxies

Abstract

1 INTRODUCTION

One of the major goals of modern galaxy formation theory is to understand the physical mechanisms that halt the star formation process, by removing, heating or preventing the infall of cold gas on to the galactic disc. X-ray observations suggest that for haloes hosting massive galaxies the majority of baryonic matter resides not only in the galaxies but also in the halo in the form of virialized hot gas (e.g. Lin, Mohr & Stanford 2003; Crain et al. 2010; Anderson & Bregman 2011). This work investigates the formation of the hot gaseous corona (also referred to as ‘hot halo’ or ‘hot atmosphere’) around galaxies, that may help reduce the rate of infall of gas on to galaxies, and has been suggested to explain the observed galaxy bimodality (Dekel & Birnboim 2006).

The hot gaseous corona is produced as a result of an important heating process, which was initially discussed by Rees & Ostriker (1977), Silk (1977), Binney (1977) and White & Rees (1978), and later in the context of the cold dark matter paradigm by White & Frenk (1991), in an attempt to explain the reduced efficiency of star formation within massive haloes. They proposed that while a dark matter halo relaxes to virial equilibrium, gas falling into it experiences a shock, and determined the cooling time of gas behind the shock. As long as the cooling time is shorter than the dynamical time, the infalling gas cools (inside the current ‘cooling radius’) and settles on to the galaxy. If, on the other hand, the cooling time exceeds the dynamical time, the gas is not able to radiate away the thermal energy that supports it. Therefore, it adjusts its density and temperature quasi-statically, forming a hot hydrostatic halo atmosphere, pressure supported against gravitational collapse. Over the past decade, the works of Birnboim & Dekel (2003) and Dekel & Birnboim (2006, hereafter DB06) investigated the stability of accretion shocks around galaxies, and concluded that a hot atmosphere forms when the compression time of shocked gas is larger than its cooling time, occurring when haloes reach a mass of about 10^11.7 M_⊙.

Numerical simulations have shown, however, that cold gas accreting through filaments does not necessarily experience a shock when crossing the virial radius, even if the spherically averaged cooling radius is smaller than the virial radius. Many groups have concluded that there are two modes of gas accretion, named as hot and cold accretion, that are able to coexist in high-mass haloes at high redshift (e.g. Kereš et al. 2005; Dekel & Birnboim 2006; Ocvirk, Pichon & Teyssier 2008; Dekel et al. 2009; Faucher-Giguère, Kereš & Ma 2011; van de Voort et al. 2011; van de Voort & Schaye 2012; Nelson et al. 2013). The hot mode of accretion refers to the accreted gas that shock-heats to the halo virial temperature. The cold mode refers to gas that flows along dark matter filaments and is accreted on to the central galaxy without being shock-heated near the virial radius. It has been found that the cold streams end up being the dominant mode of accretion on to galaxies at high redshift (e.g. Dekel et al. 2009). However, it has also been found that most of the cold gas from filaments does experience significant heating when accreted by the galaxy at radii much smaller than the virial radius (Nelson et al. 2013).

Besides the rate of gas accretion, the hot halo can be influenced by feedback mechanisms and photoionization from local sources. Feedback mechanisms can suppress cooling from the hot halo, modify the distribution of hot gas in the halo (van de Voort & Schaye 2012) and (to a limited degree) reduce the accretion rates on to haloes (Faucher-Giguère et al. 2011; van de Voort et al. 2011; Nelson et al. 2015). In this work, we investigate the impact of feedback mechanisms on the hot halo in detail and analyse whether reasonable changes to the feedback implementation result in a change to the mass scale of hot halo formation. Increasing photoionizing flux (higher star formation rate or an active nucleus) from local sources can decrease the net cooling rate of gas in the proximity of the galaxy, potentially suppressing cold gas accretion in low-mass haloes (<5 × 10¹¹ M_⊙) and decrease the mass scale for hot halo formation (Cantalupo 2010). However, the results are sensitive to the assumed escape fraction, and Vogelsberger et al. (2012) found only small effects when including local active galactic nucleus (AGN). For simplicity, we will assume the gas is only exposed to the metagalactic background radiation.

We use the suite of cosmological hydrodynamical simulations from the EAGLE project (Crain et al. 2015; Schaye et al. 2015) to investigate the physical properties of the hot gas in the halo, and their dependence on energy sources like stellar feedback and AGN feedback. Our main goal is to study the thermal properties of gas accreting on to haloes and the gas mass that remains hot in the halo (M_hot). In addition, we develop analytic and semi-analytic approaches to calculate the heating rates of gas in the halo and the mass scale of hot halo formation, which we apply in a companion work (Correa et al. in preparation, hereafter Paper II). In Paper II, we derive a physically motivated model for gas accretion on to galaxies that accounts for the hot/cold modes of accretion on to haloes and for the rate of gas cooling from the hot halo. With this model, we aim to provide some insight into the physical mechanisms that drive the gas inflow rates on to galaxies.

The outline of this paper is as follows. We describe the EAGLE simulations series used in this study and the analysis methodology in Section 2. We present our main results concerning the physical properties of hot and cold gas in the halo in Section 3 and on the modes of gas accretion in Section 4. In Section 5, we develop an analytic approach to calculate a ‘critical mass scale’, M_crit, for hot halo formation, and compare it with our numerical results and previous works. Finally, in Section 6 we summarize our conclusions.

2 SIMULATIONS

To investigate the formation and evolution of hot haloes surrounding galaxies, we use cosmological, hydrodynamical simulations from the Evolution and Assembly of GaLaxies and their Environments (EAGLE) project (Crain et al. 2015; Schaye et al. 2015). The EAGLE simulations were run using a modified version of gadget 3 (Springel 2005), a N-Body Tree-PM smoothed particle hydrodynamic (SPH) code. The EAGLE version contains a new formulation of SPH, new time stepping and new subgrid physics. Below we present a summary of the EAGLE models. For a more complete description, see Schaye et al. (2015).

The EAGLE simulations assume a Λ cold dark matter (ΛCDM) cosmology with the parameters derived from Planck-1 data (Planck Collaboration et al. 2014), Ω_m = 1 − Ω_Λ = 0.307, Ω_b = 0.04825, h = 0.6777, σ₈ = 0.8288, n_s = 0.9611. The primordial mass fractions of hydrogen and helium are X = 0.752 and Y = 0.248, respectively.

Table 1 lists the box sizes and resolutions of the simulations used in this work. We use the notation LxxxNyyyy, where xxx indicates box size (ranging from 25 to 100 comoving  Mpc) and yyyy indicates the cube root of the number of particles per species (ranging from 376³ to 1504³, with the number of baryonic particles initially equal to the number of dark mater particles). The gravitational softening was kept fixed in comoving units down to z = 2.8 and in proper units thereafter. We will refer to simulations with the mass and spatial resolution of L025N0376 as intermediate-resolution runs and to simulations with the resolution of L025N0752 as high-resolution runs.

2.1 Baryonic physics

Radiative cooling and photoheating are included as in Wiersma, Schaye & Smith (2009). The element-by-element radiative rates are computed in the presence of the cosmic microwave background (CMB), and the Haardt & Madau (2001) model for UV and X-ray background radiation from quasars and galaxies.

Black hole seeds (of mass ≈1.4 × 10⁵ M_⊙) are included in the gas particle with the highest density in haloes of mass greater than ≈1.4 × 10¹⁰ M_⊙ that do not contain black holes (Springel, Di Matteo & Hernquist 2005). Black holes can grow through mergers and gas accretion. The accretion events follow a modified Bondi–Hoyle formula that accounts for the angular momentum of the accreting gas (Rosas-Guevara et al. 2015; Schaye et al. 2015), and a free parameter that is related to the disc viscosity. AGN feedback follows the accretion of mass on to the black hole, where a fraction (0.015) of the accreted rest mass energy is released as thermal energy into the surrounding gas, and is implemented stochastically, as per the stellar feedback scheme, with a fixed free parameter heating temperature, ΔT_AGN, which is set to 10^8.5 K in the reference simulations.

When the resolution is increased, the parameters may need to be (re-)calibrated to retain the agreement with observations. The high-resolution simulation with recalibrated parameters is called Recal. In addition to Ref and Recal, we also use simulations with different feedback implementations to test the impact of feedback on the formation of the hot halo. Table 2 lists the values of the feedback parameters adopted in each simulation. In the table, the simulation identifier describes the differences in the feedback with respect to Ref. In the stellar feedback case, ‘Less/More Energetic FB’ means that in these simulations, the energy injected per mass of stars formed is lower/higher with respect to Ref. In the AGN case, ‘More Explosive AGN FB’ means that AGN feedback is more explosive and intermittent, but the energy injected per unit mass accreted by the BH does not change with respect to Ref. Additional information regarding the performance of the EAGLE simulations including an analysis of subgrid parameter variations, a study of the evolution of galaxy masses, star formation rates and sizes can be found in Crain et al. (2015), Furlong et al. (2015, 2017) and Schaye et al. (2015).

2.2 Hydrodynamics

There has been much debate regarding the systematic differences between SPH, grid codes and moving mesh grid codes when modelling fluid mixing and gas heating and cooling (e.g. Agertz et al. 2007; Vogelsberger et al. 2012; Nelson et al. 2013). It has been shown by Hutchings & Thomas (2000) and Creasey et al. (2011) that SPH simulations may not adequately resolve shocks of accreted gas. Since shocks are generally spread over several SPH kernel lengths, the heating rate is smoothed over time, potentially making it easier for radiative cooling to become important. In addition, if radiative cooling is able to limit the maximum temperature reached by the gas particle, numerical radiative losses may be enhanced.

In contrast, numerical simulations using grids do not smooth out the shocks, and are thus better at identifying shock temperatures spikes. Numerical simulations using moving mesh codes can also capture shocks accurately. However, in common with grid codes, they may suffer from numerical mixing of hot and cold gas as the fluid moves across cells. Recently, Nelson et al. (2013) compared the moving mesh code arepo (Springel 2010) with the standard SPH version of gadget, and calculated the rates of cold mode of gas accretion on to haloes and galaxies. They found that while the rates of gas accreted cold on to haloes are in very good agreement between the simulations run with gadget and arepo, the rates of gas accreted cold on to galaxies differ significantly, with galaxies in arepo having a 20 per cent lower cold fraction in 10¹¹ M_⊙ haloes. Nelson et al. (2013) concluded that most of the cold gas from filaments experiences significant heating after crossing the virial radius, implying that the numerical deficiencies inherent in different simulation codes may modify the relative contributions of hot and cold modes of accretion on to galaxies.

Some differences in the contributions of hot and cold modes of accretion on to galaxies and haloes may, however, be due to the method employed to select shock-heated gas. Previous works (e.g. Kereš et al. 2005, 2009; Faucher-Giguère et al. 2011; van de Voort et al. 2011; Nelson et al. 2013, amongst others) followed the thermal history of the gas and applied a fixed temperature cut to the distribution of the maximum past temperature (T_max), to separate hot from cold mode accretion. However, T_max is not suitable for identifying cold flows if the gas experiences a shock but cools immediately afterwards, as may happen for accretion on to galaxies. In this case, a filament that is mostly cold except at a point near the galaxy would be labelled as hot mode accretion by numerical studies using T_max, but observers would identify it as a cold flow. This practical problem may not be important for SPH simulations that suffer from ‘in shock cooling’ because they do not resolve the accretion shocks on to the galaxy, or, as in the case of van de Voort et al. (2011) and EAGLE, that impose a temperature–density relation on to high-density gas, but it may affect the conclusions inferred from moving mesh codes using the T_max statistic. To avoid this issue, we use an alternative method to identify shock-heated particles in Section 4, based on post-shock temperature values.

Hydrodynamic solvers may also produce differences in the hot/cold modes of accretion. The EAGLE version of gadget uses the hydrodynamic solver ‘Anarchy’, which greatly improves the performance on standard hydrodynamical tests, when compared to the original SPH implementation in gadget (Schaller et al. 2015, see Hu et al. 2014 for similar results). Anarchy makes use of the pressure-entropy formulation derived in Hopkins (2013), alleviating spurious jumps at contact discontinuities. It also uses an artificial viscosity switch advocated by Cullen & Dehnen (2010) that allows the viscosity limiter to be stronger when shocks and shear flows are present. In addition, Anarchy includes an artificial conduction switch (similar to that of Price 2008), the C² Wendland (1995) kernel and the time step limiters of Durier & Dalla Vecchia (2012). These changes ensure that ambient particles do not remain inactive when a shock is approaching.

Recently, Sembolini et al. (2016) compared cosmological simulations of clusters using SPH as well as mesh-based codes. They found that the modern SPH schemes (such as Anarchy) that allow entropy mixing produce gas entropy profiles that are indistinguishable from those obtained with grid-based schemes. In addition, Schaller et al. (2015) compared the EAGLE simulations with simulations run with the same subgrid physics, but using the standard gadget rather than the Anarchy hydrodynamic solver. They found that while simulations with standard SPH contain haloes with a large number of dense clumps of gas at all radii, Anarchy's ability to mix phases allows dense clumps to dissolve into the hot halo. These substantial improvements of the SPH formulation in the EAGLE simulations motivate a detailed description of the resulting predictions for hot halo formation and of hot/cold mode accretion.

2.3 Identifying haloes and galaxies

Throughout this work, we select the largest subhalo in each Friends-of-Friends (FoF) group, and use the subfind algorithm (Springel, White & Hernquist 2001; Dolag et al. 2009) to identify the substructures (subhaloes) within it. The FoF algorithm adopts a dimensionless linking length of 0.2, and the subfind algorithm calculates halo virial masses and radii via a spherical overdensity routine that centres the main subhalo from the FoF group on the minimum of the gravitational potential. We define halo masses, M₂₀₀, as the mass of all matter within the radius, R₂₀₀, for which the mean internal density is 200 times the critical density of the Universe.

To select the gas associated with the central galaxies embedded in each resolved halo, we identify the gravitationally bound cold and dense gas within R₂₀₀ that is star forming and/or has a hydrogen number density, n_H > 0.01cm^{− 3}, and temperature T < 10⁵K. We also require all particles to be contained within a sphere of radius 0.15 × R₂₀₀ in order to avoid labelling infalling cold flows (that would be included by the T − n_H cuts but are mostly at large radii) as part of the galaxy.

2.4 Measuring gas accretion

Once we have identified the haloes, we build merger trees across the simulation snapshots.¹ The standard procedure to build a halo merger tree is to link each progenitor halo with a unique ‘descendant’ halo in the subsequent output (see e.g. Fakhouri, Ma & Boylan-Kolchin 2010). To do so, we identify the main branches of the halo merger trees and compute the halo (and central galaxy) accretion rate between two consecutive snapshots. At each output redshift (snapshot), we select the most massive haloes within each FoF group and consider them to be ‘resolved’ if they contain more than 1000 dark matter particles, which corresponds to a minimum halo mass of M₂₀₀ = 10^9.8 M_⊙ (10^8.8 M_⊙) in the intermediate- (high-) resolution simulations. This limit on the number of dark matter particles results from a convergence analysis that we present in Paper II, where we find that in smaller haloes the accretion on to galaxies does not converge, indicating that the inner galaxies are not well resolved. We refer to the resolved haloes as ‘descendants’, and then link each descendant with a unique ‘progenitor’ at the previous output redshift. This is non-trivial due to halo fragmentation, in which subhaloes of a progenitor halo may have descendants that reside in more than one halo. Such fragmentation can be either spurious or due to a physical unbinding event. To account for this, we link the descendant to the progenitor that contains the majority of the descendant's 25 most bound dark matter particles (see Correa et al. 2015b for an analysis of halo mass history convergence using these criteria to connect haloes between snapshots).

We distinguish between gas accreted on to a halo and gas accreted on to a galaxy. For each descendant halo at z_i and its linked progenitor at z_j (z_j > z_i), we identify the particles that are in the descendant but not in its progenitor by performing particle ID matching. We then select particles that are new in the halo and reside within the virial radius, as particles accreted on to the halo in the redshift range of z_i ≤ z < z_j. The accretion rate on to galaxies is further explored in Paper II, where we follow the methodology described above for calculating accretion rates on to haloes, and we select the new particles within the radius 0.15 × R₂₀₀ as particles accreted on to the galaxies in the redshift range of z_i ≤ z < z_j (see Paper II, Section 2.1, for a discussion on methods for calculating gas accretion on to galaxies).

3 HOT HALO FORMATION

The simple models of galaxy formation (e.g. Rees & Ostriker 1977; White & Rees 1978) assume that as long as the cooling time, t_cool, is shorter than the dynamical time, t_dyn, the infalling gas cools (inside a ‘cooling radius’, White & Frenk 1991) and settles into the galaxy. Otherwise, the gas is unable to efficiently radiate its thermal energy and forms a hot hydrostatic atmosphere, which is pressure supported against gravitational collapse. More recent semi-analytic models of galaxy formation assume that the cooling radius expands outwards as a function of time, therefore the comparison is done between gas cooling time and a different time-scale representing the time available for cooling, like the time since the last major event or the time of halo formation (see e.g Lacey et al. 2016). In this section, we investigate when the hot hydrostatic halo forms in the EAGLE simulations, by analysing the interplay between the cooling and dynamical times of the gas particles in the halo. Throughout this work, we define hot halo gas as all gas particles that have t_cool > t_dyn and that do not form part of the galaxy, i.e. r > 0.15R₂₀₀.

Fig. 1 shows temperature profiles (left column), the logarithm of the ratio between cooling times and dynamical times (middle column) and the respective mass-weighted probability density function (PDF) of log₁₀t_cool/t_dyn (right column) for gas from haloes in the mass range of 10^11.9–10^12.1 M_⊙ (top row), 10^11.4–10^11.6 M_⊙ (middle row) and 10^10.9–10^11.1 M_⊙ (bottom row) at z = 0, taken from the Ref-L025N0752 simulation. In the right-hand panels and throughout this work, the PDFs are calculated by stacking haloes in the selected mass range and distributing the gas particles in logarithmic bins of size 0.1 dex. We then sum the mass of the gas particles in each bin and normalize the distribution by the total gas mass. In the left-hand panels, the legend lists the median values of the mass, virial temperature (defined as |$T_{\rm {vir}}=\frac{\mu m_{\rm {p}}}{2k_{\rm {B}}}V^{2}_{200}$|⁠, with |$V^{2}_{200}=GM_{200}/R_{200}$|⁠) and virial radius of haloes selected in each mass bin. The left panels also show in solid, dotted and dashed lines the median temperature per radial bin of gas from haloes taken from the simulations Ref-L025N0376, Ref-L025N0752 and Recal-L025N0752, respectively.

Figure 1.

Open in new tabDownload slide

Temperature profile (left column), logarithm of the ratio between cooling times and local dynamical times (middle column) and the mass-weighted probability density function (PDF, right column) of log₁₀t_cool/t_dyn of gas from haloes in the mass range of 10^11.9–10^12.1 M_⊙ (top row), 10^11.4–10^11.6 M_⊙ (middle row) and 10^10.9–10^11.1 M_⊙ (bottom row) at z = 0 taken from the Ref-L025N0752 simulation. The number of particles in a pixel is used for colour coding. The solid, dotted and dashed lines in the left-hand panels correspond to the median temperature per radial bin for different simulations.

The left-hand panels of Fig. 1 show that while there is very good agreement in the median temperature profiles at small (r/R₂₀₀ < 0.2) and large (r/R₂₀₀ > 0.4) radii, at intermediate radii the median temperatures from the intermediate-resolution run (Ref-L025N0376) are larger by up to 0.3 dex than those from the high-resolution runs (L025N0752). This is in agreement with the convergence analysis of Nelson et al. (2016), who concluded that the physics (different models of stellar winds or AGN feedback) has a greater impact on T(r) than resolution. We also find that in the radial range of r = [0.2 − 0.4]R₂₀₀, where the convergence with resolution is poorest, the median temperatures drop from T_gas ∼ T_vir to T_gas ∼ 10⁴ K (in agreement with Nelson et al. 2016, r_drop ≈ 0.25R₂₀₀ and van de Voort & Schaye 2012, r_drop ≈ 0.2 − 0.4R₂₀₀), because of the high densities that rapidly decrease the gas cooling times, enabling it to radiate away its thermal energy and join the ISM.

The top and middle left-hand panels of Fig. 1 show that there is relatively little gas with T ∼ 10⁵ K at small and intermediate radii, reflecting the short cooling times at these temperatures. The cooling flow in the hot halo is formed by T = 10⁶ K gas that slowly decreases in temperature as it loses hydrostatic support due to cooling. The ISM consists of T = 10⁴ K gas at r/R₂₀₀ < 0.15. For a better understanding of the hot halo forming as a function of halo mass and its effect on the infall rate of gas on to the galaxy, we next analyse the middle and right columns. The bottom middle panel shows that in haloes with masses between 10^10.9and10^11.1 M_⊙, most of the gas has low temperature (T_gas < T_vir), short cooling times (t_cool < t_dyn) and infalls towards the central galaxy, although a substantial fraction of gas has t_cool > t_dyn at r ∼ R₂₀₀. At larger halo masses, a larger fraction of the gas is unable to cool and therefore forms a hot halo. The middle panel shows that haloes between 10^11.4and10^11.6 M_⊙ are in the intermediate stage between developing a hot stable atmosphere (gas with T_gas ∼ T_vir and t_cool > t_dyn) and continuing to fuel the galaxy. The top middle panel clearly shows a stable hot halo and a reduced amount of gas infalling towards the galaxies (gas at r > 0.3R₂₀₀ and t_cool < t_dyn).

Fig. 2 is similar to Fig. 1, but shows haloes in the mass range of 10^11.9–10^12.1 M_⊙ (top panels), 10^11.4–10^11.6 M_⊙ (bottom panels) at z = 2.24. The top middle panel shows that 10¹² M_⊙ haloes develop a hot atmosphere, despite the significant fraction of cold gas at large radii that is accreted on to the halo. This cold gas forms part of the cold filamentary flows, which cross the virial radius, and are directly accreted on to the central galaxy. The cold accretion mode is best seen as the gas at T = 10⁴ K and r/R₂₀₀ > 0.4 in the bottom panel of Fig. 2, which remains cool (with T = 10⁴ K), as it is accreted on to the galaxy.

Figure 2.

Open in new tabDownload slide

Same as Fig. 1 but for haloes at z = 2.24.

The presence of cold flows produces gaseous haloes with an isothermal temperature profile of T_gas ∼ 10⁴ K at all radii. Besides analysing the cooling time profiles, Nelson et al. (2016) and van de Voort & Schaye (2012) analysed the entropy profiles of haloes at z = 2, and concluded the that while the entropy of the cold-mode gas decreases smoothly and strongly towards the centre, the entropy of the hot-mode gas decreases slightly down to 0.2R₂₀₀, after which it drops steeply. We find that the cooling time profiles of the hot (t_cool > t_dyn) and cold (t_cool < t_dyn) modes follow the entropy profiles.

Figs 1 and 2 show that as the halo mass increases, so does the hot gas mass. In the case of 10¹¹ M_⊙ haloes, the bottom left-hand panel of Fig. 1 shows that there is a large fraction of gas at r > 0.4R₂₀₀ with temperatures equal to or larger than the halo virial temperature. This seems to indicate that gas is shock-heating to the halo virial temperature when crossing R₂₀₀ and forming a hot atmosphere. However, gas with t_cool ∼ t_dyn in the outer parts of ∼10¹¹ M_⊙ haloes does not imply that the halo formed a stable hot atmosphere via gravitational accretion shocks, since the gas is affected by the extragalactic UV/X-ray background as we show in the next section.

The figures also show that as haloes are growing a hot atmosphere, the t_cool/t_dyn PDF begins to present a strong bimodal shape, with a local maximum at t_cool > t_dyn followed by a local minimum at t_cool < t_dyn (see top right-hand panel). We then conclude that the bimodality in the cooling time PDF provides a clear signature of the formation of a hot halo, and the right-hand panels indicate that the hot hydrostatic halo forms in the halo mass range of 10^11.5 M_⊙–10¹² M_⊙ with a weak dependence on redshift.

We also analyse the radial velocity distributions of gas and find that in haloes of mass 10¹² M_⊙ at z = 0, 92.2 per cent (61.3 per cent) of hot gas has an absolute radial velocity lower than 100 km s⁻¹ (50 km s⁻¹), indicating that most of it is in hydrostatic equilibrium and not accreting on to the galaxy.

3.1 The impact of photoheating in low-mass haloes

In the previous section, we analysed the dependence of the gas cooling rates on halo mass and concluded that a halo with a hot hydrostatic atmosphere should present a strongly bimodal t_cool/t_dyn PDF with a local maximum at t_cool > t_dyn. As the virial temperature decreases from 10^5.5 K to 10^5.2 K, we would naively expect the peak in the t_cool/t_dyn PDF to shift towards shorter cooling times. This is, however, not the case: we find that at z = 0, the gas in haloes less massive than 10¹¹ M_⊙ (T_vir ≈ 10^5.2 K) is affected by the extragalactic UV/X-ray background radiation, which strongly suppresses the net cooling rate of gas in the temperature range of T ∼ 10⁴–10⁵ K (Efstathiou 1992; Wiersma et al. 2009). As a result, the peak in the t_cool/t_dyn PDF remains at t_cool ∼ t_dyn. This can be seen in Fig. 3, where we show the t_cool/t_dyn PDF of gas from haloes in the mass range of 10^10.4–10^10.6 M_⊙ (olive lines), 10^10.8–10^11.0 M_⊙ (orange lines) and 10^11.2–10^11.4 M_⊙ (red lines). The solid lines correspond to the case where the cooling rates are calculated for gas exposed to the evolving UV/X-ray background from Haardt & Madau (2001), while the dashed lines correspond to the case where the cooling rates are calculated for gas in collisional ionization equilibrium (CIE) and not exposed to the background radiation. Note, however, that we apply these CIE cooling rates to simulations that were run using cooling rates that did account for photoheating, limiting the gas temperature to ∼10⁴ K.

Figure 3.

Open in new tabDownload slide

Mass-weighted PDF of the logarithm of the ratio between the net radiative cooling time and the local dynamical time for gas from haloes in the mass range of 10^10.4–10^10.6 M_⊙ (green lines), 10^10.8–10¹¹ M_⊙ (orange lines), 10^11.2–10^11.4 M_⊙ (red lines) at z = 0. The solid lines correspond to the case where the cooling rates are calculated for gas exposed to the evolving UV/X-ray background from Haardt & Madau (2001), while the dashed lines correspond to the case where the cooling rates are calculated for gas in collisional ionization equilibrium (CIE).

Fig. 3 shows that there is no large difference in the gas t_cool/t_dyn PDF for haloes more massive than 10^11.3 M_⊙, indicating that there is no strong impact of the background radiation on the cooling rates of gas from haloes with virial temperatures larger than 10⁵ K. In smaller haloes, the peak in the t_cool/t_dyn PDF curve is shifted to t_cool ∼ 0.3t_dyn in the case of no background radiation.

3.2 The impact of feedback

Feedback can affect the formation of the hot hydrostatic halo around galaxies. For example, very energetic SN activity generates large outflows and strong winds that shock against the gaseous halo. As a result, the winds can fill the halo with gas expelled from the galaxy, increasing the amount of hot gas at large radii. In this subsection, we compare the t_cool/t_dyn mass-weighted PDF at fixed halo mass obtained from simulations with different feedback implementations (see Table 2 for reference and Section 2.1 for a brief description), and determine, by analysing whether the cooling time PDF is bimodal, the mass range where the hot halo is forming.

Fig. 4 shows the mass-weighted PDF of t_cool/t_dyn of gas from haloes in the Ref-L100N1504 simulation in the mass range of 10^11.4–10^11.6 M_⊙, 10^11.9–10^12.1 M_⊙, 10^12.4–10^12.6 M_⊙ and 10^13.4–10^13.6 M_⊙ at z = 0 (left-hand panel) and z = 2.24 (right-hand panel). In the figures, the labels show the median halo mass for each mass bin. The region where t_cool > t_dyn corresponds to hot gas in the halo.

Figure 4.

Open in new tabDownload slide

Mass-weighted PDF of the logarithm of the ratio between cooling times and local dynamical times of gas from haloes in the mass range of 10^11.4–10^11.6 M_⊙, 10^11.9–10^12.1 M_⊙, 10^12.4–10^12.6 M_⊙ and 10^13.4–10^13.6 M_⊙ at z = 0 (left-hand panel) and z = 2.24 (right-hand panel).

As the halo mass increases so does the amount of hot gas (Section 3.3). The t_cool/t_dyn PDF of gas in 10¹² M_⊙ haloes at z = 0 shows a bimodal shape that becomes mostly unimodal in higher mass haloes. At z = 2.24, the bimodality persists up to the highest mass bin (10¹³ M_⊙) due to the contribution from cold flows that populate the peak at t_cool < t_dyn. We find that the presence of the bimodality in the t_cool/t_dyn PDF indicates the increasing amount of hot gas at large radii and the eventual formation of the hot halo. Then, from visual inspection, we determine that the hot hydrostatic atmosphere is forming in haloes with masses between 10^11.5 and 10¹² M_⊙ at z = 0 and z = 2.24.

The panels in Fig. 5 repeat the analysis shown in the left-hand panel of Fig. 4, but instead show t_cool/t_dyn mass-weighted PDFs for the L025N0376 simulations with different feedback prescriptions. In this case, the PDFs correspond to gas from haloes in the mass range of 10^11.4–10^11.6 M_⊙, 10^11.6–10^11.8 M_⊙, 10^11.9–10^12.1 M_⊙ and 10^12.4–10^12.6 M_⊙. In the panels, the top left legends indicate the total gas mass in the halo (M_gas).

Figure 5.

Open in new tabDownload slide

Mass-weighted PDF of the logarithm of the ratio between cooling times and local dynamical times of gas from haloes in the mass range of 10^11.4–10^11.6 M_⊙, 10^11.6–10^11.8 M_⊙, 10^11.9–10^12.1 M_⊙ and 10^12.4–10^12.6 M_⊙ at z = 0. The legends show the median mass of the haloes in each mass bin. The different panels show L025N0376 simulations with different feedback prescriptions: no AGN (top left), weak stellar feedback (top right), strong AGN (bottom left) and strong stellar feedback (bottom right).

The simulation shown in the top left-hand panel of Fig. 5 does not have AGN feedback, while the one in the bottom left-hand panel uses a more explosive AGN feedback. Both include the same feedback from star formation as in Ref. It can be seen that neither the bimodality of the t_cool/t_dyn PDF nor the amount of hot gas are strongly affected by AGN feedback in haloes with masses between 10^11.5and10^11.9 M_⊙. The right-hand panels in Fig. 5 show the t_cool/t_dyn PDFs in the less (top panel) and more energetic stellar feedback (bottom) scenarios (both including the same AGN feedback as in the Ref model). For these halo masses, stellar feedback has a strong impact on the t_cool/t_dyn PDFs. While a more energetic stellar feedback increases the fraction of hot gas, at least in the halo mass range probed by these simulations (<10¹² M_⊙) and thus limits the build-up of cold-mode gas in the halo centre (in agreement with van de Voort & Schaye 2012), a less energetic stellar feedback maintains the bimodality in the t_cool/t_dyn PDF but shifts the peak in the t_cool/t_dyn PDF of hot gas towards larger cooling times. We find that the bimodality of the t_cool/t_dyn PDF is present in 10^11.7 M_⊙ haloes with more energetic stellar feedback and in 10¹² M_⊙ haloes with less energetic stellar feedback.

In the following section, we further analyse the dependence of the total hot gas mass on halo mass, redshift and feedback. In Section 5, we derive an analytic model for the formation of a stable hot hydrostatic atmosphere. In the model, we calculate a halo mass scale for which the gravitational heating rate of the hot halo gas balances the gas cooling rate, thus keeping the gas hot and enabling the formation of a hot atmosphere. With the model, we show that the ability of a halo to develop a hot hydrostatic atmosphere depends on the amount of hot gas that the halo already contains, which we calculate in the next subsection, and on the fraction of accreted gas that shock-heats to the halo virial temperature (Section 4).

3.3 Hot gas mass

In order to better understand the build-up of the hot gas mass, M_hot (mass of gas with t_cool > t_dyn and r > 0.15R₂₀₀) as haloes evolve, in this section we look for a correlation between M_hot and the total halo mass as a function of redshift. Fig. 6 shows the median ratio of M_hot/(Ω_b/Ω_m)M₂₀₀ (with Ω_b/Ω_m = 0.146 the universal baryon fraction) taken from a range of simulations (as indicated in the legends) at z = 0 (top panel) and at z = 2.24 (the second panel from the top). In these panels, the error bars show the 1σ scatter for the Ref-L025N0752 and Ref-L100N1504 simulations. The median ratio of M_hot/(Ω_b/Ω_m)M₂₀₀ is also shown in the third panel from the top, but in this case the values are taken from the Ref-L100N1504 simulation and at various output redshifts.

Figure 6.

Open in new tabDownload slide

Fraction of hot (with t_cool > t_dyn) gas mass with respect to the total halo mass, M₂₀₀ (normalized by the universal baryon fraction), as a function of M₂₀₀ for different simulations at z = 0 (top panel), at z = 2.24 (second panel from the top), and for Ref-L100N1504 at various redshifts (from z = 0 to z = 6, third panel from the top). The error bars in the top and middle panels show the 1σ scatter. Each bin contains at least five haloes. The bottom panel shows the residual of the data points with respect to the best-fitting expression (equations 7–10).

The top panel of Fig. 6 highlights the relatively poor agreement between the intermediate- and high-resolution simulations, with the latter predicting somewhat higher hot gas fractions. Good agreement is however achieved at z = 2.24 (middle panel). Although the Ref-L100N1504 simulation is not fully converged with respect to the numerical resolution at z = 0, the convergence with box size is excellent at all redshifts. The intermediate-resolution runs show that the hot gas represents <10 per cent of the total halo gas mass for M₂₀₀ < 10^11.6 M_⊙ at z = 0. The hot gas mass fraction reaches 80–90 per cent in 10^13.6 M_⊙ haloes and remains roughly constant for higher masses. In very low-mass haloes (M₂₀₀ < 10^10.5 M_⊙), the hot gas mass fraction also remains roughly constant (M_hot/(Ω_m/Ω_b)M₂₀₀ ≈ 0.02–0.03). In these haloes, cold accretion dominates; therefore, the heating mechanism that maintains M_hot is the UV background as discussed in Section 3.1.

The third from the top panel of Fig. 6 shows the evolution of the hot gas fraction. In haloes larger than 10^11.5 M_⊙, M_hot/(Ω_b/Ω_m)M₂₀₀ remains constant over the redshift range of 3–6 and at lower redshift it increases somewhat with time. In smaller haloes (M₂₀₀ < 10^11.5 M_⊙), M_hot/(Ω_b/Ω_m)M₂₀₀ increases with time until z = 1 but decreases thereafter. We calculate the cooling rate of gas exposed to the UV background, and in the absence of it and compare the hot gas mass. We find that the hot gas mass in low-mass haloes increases due to the heating produced by the background radiation. In the case of gas not being exposed to the UV background, the total hot gas mass decreases with increasing redshift at fixed halo mass. We also find that the differences between M_hot occurs in haloes lower than 10^11.3 M_⊙.

Next, we investigate how the presence of different feedback mechanisms affect the hot gas as well as the total gas mass (M_gas) in the halo (all gas contained between 0.15and1R₂₀₀). The top panel of Fig. 7 shows the M_gas − (Ω_b/Ω_m)M₂₀₀ relation for haloes in the mass range of 10¹⁰–10¹³ M_⊙ at z = 0 for the L025N0376 simulations. The different coloured lines correspond to simulations with different feedback prescriptions.

Figure 7.

Open in new tabDownload slide

Top panel: Fraction of gas mass (M_gas(0.15R₂₀₀ < r < R₂₀₀) with respect to the total halo mass, M₂₀₀ (normalized by the universal baryon fraction), as a function of M₂₀₀ at z = 0. The different lines correspond to L025N0376 simulations with different feedback prescriptions (see Table 2 and/or Section 2). Bottom panel: Fraction of hot gas (gas with t_cool > t_dyn) with respect to M_gas as a function of M₂₀₀.

This panel shows that the impact of feedback increases with halo mass and that stellar feedback has a larger impact on the amount of gas in the halo than AGN feedback. We find that doubling the efficiency of the stellar feedback increases the gas mass fraction by a factor of 1.3 in 10¹² M_⊙ haloes relative to the Ref model, whereas halving the efficiency decreases the gas mass fraction by a factor of 2.5. No (Explosive) AGN feedback results in an increase (decrease) by a factor of 1.5 in the gas mass fraction.

While efficient stellar feedback increases the gas mass in the halo, more explosive AGN feedback decreases it. Overall, it seems that in haloes more massive than 10¹² M_⊙ there is a greater difference in the gas mass fraction between Ref and More Energetic FB than between Ref and More Explosive AGN FB. This is due to two different reasons. Physically, AGN feedback mainly ejects gas mass from the halo, or prevents it from falling into the halo, whereas stellar feedback ejects gas out of the galaxy into the inner halo. Numerically, although stellar and AGN feedback use a similar thermal implementation (Dalla Vecchia & Schaye 2012), there is a difference in the actual energetics of the processes. The energy injected per mass of stars formed changes between Ref and More Energetic FB, whereas the energy injected per unit mass accreted by the BH does not change between Ref and More Explosive AGN FB. In the latter, it is only the intermittency and the explosiveness that changes as a consequence of the change in the temperature of the AGN. In the case of Less Energetic FB, we find that the gas mass in the halo decreases because more gas is accreted by the galaxy and locked up in stars.

The bottom panel of Fig. 7 shows the variation of the ratio M_hot/M_gas with feedback. It can be seen that in haloes less massive than 10^11.5 M_⊙, the M_hot/M_gas ratio increases with decreasing halo mass, indicating that most of the halo gas is heated by the X-ray/UV background (see Section 3.1 for a discussion). In haloes more massive than 10^11.5 M_⊙, M_hot/M_gas increases with halo mass. While a more energetic stellar feedback increases the hot mass fraction by 10 per cent, no AGN feedback decreases it by 8 per cent in 10¹² M_⊙ haloes. In the case of Less Energetic FB and More Explosive AGN, M_hot/M_gas increases by 10 per cent and 3 per cent (on average) with respect to Ref, respectively, in the halo mass range of 10^11.5–10¹² M_⊙.

In this section, gas particles with long cooling times (t_cool > t_dyn) are considered hot and counted in the calculation of M_hot. Different from this work, van de Voort & Schaye (2012) separated hot and cold gas by performing a T_max cut and found that the hot fraction as a function of radius decreases not only when AGN feedback is switched on but also when stellar winds are enhanced. The reason for this is the way stellar feedback is implemented. In the more energetic stellar feedback simulation used by van de Voort & Schaye (2012), the wind velocity scales with the local sound speed, so it largely overcomes the pressure of the ISM, blowing the gas out of the galaxy and halo, thus decreasing the amount of hot gas. In our work, the efficiency of stellar feedback is regulated by the fraction of the energy budget available (f_th), which makes it more/less energetic and controls the frequency of feedback events, but the temperature increase is kept fixed.

So far we have analysed the behaviour of the hot gas mass in the halo. In the next section, we investigate the fraction of gas that is accreted via hot and cold modes as a function of halo mass and redshift.

4 HOT AND COLD MODES OF ACCRETION

Over the last decade, numerical simulations have shown that gas accretion on to haloes occurs in two different modes, gas either shock-heats to the halo virial temperature near the virial radius (the hot accretion mode), or crosses the virial radius unperturbed (the cold accretion mode). Several works have found that cold accretion dominates in low-mass haloes (M₂₀₀ < 10¹² M_⊙) and that the transition mass increases weakly with increasing redshift (e.g. Katz et al. 2003; Kereš et al. 2005, 2009; Ocvirk et al. 2008; van de Voort et al. 2011; van de Voort & Schaye 2012). The two modes coexist at high redshift in massive haloes, which develop a hot hydrostatic atmosphere despite experiencing significant cold accretion through filaments, generally referred to as ‘cold flows’ (Kereš et al. 2005; Dekel et al. 2009). Cold flows are important for galaxy formation, because even if they experience significant heating when crossing the hot atmosphere (Nelson et al. 2013), they are responsible for delivering cold, star-forming, gas deep within the halo (e.g. Dekel et al. 2009). In the following sections, we investigate different definitions that can be used to calculate the modes of accretion in the EAGLE simulations, analyse the impact of the hydrodynamic scheme and obtain best-fitting relations.

4.1 Definition of hot accretion

The contributions from the two different modes of accretion have generally been calculated from the temperature history of the accreted gas (e.g. Kereš et al. 2005; van de Voort et al. 2011; Nelson et al. 2013, amongst others), by following the maximum temperature, T_max, each gas particle has ever reached. Kereš et al. (2005) found a clear bimodality in the distribution of T_max of accreted particles, and they proposed a threshold value of T_max = 2.5 × 10⁵ K, given by the minimum in the distribution of T_max values, to determine whether gas is accreted hot (T_max ≥ 2.5 × 10⁵ K) or cold (T_max < 2.5 × 10⁵ K). This is an often used method but other approaches have also been taken. For example, Brooks et al. (2009) identified hot gas accretion based on an entropy jump criterion and concluded that their method led to a distinction between hot/cold modes in very good agreement with the selection of hot/cold gas based on the use of a constant temperature threshold. In this section, we compare T_max with other variables that can also give us some insight into whether a particle experienced a shock when crossing the virial radius.

We follow the method described in Section 2.4, and look for gas particles that crossed the virial radius between two consecutive snapshots (z_i ≤ z < z_j). Then, we calculate the mass-weighted PDFs of the gas particles T_max, temperature (T_gas) and entropy (S_gas) at redshift z_i (see Fig. A1). We find that the PDFs for the redshift interval 0.0 < z ≤ 0.1 and 2.0 < z ≤ 2.2 are bimodal, but only for haloes larger than 10¹² M_⊙, with the location of the local minimum changing with M₂₀₀. A detailed analysis of the PDFs can be found in Appendix A. We next use the minimum of the PDFs from the 10¹² M_⊙ haloes as a threshold value to separate particles accreted hot and cold. For T_gas and T_max PDFs, we select the threshold to be T_min = 10^5.5 K, and for S_gas, we use S_min = 10^7.2 K cm².

Fig. 8 shows the fraction of gas particles accreted hot, f_acc,hot, in the redshift range of 0–0.1 as a function of halo mass. The different lines correspond to f_acc,hot calculated using different definitions. We first calculated f_acc,hot requiring that the gas particles at redshift 0 have temperatures higher than 10^5.5 K (olive solid line), or that the gas particles have temperatures higher than the host halo virial temperature (orange dotted line). Then, we calculated f_acc,hot requiring that the gas particle maximum past temperature is higher than 10^5.5 K (red dashed line) or higher than the host halo virial temperature (dark red dot–dashed line). Finally, we calculated f_acc,hot requiring that the entropy of the gas particles is larger than 10^7.2 K cm² (purple dot-dot-dot–dashed line) or that the gas particles cooling time is larger than the local dynamical time (blue dashed line).

Figure 8.

Open in new tabDownload slide

Fraction of gas accreted hot during the redshift range of 0 ≤ z < 0.1 as a function of halo mass. The various curves correspond to the hot fraction calculated using different methods as indicated in the legends.

As previous works have shown (e.g. van de Voort et al. 2011; Nelson et al. 2013), the hot fraction depends very much on the definition. The fixed temperature cut, T_gas ≥ 10^5.5 K (T_max ≥ 10^5.5 K), gives a hot fraction that increases from 0.2 (0.5) in 10^11.5 M_⊙ haloes to 0.7 (0.9) in 10¹² M_⊙ haloes, thus showing a smooth transition from cold to hot accretion in the halo mass range of 10¹¹ to 10¹³ M_⊙. On the other hand, the criterion T_gas ≥ T_vir (T_max ≥ T_vir) indicates that for the most massive haloes the hot mode accounts for only 20 per cent (60 per cent) of the accreted gas particles. For the low-mass haloes, T_gas ≥ T_vir (T_max ≥ T_vir) gives a hot fraction that increases from 0.1 (0.35) in 10¹¹ M_⊙ haloes to 0.35 (1.0) in 10¹⁰ M_⊙ haloes. This seems to indicate that there is another mechanism, such as shocks driven by winds or heating by the extragalactic UV/X-ray background radiation, that makes gas reach temperatures higher than the halo's T_vir. In the case of T_max, stellar feedback events certainly increase the hot fraction, since we obtain a peak in the T_max PDFs at 10^7.5 K for all halo mass, which disappears when stellar feedback is switched off.

Fig. 8 also shows that f_{acc, hot} calculated using S_gas ≥ 10^7.2 K cm² and t_cooling ≥ t_dyn decreases in the halo mass range of 10¹⁰–10¹² M_⊙ and increases for higher halo masses, being in agreement with f_acc,hot from T_gas ≥ 10^5.5 K in haloes larger than 10¹² M_⊙. From this upturn, we conclude that for haloes with T_vir ≫ 10⁵ K (M₂₀₀ > 10¹² M_⊙), a large fraction of the accreted gas goes through a virial shock, and thus we can safely separate hot and cold accretion using T_gas ≥ 10^5.5 K or T_max ≥ 10^5.5 K. In lower mass haloes (M₂₀₀ < 10¹² M_⊙), separating hot and cold accretion is not so easy. Although the hot halo is not expected to form (DB06) and the T_gas, T_max and S_gas PDFs are unimodal, UV background radiation (as discussed in Section 3.1) can significantly increase the cooling time and entropy of accreted gas.

Throughout this work, we have investigated how gas heated by either stellar or AGN feedback, UV/X-ray background radiation or accretion shocks, evolves with halo mass and redshift. We next aim to identify gas that is mostly heated by accretion shocks when crossing the virial radius. While stellar or AGN feedback does not strongly impact gas falling on to the halo (van de Voort et al. 2011, it does strongly affect gas falling on to the galaxy, van de Voort et al. 2011; Paper II), UV radiation does for low-mass haloes (<10^11.5 M_⊙). This can be seen in Fig. 8 through the unexpected increase in f_acc,hot towards very low halo masses using S_gas ≥ 10^7.2 K cm², t_cool ≥ t_dyn, T_gas ≥ T_vir and T_max ≥ T_vir. These methods clearly indicate that gas is hot after falling on to halo, but not necessarily due to shock-heating. For that reason, we decide to use a fixed temperature cut to calculate the hot mode of accretion. Note that T_max is updated whenever the gas particle reaches a higher temperature. However, if the gas particle is star-forming, T_max is not updated, because we impose a lower limit on the temperature of such gas. As in van de Voort et al. (2011), ignoring shocks in the ISM is appropriate because we are interested in the T_max the gas reached before accreting on to the galaxy.

Because the gas temperature after shock-heating can be slightly higher or lower than T_vir, we analyse the dependence of f_acc,hot on the fraction of T_vir used in the definition of hot accretion. As expected, at fixed halo mass f_acc,hot increases with decreasing fraction of T_vir, reaching f_acc,hot = 0.3 and 0.6 in 10¹² M_⊙ haloes for T_gas ≥ 0.8T_vir and T_gas ≥ 0.5T_vir, respectively. In addition, the virial shock can be located slightly inwards or outwards of R₂₀₀. We therefore calculated f_acc,hot (using T_gas ≥ 10^5.5 K) for gas crossing 0.8R₂₀₀ and 1.2R₂₀₀. We found that for haloes more massive than 10¹² M_⊙, f_acc,hot is insensitive to the precise value of the radius. In lower mass haloes, the difference between f_acc,hot for gas crossing 0.8R₂₀₀ and 1.2R₂₀₀ can be as high as 0.2dex in 10^11.5 M_⊙ haloes.

Fig. 9 shows f_acc,hot calculated using T_max ≥ 10^5.5 K (solid lines) and T_gas ≥ 10^5.5 K (dashed lines) at z = 0 (green lines), z = 2 (red lines) and z = 4 (thick blue lines). For 10¹² M_⊙ haloes, the fraction of shock-heated gas particles with T_gas ≥ 10^5.5 K is a factor of 1.25 lower than f_acc,hot given by T_max ≥ 10^5.5 K at z = 0, a factor of 1.6 at z = 2 and a factor of 2 at z = 4. Although changing the threshold value can bring the f_acc,hot curves into better agreement, some disagreement is expected because gas that was heated once may cool later.

Figure 9.

Open in new tabDownload slide

Fraction of gas accreted hot during the redshift ranges 0 ≤ z < 0.1 (green lines), 2.0 ≤ z < 2.2 (red lines) and 4.0 ≤ z < 4.49 (thick blue lines), against halo mass. The solid curves correspond to the hot fraction calculated using T_max ≥ 10^5.5 K, whereas the dashed curves correspond to T_gas ≥ 10^5.5 K. The error bars show the 1σ scatter of the fractions.

T_gas and T_max are able to identify the gas particles that are shock-heated when crossing the virial radius; however since we have found that T_max is affected by stellar feedback, we decide to use the gas particle temperature, T_gas, after accretion. We find that a lower limit on T_gas is the most suitable method to calculate f_acc,hot, since it also does not include gas that goes through a shock but cools immediately afterwards and therefore does not contribute to the hot halo formation process.

4.1.1 gadget and Anarchy

In this section, we extend the discussion presented in Section 2.2 and analyse the differences in the hot/cold modes of accretion on to haloes when the formulation of the hydrodynamic scheme is changed. We compare two L050N0752 simulations that use the same subgrid models, one employs the standard SPH code gadget, while the other employs the Anarchy hydrodynamic solver used in the fiducial EAGLE runs. Fig. 10 shows the same as Fig. 9 for z = 0 (top panel) and z = 2 (bottom panel), but the lines correspond to the standard gadget (blue lines) and Anarchy (orange lines) simulations.

Figure 10.

Open in new tabDownload slide

Fraction of gas accreted hot during the redshift ranges 0 ≤ z < 0.1 (top) and 2.0 ≤ z < 2.2 (bottom) against halo mass. The panels show the same as Fig. 9, but in this case, the curves correspond to the hot fraction calculated from the Ref-L050N0752 Anarchy (orange lines) and gadget (blue lines) simulations. The error bars show the 1σ scatter of the fractions. The symbols correspond to the hot fraction estimates of van de Voort et al (2011, blue triangles), Faucher-Giguere et al. (2011, open diamonds) and Nelson et al. (2013, open circles).

The top and bottom panels show excellent agreement between gadget and Anarchy in haloes less massive than 10^11.5 M_⊙ and 10¹² M_⊙, respectively, and modest differences in larger haloes, irrespective of whether we use the T_max hot gas accretion fractions (solid lines) or the T_gas hot gas accretion fractions (dashed lines). The Anarchy simulation exhibits a somewhat larger fraction of hot accretion on to massive haloes than its gadget counterpart. This is expected, since the spurious surface tension appearing in the gadget formulation of SPH prevents the cold dense clumps of gas from being disrupted, mixed and heated when crossing the virial shock (e.g. Schaller et al. 2015).

The top panel of Fig. 10 shows that the T_max hot gas accretion fractions taken from the gadget simulation are in agreement with van de Voort et al. (2011) analysis of the OWLS simulations (Schaye et al. 2010). The bottom panel also shows broad agreement with van de Voort et al. (2011), but substantial differences with Nelson et al. (2013). While van de Voort et al. (2011) used the standard gadget hydrodynamic solver in their simulations, Nelson et al. (2013) analysed two simulation series that employed either the moving mesh code arepo (Springel 2010) or standard gadget, both without stellar, AGN feedback or metal cooling. Nelson et al. (2013) traced the evolution of the gas properties using a Monte Carlo tracer particle technique that enable them to compute T_max. They did not find large differences between gadget and arepo in the cold mode of accretion on to haloes, and concluded that the cold fraction on to haloes mainly depends on the manner (either with T_max or other cut-off temperature) in which it is measured.

The large differences between our work and that of Nelson et al. (2013) is intriguing. Nelson et al. (2013) calculated the accretion rates and so the hot and cold fractions considering only smooth accretion (without including the merger contribution) and over an accretion time window of 1 Gyr. In this work, we did not separate gas accreted smoothly or through mergers and used smaller time windows; however, we find that considering only smooth accretion and/or larger time window does not significantly change the fraction of gas accreted hot (but increases the total rate of gas accretion). In addition, Nelson et al. (2013) used simulations without metal cooling, stellar feedback or AGN feedback. We compared the fraction of hot gas accretion between our reference model and a model without feedback and found that in the simulation without feedback the hot fraction increases from 10 per cent to 20 per cent in the mass range of 10¹² M_⊙ to 10^12.5 M_⊙, and agree in 10¹³ M_⊙. However, this increase is not enough to explain the large differences we find with Nelson et al. (2013). Unfortunately, we do not have a model without metal cooling to test whether this plausible explanation is sufficient to account for the remaining differences.

We also compare our results with the work of Faucher-Giguère et al. (2011), who used a series of cosmological simulations run with the standard SPH code gadget and including stellar feedback and metal cooling but no AGN feedback. They calculated the rates of gas crossing a virial shell, and similar to this work, they used an instantaneous temperature (T_gas > 2.5 × 10⁵ K rather than T_max) to separate hot from cold accretion. We find good agreement at z = 0, but at z = 2, we find large differences. The Faucher-Giguère et al. (2011) transition mass (i.e. the halo mass where the hot mode of accretion equals the cold mode) at z = 2 is between 10^11.3 M_⊙ and 10^11.5 M_⊙ (depending on the stellar feedback).

4.2 Hot/cold fraction

The final ingredient for our model of hot halo formation, which we will present in Section 5, is the fraction of gas accreted on to haloes in the hot and cold modes. Throughout this work, we calculate the fraction of hot mode gas accretion, f_acc,hot(M, z), using T_gas ≥ 10^5.5 K. f_acc,hot can be considered as an indirect measure of the presence of hot gas in the halo, since large values of f_acc,hot imply large values of M_hot/M₂₀₀. Fig. 11 shows f_acc,hot(M, z) at z = 0–0.1 (top left-hand panel), z = 1.0–1.26 (top right-hand panel) and z = 2.0–2.24 (bottom left-hand panel). For each redshift range, we find excellent agreement between the f_acc,hot(M, z) curves taken from simulations with different resolution and box size. We find that f_acc,hot(M, z) increases smoothly with halo mass and decreasing redshift.

Figure 11.

Open in new tabDownload slide

Fraction of hot mode gas accretion during 0 ≤ z < 0.1 (top left-hand panel), 1.0 ≤ z < 1.26 (top right-hand panel), 2.0 ≤ z < 2.24 (bottom left-hand panel) against halo mass. In these panels, the different lines correspond to simulations with different resolutions and box sizes and the error bars show the 1σ scatter. The bottom right-hand panel compares different redshifts (symbols) to the best-fitting expression (lines).

We have also investigated how f_acc,hot(M, z) is affected when we vary the feedback mechanisms. We find that although the total gas accretion rate on to the halo remains nearly unchanged under varying feedback mechanisms (in agreement with van de Voort et al. 2011), f_acc,hot(M, z) increases somewhat at fixed halo mass for the strong stellar feedback and no AGN feedback scenarios. The impact of stellar feedback is largest. For example, strong stellar feedback increases f_acc,hot(M, z) by a factor of 1.2 in 10¹² M_⊙ haloes, whereas weak stellar feedback decreases it by a factor of 1.26. Strong AGN feedback decreases f_acc,hot(M, z) but only in high-mass haloes and by up to a factor of 1.1. As in van de Voort et al. (2011), we find that the impact of feedback mechanisms on the fraction of hot mode gas accretion is small.

In the next section, we derive an analytic model for hot halo formation. In the model, we assume that the halo develops a hot atmosphere depending on the fraction of hot mode gas accretion and on the amount of hot gas already in the halo.

5 TOY MODEL

In Section 3, we investigated the formation of the hot halo in the EAGLE simulations. We found that the development of a strong bimodality in the cooling time PDF of the halo gas provides a clear signature of hot halo formation, which occurs in the halo mass range of 10^11.5–10¹² M_⊙. We noticed however that, even when a stable hot atmosphere has not yet been formed, there is already some hot gas in haloes less massive than 10^11.5 M_⊙. This is because gas can be heated by the extragalactic UV/X-ray background or by shocks with stellar or AGN outflows. In this section, we aim to determine the heating rate of gas produced by accretion shocks only, and the halo mass at which this heating overcomes the cooling. To do so, we present an analytic model for the shock-heating rate that takes into account the hot gas mass already in the halo, and the fraction of gas accretion occurring in the hot mode. With the model, we aim to assess the impact of feedback mechanisms (that change the hot gas mass), as well as filamentary cold accretion (that decreases the hot mode fraction of gas accretion), on the formation of a stable hot halo. The model assumes that the hot halo forms when the heating rate produced by accretion shocks is able to balance the radiative cooling rate.

5.1 Virial heating rate and accretion history

5.2 Cooling rate

5.3 Critical halo mass for the formation of a hot halo

5.3.1 Analytic estimate

In this subsection, we use the energy condition of post-shock gas given by equation (16) and calculate the critical halo mass, M_crit, for which the mechanical heating rate, Γ_heat (equation 19), equals the gas cooling rate, Γ_cool (equation 22). M_crit is the halo mass above which the heating rate exceeds the cooling rate, and as a result the halo develops a stable hot hydrostatic atmosphere. To calculate M_crit, we assume values for the fraction of hot mode gas accretion (f_acc,hot), as well as the fraction of hot gas mass in the halo (f_hot).

The top panel of Fig. 12 shows M_crit as a function of f_hot and f_acc,hot for z = 0 (solid lines) and z = 2 (dashed lines). From the panel, it can be seen that for fixed f_acc,hot, M_crit increases with increasing f_hot. This is because as f_hot increases, so does Γ_cool, and therefore for Γ_heat(⁠|$\propto M_{200}^{2/3}\dot{M}_{200}$|⁠) to be able to balance Γ_cool, M₂₀₀ needs to be larger.

Figure 12.

Open in new tabDownload slide

Top panel: Halo mass obtained by equating Γ_heat and Γ_cool for redshift z = 0 (solid lines) and z = 2 (dashed lines) as a function of f_hot = M_hot/(Ω_b/Ω_m)M₂₀₀. The different colour lines correspond to M_crit calculated assuming fixed values for the fraction of the hot mode gas accretion (f_acc,hot, as indicated in the legends). Bottom panel: Same as the top panel, but in this case M_crit is calculated assuming f_hot = 0.1 (solid lines) and f_hot = 0.5 (dashed lines) as a function of redshift.

For fixed f_hot, M_crit increases with decreasing f_acc,hot. In this case, when f_acc,hot decreases, the heating rate is able to balance the cooling rate only if the accretion rate is large (⁠|$\Gamma _{\rm {heat}}\propto \dot{M}_{200}T_{\rm {vir}}$|⁠), and since the accretion rate increases with halo mass, the halo needs to grow in mass in order to develop a heating rate large enough to keep the gas hot. The top panel of Fig. 12 also shows that for fixed f_acc,hot and f_hot, M_crit at z = 2 is lower than M_crit at z = 0. This can also be explained in terms of the halo accretion rate. If f_acc,hot and f_hot do not change, the heating rate for fixed halo mass still increases because |$\dot{M}_{200}$| and T_vir increase with increasing redshift. As a result, lower mass haloes are able to produce a heating rate that balances the gas cooling rate.

It is challenging to calculate analytically the hot gas mass in the halo, and the fraction of gas accreted hot, as a function of halo mass and redshift. For that reason, we follow our analysis from Section 3 and make the ansatz that a halo develops a hot atmosphere when the hot gas mass is 10 per cent (which is roughly the hot fraction in haloes with masses between 10¹¹–10¹² M_⊙ in the redshift range of 0–6, see the bottom panel from Fig. 7). In the case of the fraction of hot gas accretion, it is known that for fixed halo mass, f_acc,hot is large at low redshift, and it decreases with increasing redshift due to the presence of cold gas accretion from filaments (as shown in Section 4.2). We assume that at z = 0 f_acc,hot ∼ 0.5–1, and obtain that the mass scale of hot halo formation is between 10^11.4and10^11.7 M_⊙. This is in agreement with the analysis from Section 3, where, by visual inspection of the gas cooling time PDF, we concluded that the hot halo forms between 10^11.5and10¹² M_⊙.

We next analyse how M_crit changes with redshift. The bottom panel of Fig. 12 shows M_crit as a function of redshift for f_hot = 0.1 (solid lines) and f_hot = 0.5 (dashed lines). The different colour lines correspond to M_crit calculated assuming fixed values for f_acc,hot (as indicated in the legend in the top panel). It can be seen from the panel that, for any redshift, the higher f_hot and lower f_acc,hot, the higher M_crit. It can also be seen that for any f_hot and f_acc,hot values, M_crit remains roughly constant in the redshift range of 6–2, and then increases. This is possibly due to the rapid drop of the accretion rate (⁠|$\dot{M}_{200}$|⁠, hence Γ_heat) in the redshift range of 0–1 caused by the accelerated expansion of the Universe.

Ocvirk et al. (2008), along with DB06, argued that chemical enrichment has a crucial impact on shock stability, since metallicity, as well as gas density, determines the cooling rate. DB06 and Ocvirk et al. (2008) found that increasing the metallicity increases the critical halo mass for shock stability. We analyse the impact of metallicity on M_crit in Appendix A, where we show that increasing metallicity increases M_crit, but if the hot gas metallicity is lower than 0.1Z_⊙, it does not strongly impact the normalization of M_crit. This is expected, since metal cooling only becomes important for Z ≳ 0.1Z_⊙ (e.g. Wiersma et al. 2009).

5.3.2 Semi-analytic estimate

In the previous subsection, we used the analytic model to analyse how the mass scale of hot halo formation changes with the fraction of hot gas mass in the halo, f_hot, and the fraction of gas that shock-heats when crossing the virial radius, f_acc,hot. Since it is challenging to derive analytically f_hot and f_acc,hot as a function of halo mass and redshift, we assumed typical values and concluded that the hot halo forms in the halo mass range of 10^11.4–10^11.7 M_⊙ at z = 0 and remains roughly constant with redshift. In this section, we make a ‘semi-analytic’ estimate of M_crit (hereafter M_{crit, sa}) by using the best-fitting relations from our simulations for f_hot(M₂₀₀, z) (equations 7–10), f_acc,hot(M₂₀₀, z) (equations 12–15). Note that f_acc,hot(M₂₀₀, z) relation may underestimate the fraction of hot gas accretion for haloes with virial temperatures lower than 10^5.5K (see Section 4.1 for a discussion). In this section, however, we do not intend to predict a mass scale for hot halo formation, but compare the result with the analysis done in the previous subsection.

Fig. 13 shows the semi-analytic estimate of M_crit in black solid lines and the analytic estimates calculated in the previous subsection in coloured dashed lines. By using the best-fitting relations, we obtain a critical mass scale of 10^11.75 M_⊙ at z = 0 that remains roughly constant with redshift in agreement with the analysis done in the previous section.

Figure 13.

Open in new tabDownload slide

Halo mass obtained by equating Γ_heat and Γ_cool, for constant values of f_hot and f_acc,hot (analytic M_crit, coloured dashed lines) and for f_hot and f_acc,hot as a function of halo mass and redshift (semi-analytic M_crit using equations 7–10 and 12–15, black solid line).

We find that f_acc,hot(M_{crit, sa}) ≈ 0.3 at z = 0 and f_acc,hot(M_{crit, sa}) < 0.3 at z > 0, from which we conclude that M_crit does not necessarily correspond to the mass scale where the hot and cold modes of accretion contribute equally (f_acc,hot = 0.5), because equal hot/cold modes of accretion do not imply the existence of a stable hot atmosphere or lack thereof. In fact, f_acc,hot(M_{crit, sa}) < 0.3 at z > 0 implies that massive haloes are able to develop a hot atmosphere (and hence virial shocks), even when they are accreting the majority of gas in the cold mode.

We also find that M_{crit, sa} is somewhat lower than the analytic halo mass calculated by DB06 at r = R₂₀₀. In their work, DB06 found that the critical halo mass of shock-heating occurring in the inner halo at r = 0.1R₂₀₀ is 6 × 10¹¹ M_⊙ and 2 × 10¹² M_⊙ at r = R₂₀₀, in both cases the critical halo mass was nearly constant with redshift. In our case, we focus on shocks occurring at R₂₀₀ and find that M_{crit, sa} is a factor of 3.5 lower and changes slightly with redshift. The small change of M_{crit, sa} with redshift is driven by the interplay between the accretion rate, |$\dot{M}_{200}$|⁠, which increases with redshift, and the fraction of hot gas accretion, f_acc,hot, which decreases with redshift, but it is also due to the fact that we assume a fixed value of 10^0.6ρ_crit for the hot gas density. In Appendix D, we do a more detailed comparison with the work of DB06, and in Appendix C, we discuss how the density of the hot gas in the halo changes with redshift.

5.4 Comparison between t_cool/t_dyn and Γ_heat/Γ_cool

In this subsection, we aim to show that the analytical model gives a better prediction for M_crit compared to the ratio of the halo dynamical and cooling time-scales. To do so, we compute Γ_heat and Γ_cool, as well as t_cool and t_dyn, for each halo in the simulations.

In the calculation of Γ_heat, we compute the individual gas accretion rates for each halo, as well as the fraction of particles that shock-heat, and define hot gas as all particles with temperatures larger than the host halo's T_vir (to avoid using t_cool/t_dyn as we do in section 3.3). In the case of t_cool and t_dyn, we assume that the hot halo is formed when the cooling time at the virial radius is larger than the dynamical time, and to calculate them we only use gas between (0.8 − 1) × R₂₀₀.

The top panel of Fig. 14 shows the ratio between Γ_heat and Γ_cool (olive lines) and t_cool and t_dyn (blue lines). We find that the halo mass where Γ_heat = Γ_cool at z = 0 is in very good agreement with the semi-analytic prediction (M_{crit, sa} = 10^11.8). At z = 2, the halo mass is larger than M_{crit, sa} due to the hot gas density being slightly lower than the fiducial value adopted in Section 5.2. We find that at z = 0, the median gas t_cool is always larger than t_dyn in the halo mass range of 10¹¹–10^12.5 M_⊙, indicating that the hot halo should form, in contradiction with the results from Section 3. Such long cooling times at all masses are due to the presence of additional heating mechanisms (like UV/X-ray background). At z = 2, the gas median cooling times at large radii are shorter, due to the cold, dense filamentary gas, and t_cool overcomes t_dyn only in haloes larger than 10^12.2 M_⊙.

Figure 14.

Open in new tabDownload slide

Top panel: median logarithmic ratio between Γ_heat and Γ_cool (green lines), and between t_cool and t_dyn (blue lines) as a function of halo mass (M₂₀₀) for gas in haloes at z = 0 (solid lines) and z = 2 (dashed lines). Bottom panel: median logarithmic ratio of Γ_heat and Γ_cool calculated from simulations with different feedback prescriptions as indicated in the legend.

We believe that Γ_heat = Γ_cool is a better method to determine when the hot halo forms because, unlike t_cool = t_dyn, Γ_heat by definition only considers the heating due to halo growth and accretion shocks.

5.4.1 Feedback variations

The model presented in this section assumes that the formation of the hot halo is only driven by the heating from gravitational accretion shocks. However, in the presence of other energy sources, like stars or AGN, the heating rate should increase, and therefore extra terms (like Γ_stellar or Γ_AGN) should be added to equation (16). Although we do not include extra heating sources accounting for the presence of feedback, we still find good agreement between the analytical results and the numerical analysis. Also, varying feedback may change the metallicity and the hot gas density (e.g. Crain et al. 2013) and hence Γ_cool. The sign of the effect is however difficult to predict. On the one hand, a more efficient feedback will reduce the stellar and hence the total metal mass. On the other hand, a greater fraction of the metals may reside in the hot halo gas.

We next investigate how the critical halo mass is affected by different feedback implementations. To do so, we calculate Γ_heat and Γ_cool from simulations with feedback variations and show these in the bottom panel of Fig. 14. We find that at z = 0, in the less energetic stellar feedback scenario |$M^{^{\prime }}_{\rm {crit}}$| (at which Γ_heat = Γ_cool) is ∼10^11.8 M_⊙, in the more explosive AGN FB |$M^{^{\prime }}_{\rm {crit}}\sim 10^{11.9}\, {\rm M_{{\odot }}}$|⁠, in the no AGN FB (but moderate stellar feedback) |$M^{^{\prime }}_{\rm {crit}}\sim 10^{11.9}\, {\rm M_{{\odot }}}$|⁠, and in the more energetic stellar feedback case |$M^{^{\prime }}_{\rm {crit}}\sim 10^{12}\, {\rm M_{{\odot }}}$|⁠. This is in excellent agreement with the analytic model, which predicts that M_crit increases when f_hot increases, which occurs when stellar feedback is more energetic and/or there is no AGN feedback. The result is also in good agreement with Section 3.2, where we concluded that the hot halo formation is not as strongly affected by AGN feedback as it is by stellar feedback.

5.5 Comparison with simulated t_cool/t_dyn

In this subsection, we compare the semi-analytic predictions for M_{crit, sa} presented in the previous section with results from our simulations. To do so, we investigate the halo mass range for which haloes develop a hot atmosphere. The model predicts that haloes with masses M₂₀₀ > M_{crit, sa} form a hot atmosphere. We test this by analysing the PDF of t_cool/t_dyn for haloes with masses M₂₀₀ ≪ M_{crit, sa} (the left-hand panel of Fig. 15), M₂₀₀ ∼ M_{crit, sa} (middle panel) and M₂₀₀ ≫ M_{crit, sa} (right-hand panel) at z = 0–4. We select gas particles in the halo (located between r = [0.15 − 1]R₂₀₀) that are not star forming, and use the Ref-L100N1504 simulation.

Figure 15.

Open in new tabDownload slide

Mass-weighted PDF of the ratio of the cooling and local dynamical time for gas in haloes with masses M₂₀₀ ≪ M_crit (left-hand panel), M₂₀₀ ∼ M_crit (middle panel) and M₂₀₀ ≫ M_crit (right-hand panel) at z = 0, 1, 2, 3, 4 and 5, taken from the Ref-L100N1504 simulation.

The left-hand panel of Fig. 15 shows that in haloes with masses M₂₀₀ ∼ 0.1M_{crit, sa} most of the gas has t_cool < t_dyn and is thus able to cool effectively and accrete on to the galactic disc. In larger haloes (M₂₀₀ ∼ M_{crit, sa}, middle panel), the PDF is bimodal with a peak in either side of t_cool/t_dyn = 1, indicating that a hot atmosphere has been formed in these haloes at each redshift, despite the increasing contribution of cold gas from filaments at higher redshifts. In haloes with masses M₂₀₀ ∼ 10M_{crit, sa} (right-hand panel), the increase in shock-heated hot gas enhances the peak at t_cool ≫ t_dyn at the expense of the peak at t_cool ≪ t_dyn. We conclude that the semi-analytic model for hot halo formation captures the mass scale of hot halo formation in the simulations.

6 SUMMARY AND CONCLUSIONS

We have studied the formation of the hot hydrostatic halo, its dependence on feedback mechanisms, on the hot gas mass that is already in the halo, and on the fraction of gas accreted hot using the EAGLE suite of hydrodynamical simulations, as well as analytic calculations.

We began by analysing the PDF of the ratio of the radiative cooling time and the dynamical time for gas in the halo and found that when the hot halo is formed, it produces a strong bimodality in the PDF (Figs 1 and 2, top right-hand panels). By inspection of cooling time PDFs, we found that the mass scale for hot halo formation is 10^11.5–10¹² M_⊙ at z = 0–4 (Fig. 4).

We found, however, that the cooling time PDF is strongly affected by the extragalactic UV/X-ray background radiation and stellar feedback. The UV/X-ray background radiation suppresses the net cooling rate of gas in the temperature range of T ∼ 10⁴–10⁵ K; therefore, the peak of the cooling time PDF shifts towards larger cooling times as the halo's virial temperature decreases to these values (Fig. 3).

In the case of stellar feedback, galactic winds expel gas from the galaxy into the halo, changing the distribution of the hot gas. As a result, more energetic stellar feedback increases the fraction of hot gas. We also analysed the build-up of the total gas mass, M_gas, as well as the hot gas mass, M_hot, as haloes evolve, and concluded that stellar feedback has a large impact on the amount of gas in the halo. For example, doubling the strength of the stellar feedback increases the gas mass fraction by a factor of 1.3 in 10¹² M_⊙ haloes relative to the Ref model, whereas halving the strength of the stellar feedback decreases the gas mass fraction by a factor of 2.5 (Fig. 7).

In the case of AGN feedback, neither the bimodality of the cooling time PDF nor the hot gas mass are strongly affected in haloes smaller than 10¹² M_⊙ since they do not form massive black holes. However, the PDFs and hot gas mass do change at higher halo masses and in a manner opposite from that of stellar feedback. While efficient stellar feedback increases the gas mass in the halo, more explosive AGN feedback decreases it. For example, with strong (without) AGN feedback the total gas mass in the halo decreases (increases) by a factor of 1.5 (Fig. 7). The effect on M_hot is slightly different, efficient stellar feedback increases the hot mass fraction by 10 per cent relative to Ref, but no AGN feedback decreases it by 8 per cent in 10¹² M_⊙ haloes. In the case of less energetic stellar feedback and more explosive AGN feedback, the ratio M_hot/M_gas increases by 10 per cent and 3 per cent (on average) with respect to Ref, respectively, in the halo mass range of 10^11.5–10¹² M_⊙.

In addition to the hot gas in the halo, we calculated the fraction of gas accretion occurring in the hot mode (Fig. 11), f_acc,hot. Rather than using a lower limit on the maximum past temperature to select shock-heated gas particles, as done in most previous works, we used the gas temperature after accretion, T_gas, which yields lower f_acc,hot values than the maximum past temperature. We believe that a lower limit on T_gas is a better method to select hot gas accretion, because it excludes the gas that goes through a shock but cools immediately afterwards or that has not passed through an accretion shock but was heated in the past by stellar feedback and has since cooled, and therefore does not contribute to the formation of a hot halo.

We derived an analytic model for hot halo formation that depends on the fraction of hot gas mass that is already in the halo, f_hot, as well as on the accretion rates and the fraction of gas accretion occurring in the hot mode, f_acc,hot. We assumed that a hot halo develops when the heating rate from accretion shocks balances the radiative cooling rate. We computed M_crit, the critical mass scale above which the hot halo forms, as a function of f_hot and f_acc,hot. We found that M_crit increases with increasing f_hot and decreasing f_acc,hot (Fig. 12). The analytic model yields a mass estimate of M_crit ≈ 10¹²–10^11.5 M_⊙ at z = 0, which agrees with the simulation results.

Because estimating f_hot and f_acc,hot analytically is very challenging, we combined the analytic model with fits to the hot gas mass, and hot mode accretion rates as a function of mass and redshift. We computed a semi-analytic critical mass, M_{crit, sa}, and found that M_{crit, sa} = 10^11.75 M_⊙ at z = 0. At higher redshift, M_{crit, sa} remains roughly constant (Fig. 13). We tested the M_{crit, sa} values by inspecting the cooling time PDF of hot gas in haloes with mass M_{crit, sa} and confirmed that at all redshifts, the PDF has a clear bimodal shape (Fig. 15). Note that because the semi-analytic model uses our simulation results as input, its prediction for the mass scale for hot halo formation cannot be tested using the same simulations. We can, however, use it to compare with the analytic analysis.

We compared the ratio of the heating due to accretion (Γ_heat, derived in the analytic model) and the radiative cooling (Γ_cool) rates of hot gas, with the ratio of the cooling (t_cool) and dynamical (t_dyn) times of gas at the virial radius, and found that unlike Γ_heat/Γ_cool, the median t_cool/t_dyn of gas is always greater than unity in the halo mass range of 10¹¹–10^12.5 M_⊙ (Fig. 14). On the contrary, Γ_heat/Γ_cool is only greater than unity for haloes more massive than 10^11.8 M_⊙, indicating that this ratio better captures the heating due to halo growth and accretion shocks. We believe that compared with the ratio of cooling and dynamical times, the analytic model of hot halo formation is better at indicating when a hot halo forms.

Finally, we investigated how feedback impacts the hot halo formation mass scale. We calculated Γ_heat/Γ_cool using simulations with different feedback prescriptions, and concluded that while a hot hydrostatic atmosphere forms in more (less) massive haloes in scenarios with more (less) energetic stellar feedback, the mass scale of hot halo formation is not strongly affected by AGN feedback. This result is driven by the dependence of Γ_heat and Γ_cool on the hot gas mass fraction, f_hot. When f_hot increases (i.e. in the strong stellar feedback scenario), so does the rate of cooling and therefore the halo needs to grow in mass in order to develop a heating rate that overcomes the cooling rate.

Cosmological hydrodynamical simulations have shown that the manner in which galaxies accrete gas depends on the complex interaction between the hot halo, AGN and stellar feedback (see e.g. van de Voort et al. 2011; Faucher-Giguère et al. 2011; Nelson et al. 2015). In Paper II, we will make use of the semi-analytic calculations presented in this work and derive a semi-analytic model for gas accretion on to galaxies that accounts for the hot/cold modes of gas accretion on to haloes and for the rate of gas cooling from the hot halo. By doing so, we aim to provide some insight into the physical mechanisms that drive the gas inflow rates on to galaxies.

Acknowledgements

We are grateful to the referee for fruitful comments that substantially improved the original manuscript. This work used the DiRAC Data Centric system at Durham University, operated by the Institute for Computational Cosmology on behalf of the STFC DiRAC HPC Facility (www.dirac.ac.uk). This equipment was funded by BIS National E-infrastructure capital grant ST/K00042X/1, STFC capital grant ST/H008519/1, and STFC DiRAC Operations grant ST/K003267/1 and Durham University. DiRAC is part of the National E-Infrastructure. The EAGLE simulations were performed using the DiRAC-2 facility at Durham, managed by the ICC, and the PRACE facility in France at TGCC, CEA, Bruyeres-le-Chatel. This work was supported by the European Research Council under the European Union's Seventh Framework Programme (FP7/2007-2013)/ERC Grant agreement 278594-GasAroundGalaxies and by the Netherlands Organisation for Scientific Research (NWO) through VICI grant 639.043.409. RAC is a Royal Society University Research Fellow. We also acknowledge support from STFC (ST/L00075X/1).

REFERENCES

APPENDIX A: SHOCK ANALYSIS

In this section, we analyse the mass-weighted PDF of T_max, T_gas, S_gas. Fig. A1 shows the PDFs for the redshift interval 0.0 < z ≤ 0.1 (top left-hand panel) and 2.0 < z ≤ 2.2 (top right-hand panel). The curves are coloured according to the colour bars at the top of the figure, which indicate the halo mass (and virial temperature) of the halo that gas is accreted on to. It can be seen that the T_max PDF varies with M₂₀₀, being unimodal in low-mass haloes and bimodal in high-mass haloes. The location of the local minimum of the bimodal distribution also changes with M₂₀₀, going from |$T_{\rm {max,min}}\sim 10^{5.5}\hspace{2.84526pt}\rm {K}$| in 10¹² M_⊙ haloes to |$T_{\rm {max,min}}\sim 10^{6}\hspace{2.84526pt}\rm {K}$| in 10¹⁴ M_⊙ haloes. Besides this local minimum, there is a local maximum at |$10^{7.5}\hspace{2.84526pt}\rm {K}$| for all halo masses. This peak is produced by stellar feedback instead of accretion shocks. Some of the gas that is ejected out of the halo due to stellar feedback is eventually re-accreted. However, if it does not reach a temperature larger than |$10^{7.5}\hspace{2.84526pt}\rm {K}$| when crossing R₂₀₀, T_max is not updated, and the gas will be considered hot mode accretion by the maximum temperature criterion. When applying a T_max criterion to separate hot from cold accretion, rather than calculating a T_{max, min} threshold value that changes with M₂₀₀, we follow previous works from the literature (e.g. van de Voort et al. 2011; Nelson et al. 2013) and use |$T_{\rm {max}}=10^{5.5}\hspace{2.84526pt}\rm {K}$|⁠.

Figure A1.

Open in new tabDownload slide

PDF of maximum past temperature (top panels), temperature (middle panels) and entropy (bottom panels) of gas accreted on to haloes in the redshift ranges z = 0–0.1 (left-hand panels) and z = 2.0–2.2 (right-hand panels). The curves are coloured according to the colour bars at the top of the figure, which indicate the halo mass and virial temperature.

As per T_max, the middle and bottom panels of Fig. A1 show that the T_gas and S_gas PDFs have a bimodal shape. To identify the accreted gas that does not cool immediately after the shock, we analyse the post-shock gas temperature and entropy. We select the gas particles that were accreted during the redshift interval z_i − z_j (z_i < z_j), and are hot using a temperature and entropy threshold value. We calculate the gas temperature and entropy mass-weighted PDFs at redshift z_i. We use the local minima of the bimodal distribution as the threshold values (⁠|$T_{\rm {min}}=10^{5.5}\hspace{2.84526pt}\rm {K}$| and |$S_{\rm {min}}=10^{7.2}\hspace{2.84526pt}\rm {K}\hspace{2.84526pt}\rm {cm}^{2}$|⁠) to calculate the fraction of the gas accreted hot. Our motivation for using the gas entropy (besides temperature) to identify shocked gas is based on the fact that gas generally undergoes a large entropy increase when it encounters a shock (Brooks et al. 2009).

APPENDIX B: METALLICITY

In Section 5, we derived an analytic model for hot halo formation. The model predicts that a hot halo develops when the heating rate from accretion shocks balances the radiative cooling rate. To calculate the gas cooling rate, we assumed a constant hot gas metallicity, |$Z_{\rm {hot}\hspace{1.42271pt}\rm {gas}}$|⁠, of 0.1 Z_⊙. However, it has been argued that chemical enrichment has a crucial impact on shock stability (e.g. DB06; Ocvirk et al. 2008). Therefore, in this section, we investigate how the hot gas metallicity in haloes changes with halo mass and redshift, and analyse how the halo mass scale, M_crit, at which the heating rate balances the gas cooling rate, depends on metallicity.

Using the Ref-L100N1504 simulation, we define hot gas as all gas that is within R₂₀₀ and that has a cooling time longer than the local dynamical time, and calculate the mass-weighted median metallicity per halo mass bin. The left-hand panel of Fig. B1 shows the |$Z_{\rm {hot}\hspace{1.42271pt}\rm {gas}}{\rm -}M_{200}$| relation for various output redshifts (here |$Z_{\rm {hot}\hspace{1.42271pt}\rm {gas}}$| is normalized by solar metallicity, which we assume to be Z_⊙ = 0.0129). The error bars show the 16–84th percentiles and each bin contains at least five haloes. Interestingly, at fixed halo mass the mass-weighted median metallicity of the hot gas slightly increases with redshift (i.e. by up to a factor of 1.6 in 10¹² M_⊙ haloes). Many studies of metal abundance have shown that galaxies tend to have lower metallicities at higher redshift (see e.g. Prochaska et al. 2003; Nagamine et al. 2004; Savaglio 2006; Kulkarni et al. 2007; Péroux et al. 2007). It is important to note that the left-hand panel does not show the median metallicity of the ISM (Z_ISM), but of the hot mostly ionized gas in the halo. In the case of the ISM, we obtain Z_ISM–M_* relations (with M_* stellar mass) in agreement with the galaxy mass–metallicity relation from Andrews & Martini (2013) and Zahid et al. (2013), as recently shown by Somerville & Davé (2015). As expected, Z_ISM at fixed stellar mass decreases with increasing redshift.

Figure B1.

Open in new tabDownload slide

Left-hand panel: Mass-weighted median metallicity of the hot gas in the halo as a function of halo mass for various output redshifts, as indicated in the legends. Hot gas is defined as gas within R₂₀₀ that has a cooling time longer than the local dynamical time. The error bars in the panel correspond to the 16–84th percentiles (i.e. 1σ scatter) and each bin contains at least five haloes. Right-hand panel: Halo mass obtained by equating Γ_heat and Γ_cool for redshift z = 0 as a function of f_hot = M_hot/(Ω_b/Ω_m)M₂₀₀. The different colour lines and types correspond to M_crit calculated assuming fixed values for the fraction of the hot mode gas accretion (f_acc,hot) and metallicity, respectively, as indicated in the legends.

We find that for all halo masses, the median mass-weighted metallicity of the gas in the halo increases towards the halo centre (in agreement with Ocvirk et al. 2008; van de Voort & Schaye 2012). Interestingly, we find that the diffuse hot gas in the halo has lower median mass-weighted metallicity (by up to 0.5 dex in haloes larger than 10¹¹ M_⊙) than the cold gas. However, if we define hot gas as all gas within R₂₀₀ that has a maximum past temperature lower than 10^5.5 K, we obtain that the median mass-weighted metallicity of the hot gas is higher (also by to 0.5 dex in haloes larger than 10¹¹ M_⊙), as in van de Voort & Schaye (2012).

We next analyse how M_crit changes with metallicity. The right-hand panel of Fig. B1 shows M_crit at z = 0, calculated assuming constant values of the fraction of hot mode gas accretion, f_acc,hot, and metallicity (as indicated in the legends), as a function of the fraction of hot gas mass in the halo, f_hot. The panel shows that increasing metallicity increases M_crit, but if the hot gas metallicity is lower than 10⁻¹Z_⊙, it does not strongly impact on the normalization of M_crit. This is expected, since metal cooling only becomes important for Z ≳ 0.1Z_⊙ (e.g. Wiersma et al. 2009). The left-hand panel of Fig. B1 shows that the typical metallicity reached by hot gas in the redshift range of 0–4 is ∼0.1 Z_⊙ in ∼10¹² M_⊙ haloes, we then find that assuming |$Z_{\rm {hot}\hspace{1.42271pt}\rm {gas}}\sim 0.1\, {\rm Z}_{{\odot }}$| in the analysis of Section 5.3 is a good approximation.

APPENDIX C: DENSITY

In Section 5, we derived an analytic model for hot halo formation, which considers the heating and cooling rates of gas in the halo. When calculating the cooling rates, we assumed that the hot gas density is about 10^0.6ρ_crit. In this section, we analyse the density and temperature of the hot gas in 10^11.5–10¹² M_⊙ haloes.

Fig. C1 shows the temperature versus density of gas in 10¹² M_⊙ (top left-hand panel) and 10^11.5 M_⊙ haloes (top right-hand panel) at z = 0, and in 10¹² M_⊙ (bottom left-hand panel) and 10^11.5 M_⊙ haloes (bottom right-hand panel) at z = 2.2. The red contours indicate the distribution of the hot gas (i.e. with t_cool ≥ t_dyn), while the blue contours indicate the cold gas (i.e. t_cool < t_dyn). It is interesting to see that there is gas with temperatures below 10⁵ K (and a large spread in density) that has a net cooling time longer than the local dynamical time. This gas is close to the equilibrium temperature obtained when photoheating balances radiative cooling rate. We define the hot halo gas as all gas with densities below 10²ρ_crit and temperature above 10⁵ K. At z = 0, we find that the mean density of the hot gas in the halo is around 10^0.6ρ_crit in both 10¹² M_⊙ and 10^11.5 M_⊙ haloes. At z = 2.2, the mean density increases slightly to 10^0.8ρ_crit for both halo masses. We find that the hot gas density does not change significantly with the host halo mass, nor with redshift (we also did the analysis for z = 3 and 4), we then assume that the mean density of the hot gas is 10^0.6ρ_crit and use this value to calculate of the mass scale for the hot halo formation.

Figure C1.

Open in new tabDownload slide

Temperature versus density of gas in 10¹² M_⊙ (top left-hand panel) and 10^11.5 M_⊙ haloes (top right-hand panel) at z = 0, and in 10¹² M_⊙ (bottom left-hand panel) and 10^11.5 M_⊙ haloes (bottom right-hand panel) at z = 2.2. The contours in the panels enclose 25, 50 and 75 per cent of the distribution and the different colours correspond to the distribution of hot gas (i.e. with t_cool ≥ t_dyn, red contours) and cold gas (blue contours). Note that the cooling time is computed from the net cooling rate, i.e. the difference between the photoheating and radiative cooling rate, and depends on metallicity.

APPENDIX D: COMPARISON WITH THE DEKEL AND BIRNBOIM MODEL

We find that our model is not equivalent to that of DB06 due to the different redshift dependence of t_comp and t_heat. While |$t^{-1}_{\rm {comp}}\propto [\Omega _{\rm {m}}(1+z)^{3}+\Omega _{\Lambda }]^{1/2}$|⁠, |$t^{-1}_{\rm {heat}}\propto (1+z)[\Omega _{\rm {m}}(1+z)^{3}+\Omega _{\Lambda }]^{1/2}$|⁠. However, we obtain similar results at all redshifts when comparing to DB06's critical mass for the shock at 0.1R₂₀₀ (6 × 10¹¹ M_⊙). We believe that we have improved upon the DB06 model by introducing a dependence on the amount of shock-heated gas, which we find to decrease with increasing redshift at fixed halo mass due to the presence of cold filaments (see panels in Fig. 9). To include the impact of cold filaments in their calculations, DB06 had to modify the gas density and assume ρ_stream/ρ_vir ∼ (3M_*/M)^−2/3 (with ρ_stream the filamentary gas density, ρ_vir the gas density at the virial radius and M_* the non-linear clustering mass scale). By doing so they obtained an upper limit for cold streams that increases with increasing redshift, in agreement with our results, but reaches 10¹⁴ M_⊙ at z = 3.5.