The origin of cold gas in the Circumgalactic Medium

The presence of cold ($T \lesssim 10^4$ K) gas in the circumgalactic medium (CGM) of galaxies has been confirmed both in observations and high-resolution simulations, but its origin still represents a puzzle. Possible mechanisms are cold accretion from the intergalactic medium (IGM), clumps embedded in outflows and transported from the disk, gas detaching from the hot CGM phase via thermal instabilities. In this work, we aim at characterizing the history of cold CGM gas, in order to identify the dominant origin channels at different evolutionary stages of the main galaxy. To this goal, we track gas particles in different snapshots of the SPH cosmological zoom-in simulation Eris2k. We perform a backward tracking of cold gas, starting from different redshifts, until we identify one of the followings origins for the particle: cold inflow, ejected from the disk, cooling down in-situ or stripped from a satellite. We also perform a forward tracking of gas in different components of the galaxy (such as the disk and outflows). We find a clear transition between two epochs. For $z>2$, most cold gas (up to 80%) in the CGM comes from cold accretion streams as the galaxy is accreting in the"cold mode"from the IGM. At lower $z$, gas either cools down in-situ after several recycles (with 10-20% of the gas cooling in outflow), or it is ejected directly from the disk (up to 30%). Outflows have a major contribution to the cold CGM gas budget at $z<1$, with almost 50% of hot gas cooling in outflow. Finally, we discuss possible mechanisms for CGM cooling, showing that the thermally unstable gas with $t_{\rm cool}/t_{\rm ff}<1$ (precipitation-regulated feedback) is abundant up to $r\sim 100$ kpc and cooling times are shorter than 50 Myr for densities $n>10^{-2}\,{\rm cm}^{-3}$.


Introduction
The circumgalactic medium (CGM) of a galaxy includes all the gas in the potential well of the dark matter halo that extends outside the disk and up to a few times the virial radius.There is evidence that this gas has a complex multi-phase structure, including hot virialized gas at T ∼ 10 6−7 K (Allen et al. 2008;Dai et al. 2010) and cold gas with hydrogen in atomic and molecular form (T < 10 5 K, Prochaska et al. 2013;Battaia et al. 2016;Meiksin et al. 2017;Fossati et al. 2021).
The presence of cold gas when the virial temperature exceeds 10 6 K is currently a puzzle.The cooling times of low-density gas are generally very long, and hence mechanisms to locally increase the density are needed.In principle, gas could detach from the hot phase because of density perturbations triggering thermal instability (Field 1965;McCourt et al. 2012).Density can also increase due to convective instability (Balbus & Soker 1989;Malagoli et al. 1990), wind-driven shocks (Zubovas & King 2014;Ferrara & Scannapieco 2016;Richings & Faucher-Giguere 2018), or a condensation-precipitation mechanism (Voit et al. 2015(Voit et al. , 2017)).Another possible scenario is that cold gas is lifted directly from the disk and entrained into an outflow, which could be an explanation for the amount of outflowing molecular gas in active galactic nuclei (Cicone et al. 2014).This gas is subject to strong shocks due to the extremely turbulent environment and surrounded by a high-velocity hot medium, so clouds can be disrupted because of thermal conduction and Kelvin-Helmholtz instabilities (McKee & Begelman 1990;Ferrara & Scannapieco 2016).Furthermore, the clouds are exposed to the intense ultraviolet radiation field from the center of the galaxy that is able to dissociate molecular gas, making it challenging for clouds to survive and travel up to distances of a few kiloparsecs (Decataldo et al. 2017(Decataldo et al. , 2019)).Finally, cold gas can be found inside accretion filaments that fuel the galaxy disk.Observations (Rubin et al. 2010(Rubin et al. , 2015;;Lehner et al. 2014) and simulations (Fardal et al. 2001;Kereš et al. 2005;van de Voort et al. 2012;Nelson et al. 2013) agree that a large fraction of gas accretes into galaxies via a "cold mode," especially in smaller galaxies where accretion shocks might not be able to heat up the rapidly cooling gas in dense filaments.
Numerical simulations are a powerful tool for improving our understanding of the origin of cold gas.Simulations can be done in two ways.First, it is possible to study the gas on small scales to understand the evolution of thermal instabilities and the dynamics of cold gas in different environments, for example simulations of stratified media (Joung et al. 2012), filamentary streams (Mandelker et al. 2019), and the cloud-crushing problem (Scannapieco & Brüggen 2015;Schneider et al. 2018;Gronke & Oh 2018).The second way is by using cosmological simulations, which allow us to visualize the large-scale structure, showing where the cold gas is more abundant and where it comes from; one approach is to track the gas trajectory, and hence determining if the CGM gas cools down in situ or if it is carried already in cold phase within accretion streams or outflows.By analyzing simulations of large cosmological boxes (Oppenheimer et al. 2018;Hafen et al. 2019;Nelson et al. 2020) it is possible to obtain statistics over many halos of different masses.On the other hand, zoom-in simulations allow us to reach a greater level of detail on a single halo, with the inclusion of more detailed physics (Peeples et al. 2019;Suresh et al. 2019).
Nevertheless, resolving the cold gas in the CGM is a challenging task, even for the most modern simulations.In particular, the cold gas mass, and the abundance of low ionization potential ions as a consequence, are typically under-reproduced in comparison to the observational results (e.g., the COS-Halos survey; Tumlinson et al. 2011).It has been proposed (Hummels et al. 2019) that this is due to the lack of resolution since cold gas is usually organized in dense clumps and clouds that cannot be resolved (Rauch et al. 1999;Crighton et al. 2015;Stern et al. 2016;Rubin et al. 2018), ending up mixed with the hot phase.Furthermore, since Lagrangian methods and common refinement schemes in Eulerian methods focus on resolving the dense gas, the diffuse phase of the CGM is generally poorly resolved.Some simulations manage to find a compromise between the improved resolution and the increase in computational cost, adopting refinement criteria that explicitly target the CGM, for example additional spatial refinement up to the virial radius (Hummels et al. 2019;van de Voort et al. 2019;Peeples et al. 2019) and super-Lagrangian zoom schemes (Suresh et al. 2019).
Most of past theoretical works discussed the origin of the CGM gas, without focusing on the cold phase.For example, Ford et al. (2014) track particles backward and forward from a simulation snapshot at z = 0.25 to distinguish between gas in the process of accretion, gas ejected at a previous time, and ambient gas that always belonged to the CGM.They find ambient gas to be the dominant component (up to 85%), with outflows becoming more important (up to 43%) for low-mass halos (< 10 11.5 M ⊙ ).A fraction of outflowing gas is later re-accreted, and this cycle can happen multiple times (gas recycling; Oppenheimer et al. 2010;Anglés-Alcázar et al. 2017).In a similar analysis over a set of FIRE-2 simulations, Hafen et al. (2019) classify 60-80% of gas as intergalactic medium (IGM) accretion and ∼ 10-30% from winds, but the IGM accretion category includes any gas particles not belonging to other galaxies, hence also including the ambient gas in Ford et al. (2014).An analysis focused on cold gas was presented in Suresh et al. (2019), where galactic fountains were identified as the major source, with most gas recycled more than once.Interestingly, a comparison between two simulations with and without winds shows that winds are able to trigger gas cooling, increasing by ∼ 50% the HI covering fraction.Finally,  2020) identify cold gas clouds in the TNG50 cosmological box (Pillepich et al. 2019) at redshift z ∼ 5, and follow their assembling history.They find that clouds form by accretion of gas into dense pre-existing seeds, which mostly originate from density perturbations induced by mergers with infalling satellite galaxies.
In this work we performed an analysis on the origin and fate of the cold CGM in gas in the cosmological zoom-in simulation Eris2k.Differently from most previous works that typically track the cold gas back in time until some location criteria are met, we consider simultaneously the position and temperature evolution of gas particles in order to understand when and where gas cooling happens.In particular, we distinguish between gas cooling down directly in the CGM (while in outflow, in inflow, or static) and gas entering the CGM already in a cold state (from the IGM, from the disk of the main galaxy, or from a satellite).This allows us to investigate the mechanisms responsible for gas cooling, and the survival of cold gas clouds transported into the CGM.
The paper is organized as follows.In Sect. 2 we present the simulation Eris2k, and describe our algorithm to track smoothed particle hydrodynamic (SPH) particles.In Sect. 3 we detail the properties of cold gas in Eris2k (Sect.3.1), we discuss its origin (Sect.3.2), and we analyze in more detail the gas temperature evolution by forward tracking gas in the disk (Sect.3.3) and in outflow (Sect.3.4).A discussion of the cooling physics and cooling times is presented in Sect.3.5.We summarize our conclusions in Sect. 4.

The Eris2k simulation
In this work we analyze the cosmological zoom-in simulations of a Milky Way-like spiral galaxy, Eris2k.The simulation was run with the SPH code Gasoline (Wadsley et al. 2004), assuming flat WMAP-3 cosmology (Ω M = 0.24, Ω b = 0.042, H 0 = 73 km s −1 Mpc −1 , n = 0.96, σ 8 = 0.76).The mass resolution is m dm = 9.8 × 10 4 M ⊙ for dark matter and m g = 2 × 10 4 M ⊙ for gas particles.The ionization fractions of H and He were computed taking into account photoionization by a uniform UV background (Haardt & Madau 2012), recombination, and collisional ionization, with rates from Abel et al. (1997).Metal-line cooling rates were pre-computed with the code Cloudy (Ferland et al. 1998) on the assumption of ionization equilibrium (Shen et al. 2010).Star particles, representing a population of stars with a Kroupa (2001) initial mass function (IMF), are formed stochastically by converting dense (n > 100 cm −3 ) and cold (T < 10 4 K) gas particles.Supernova feedback is included by injecting thermal energy (with the delayed cooling approximation) and metals into the gas.
At z ≃ 0.3, the main halo in Erisk2k (total mass M = 7.2 × 10 11 M ⊙ , virial radius R vir = 206 kpc) hosts a Milky Way-like galaxy with 4.9 × 10 10 M ⊙ in gas and 3.9 × 10 10 M ⊙ in stars.The galaxy has a quiet merging history, with no major mergers (ratio > 1 : 10) happening for z ≲ 3. The star formation rate is quite constant (∼ 10 M ⊙ yr −1 ) from z ≃ 5 to z ≃ 2, and then drops to a value of about 2 M ⊙ yr −1 at z ≃ 0.3.
In Fig. 1 we show the spatial resolution (i.e., the smoothing length of the particles) in Eris2k, compared to other simulations that have been adopted in other CGM studies.The median spatial resolution is computed over radial bins (blue curve), with the variation from the 10th to the 90th percentile shown as a shaded area.As is natural for SPH simulations, the highest resolution is obtained where the most particles are: in the disk (r < 10 − 30 kpc).The resolution degrades as we move farther from the center, to a value λ ∼ 30 − 100 pc at a radius R/R vir ≃ 0.1 − 0.2, where most of the cold gas is found (see Sect. 3.1), which is approximately the size of a giant molecular cloud.This resolution is comparable to that of other works focusing on the CGM, using different codes and numerical techniques, such as the superrefinement zoom-in simulations of Suresh et al. (2019) and the cosmological box Illustris TNG50 (Genel et al. 2014;Nelson et al. 2020).It is possible that cold gas clumps and filaments in the CGM could be smaller, implying that our analysis on the origin of cold gas is limited to the scale we are able to resolve.Further convergence studies are needed to investigate the importance of smaller-scale mechanisms (e.g., Kelvin-Helmholtz and Rayleigh-Taylor instabilities) on gas condensation.

Ion abundances and column density maps
To make a comparison with observations of cold gas around lowredshift galaxies, we computed the ion abundances for two common tracers, HI and Mg II, from our simulation, and their column density distributions.The HI abundance in the simulation is computed on-the-fly through a primordial chemical network including six species (e − , HI, HII, HeI, HeII, HeIII), accounting for collisional ionization, recombination, and photoionization by the Haardt & Madau (2012) UV background, with corresponding heating and cooling rates (Abel et al. 1997).Gas self-shielding is not taken into account in the simulation.The abundance of MgII is obtained interpolating from Cloudy tables, in the same way as in Shen et al. (2010Shen et al. ( , 2013)).Tables were obtained by running single-zone models (i.e., thin slabs) for a range of densities and temperatures, assuming photoionization equilibrium under the Haardt & Madau (2012) uniform UV background (the same background used for HI and MgII abundance calculations).

Particle tracking
SPH particles in the simulations are assigned a unique ID, so that it is possible to track them across snapshots.In this way Fig. 2.After selecting cold gas particles in the CGM (hence excluding gas in the disk of the main galaxy and gas in satellites) at a chosen redshift z, their origin was discovered by tracking back in the previous snapshots of the simulation.The green boxes represent location criteria: particles coming from the disk, coming from a satellite, or always present in the CGM.In the red boxes are shown particles in the CGM that have always been cold (hence in cold inflow since they must have been accreted directly in cold phase at some moment in the galaxy evolution) and particles that cooled down from above the threshold temperature (T = 3 × 10 4 K).Finally shown are static particles (radial velocity |v r | < 50 km s −1 ) and particles with inflowing or outflowing velocity.
we are able to check the variation in the physical properties of parcels of gas (with the mass equal to the mass resolution of the simulation, i.e., around 2 × 10 4 M ⊙ ) across the evolution of the halo.In particular, we are interested in the temperature evolution to understand how cold gas in the CGM forms, and we follow its trajectory in the halo as well.
We only analyze the main halo of the simulation, which is identified by using the Amiga Halo Finder (AHF) (Knollmann & Knebe 2009).Then we find the center of mass of the halo, and we define the disk and CGM in the following way.For simplicity, we classify as "disk" all the gas within 10 kpc of the center, which is approximately twice the half-mass radius of the stellar component (evolving from ∼ 7 kpc at z = 2 to ∼ 13 kpc at z = 0.3).The CGM includes all the gas at a distance 10 kpc < r < 1.5R vir from the center (i.e., up to a radius of ∼ 110 kpc at z = 2 and ∼ 310 kpc at z = 0.3).Throughout the paper, we refer to cold gas as gas with temperature T < 3 × 10 4 K (as in Sokołowska et al. (2018)) in order to include all neutral gas and gas bearing low ionization potential ions (e.g., HI, MgII, SiII).
To understand the origins of cold gas at a given redshift z, each SPH particle that belongs to the cold CGM (T < 3 × 10 4 K, r > 10 kpc) at that redshift is tracked back in the previous snapshots until one of the possible origin channels is identified (see Fig. 2).First (green boxes), we identify whether the particle comes from the disk or from a satellite (i.e., belongs to a satellite of the main halo in the AHF catalog).If not, then the particle could cool down in the CGM from a hot phase or be accreted di-rectly in cold phase (red boxes).Finally, among particles cooling down, we distinguish those with an outflowing velocity (radial velocity v r > 50 km s −1 ) or inflowing velocity (v r < 50 km s −1 ), or that are static (−50 km s −1 < v r < 50 km s −1 ).
Moreover, to understand the trajectory and thermal evolution of gas in the CGM, we track forward gas particles from specific regions of the galaxy.In particular, we consider three setups: forward tracking of cold disk gas, to understand to what extent it contributes to the cold gas in the CGM, by ejection from the disk; forward tracking of cold outflowing gas (T < 3 × 10 4 K, v r > 50 km s −1 ), to understand if this gas survives and reaches the CGM without being heated; and forward tracking of hot outflowing gas (T > 3 × 10 4 K, v r > 50 km s −1 ), to understand if cooling in the outflow occurs.

Overview of the cold gas in Eris2k
We begin by presenting the main thermal properties of gas in the largest halo in the simulation Eris2k.A map of gas temperature, for a thin slice passing through the center of the box, is shown in the upper left panel of Fig. 3.The upper right panel of the same figure shows only the cold gas (i.e., with temperature lower than 3 × 10 4 K in our definition).Cold gas is present in the form of numerous small clumps distributed mainly in the region within 100 kpc of the center.These clumps are not easily seen in the left panel, due to the large amount hot gas in which they are embedded.The temperature can go as low as 10 2 K in the clumps nearer to the center (i.e., clumps in the disk of the galaxy).The mass of the cold gas in the CGM is 7.6 × 10 9 M ⊙ , representing ∼ 15.7% of the total CGM gas mass; instead, cold gas makes up ∼ 40% of the total gas in the disk.
The lower panels of Fig. 3 show the mass-weighted distribution of gas in the density-temperature plane (lower left panel) and the radius-temperature plane (lower right panel).Cold gas is present with higher mass in the region at T ≃ 10 4 K, corresponding to gas in photoionization equilibrium by the UV background, and in the tail of dense gas going down to densities n > 1 cm −3 .Gas with high density is favored to cool down because of shorter cooling times (see Sect. 3.5).In the lower right panel of Fig. 3, we can see how most of the cold gas is located in the central region, as was already observed from the temperature slice.Nevertheless, we note that there is a substantial amount of gas below 10 4 K within a relatively large radius, beyond the disk region, up to ∼ 50 kpc.Cold gas is also present at larger distances, without showing any specific preference at any radius.
Figure 4 shows the column densities of two main cold gas tracers, HI and MgII, computed for the main halo in Eris2k at z ≃ 0.3.Most of cold gas is clearly seen in the central region of the galaxy (r < 20 kpc), associated with the disk, with column densities N HI ≃ 10 21−22 cm −2 and N MgII ≃ 10 17−18 cm −2 .Nevertheless, co-spatial clumps of HI and MgII are also seen up to r ∼ 100 kpc, in the CGM of the galaxy, with sporadic clumps also visible at the virial radius and beyond.In Fig. 5 we show the radial profiles of N HI and N MgII , averaging over annuli of radial separations of 0.5 kpc.The profiles are flat in the central region, corresponding to the disk, and decline quickly for r > 0.1 R vir .The N HI profile presents a plateau for r > R vir , which is not seen for the N MgII profile, declining below 10 7 cm −2 at R vir .This points to a more extended distribution of HI with respect to MgII, as also found in other analyses of simulated halos of similar mass (Oppenheimer et al. 2018;Nelson et al. 2020).Points representing the data from L * galaxies in the COS-Halos survey are shown in the plot (Werk et al. 2013;Prochaska et al. 2017).The data present a large scatter and a comparison with Eris2k is not straightforward since the observed objects span a wide range of properties, such as different halo masses and star formation rates.Nevertheless, by scaling the impact parameter of the observation with the virial radius, we find many galaxies with column densities compatible to those obtained from the simulation.
Cold gas extending up to the virial radius is also found in other works focusing on the thermal properties of the CGM, for simulations run with different codes and different refinement schemes for the CGM gas.Hummels et al. (2019) pointed out that the abundances of HI and low ionization potential ions increase, at the expense of high ionization potential (as OVI), when the CGM is better resolved.Cold HI-bearing gas is found in small and dense clumps, that end up mixed with the hot phase in the low-resolution runs.In a similar fashion, Peeples et al. (2019) showed that resolving smaller structures in the CGM generally increases the HI column density.Suresh et al. (2019) analyzed two zoom-in simulations, finding that the majority of dense (n > 10 −3 cm −3 ) gas in the galaxy is cold (T < 10 5 K) up to 0.6 − 0.7 R vir , with a significant fraction extending up to R vir .Nelson et al. (2020) resolved individual cold clumps by producing HI and MgII maps for massive halos (M > 10 13 M ⊙ ) in the simulation box TNG50, though the CGM is not explicitly targeted in the refinement criteria.

The origin of cold gas
Starting from different redshifts, we tracked back the cold SPH gas particles in the CGM to determine one of their possible origins.We identified different possible origins, which are summarized in Fig. 2. First, we adopted a criterion based on the location of the particle to understand if it comes cold from the disk or from a satellite.Otherwise, the particle could have been accreted directly in cold phase or cooled down in the CGM from a hot phase.In the latter case, we also distinguished whether the particle cools down with inward velocity (inflowing, v r < −50 km s −1 ) or outward velocity (outflowing, v r > 50 km s −1 ), or is static (−50 km s −1 < v r < 50 km s −1 ).
After applying this algorithm to the cold CGM gas at z ≃ 0.3, we found a selection of particles with different origins, and shown in Fig. 6.In each row, there is a selection of 30 SPH particles labeled with the origins "cooling CGM, outflowing," "cooling CGM, inflowing," "from disk," and "from satellites."In the left column we show the trajectories projected in the 2D plane from z ≃ 0.5 to z ≃ 0.3, with dots becoming progressively smaller and fainter as we go to higher redshifts.Particles are shown in blue when they are cold and in red when they are hot.In the right column the distance of each particle from the center is shown as a function of redshift.We note that the behavior of particles cooling in outflow and coming from the disk is similar: all of them are located very close to the disk (reaching at most a distance r ∼ 25 kpc from the center), and during their history they are periodically ejected from the disk and then fall back.This means that the only distinction between the two criteria is the moment in the velocity and the distance from the center at which they cool down.We also note that these particles are heated up and cool down many times during their motions.This multiple recycling of gas has also been observed in Suresh et al. (2019), finding that 60% of the gas coming from the disk undergoes at least two cycles.A similar behavior is found in the FIRE simulation (Anglés-Alcázar et al. 2017) for the different types of galaxies analyzed (a massive dwarf spheroidal, Milky-Way mass  systems, and an early-type galaxy): the recycled-to-ejected mass ratio does not show a clear dependence on halo mass and on redshift ( f rec ≃ 50 − 95%), while the average number of recycles is higher for low-mass systems (⟨N rec ⟩ ≃ 3 − 4 for Milky Way-like galaxies).
Particles cooling while inflowing show a more regular pattern, with many of them inflowing together as a larger stream, and hence showing a coherent position as a function of redshift r(z), and also cooling down at a similar moment (around z ≃ 0.4 for the particles shown).Similarly, also particles coming from satellites have a coherent behavior, being bonded by the gravitational well of the subhalo, and most of them are always cold during their motion.
Figure 7 shows the cumulative result of the origin tracking.From each chosen redshift, we track back the cold CGM gas until we identify the origin.The total mass of gas in each identified origin is shown in the left panel, while we show the percentage with respect to the total cold CGM gas in the right panel.
Article number, page 5 of 14 A&A proofs: manuscript no.paper When tracking back from low redshift (z ≃ 0.3), most of gas cools down directly in the CGM, with a little contribution from the disk or inflowing streams.The fraction of gas cooling in the CGM decreases at high redshift, making space for gas that is accreted cold from outside the galaxy.
An interpretation for this transition, occurring between z = 1.5 and z = 2, can be found considering the accretion history of the halo.In Fig. 7 we show the total hot (red shaded area) and cold (blue shaded area) inflowing gas, and the total outflowing gas (orange dashed line), in solar masses in the left panel and as percentage over the total inflowing plus outflowing gas in the right panel.At z ≃ 2.2, when the halo mass is around 3×10 11 M ⊙ , we have the transition between accretion via cold mode and hot mode, as discussed extensively in Kereš et al. (2005).Hence, when the galaxy is in the cold mode accretion phase, part of this accreted gas is in the CGM, and explains the gas classified as "cold inflow".On the other hand, at lower z most accreted gas is in the hot phase, so the main channel for cold gas in the CGM becomes gas cooling in situ ("cooling CGM", either "inflowing", "outflowing", or "static" depending on its radial velocity at the moment of the cooling).In addition, the peak of outflow activity is reached between z = 1 and z = 1.5, implying that there will be a large contribution of recycled gas to the CGM.At this redshift, we have an increase in "cooling CGM, outflowing" and gas coming cold from the disk (i.e., cold gas entrained in outflows).To be more quantitative, the majority of cold gas in the CGM at lower redshift (z < 2) cooled down in the CGM itself (60-80%), with a nonnegligible fraction coming from the disk (10-30%).
Gas metallicity can also be used to get insight into the cold CGM gas origin since gas from the disk is expected to be metal enriched, while gas populating the CGM directly from the IGM is expected to be pristine.In Fig. 8 we plot the PDF of gas metallicity for cold gas particles that have been classified as originating "from disk," "cooling CGM" (without distinguishing among outflowing, inflowing, and static), and "cold inflow."At both redshifts shown (z = 1 and z = 0.3), gas coming from the starforming disk has solar or super-solar metallicity, while gas coming from cold inflows is mainly metal-free.Only a small fraction of the "cold inflow" gas is enriched at z = 0.3, and we verified that this gas is located (or was recently located) close to the disk at r = 10 − 15 kpc, so it is not classified as "from disk" because of our geometrical definition of disk, but it still displays properties similar to the cold CGM gas coming from the disk.On the other hand, the gas classified as "cooling CGM" has a bimodal metallicity, showing that even metal-poor gas is able to cool down from a hot phase in the CGM.It is indeed possible for gas with metallicity as low as Z = 10 −3 to cool down from temperatures T = 10 5−6 K when the density is high enough (n > 10 −2 cm −3 ), as explored in more detail in Sect.3.5.We note that the Eris2k simulation includes metal diffusion due to turbulent mixing.However, as shown in (Shen et al. 2010) (e.g., their Figure 9), this model only acts on a subgrid scale, and does not appear to significantly alter large-scale mixing between pristine gas in accretion streams and gas outflowing from the disk.
To make a comparison with previous works on the origin of CGM gas, we have to consider that gas is usually tracked back without stopping the tracking when the particle is heated; in other words, the studies focus on understanding where particles come from, regardless of when it cools down.Tracking back farther in our simulation the particles labeled as "cooling CGM, static" from z ≃ 0.3, we find that ∼ 55% come from the disk (even if cooling down later in the CGM) and ∼ 33% come from r > R vir .In this way, we obtain that ∼ 65% of cold CGM gas at z ≃ 0.3 comes from the disk, and ∼ 30% is accreted from the IGM.The importance of galactic winds in feeding the CGM has already been assessed in previous works (Sokołowska et al. 2018), with a contribution spanning from 30% (Ford et al. 2014;Hafen et al. 2019) to 70% (Suresh et al. 2019).Nevertheless, the relative importance among different origins (IGM accretion, outflows, and satellites) appears to be strongly dependent on the halo mass and the redshift from which the gas is tracked back.
The contribution of different components to the cold CGM gas also reflects its dynamical properties.In particular, observational evidence of the contribution of galactic fountains to the cold CGM gas budget lies in the measured vertical gradient of the rotational velocity of gas in planes above and below the midplane of the galaxy (known as extra-planar gas).This has been observed in nearby spiral galaxies (e.g., Boomsma et al. 2004;Fraternali et al. 2004;Oosterloo et al. 2007;Voigtländer et al. 2013;Li et al. 2021), and is interpreted as a loss of angular mo-Article number, page 6 of 14 Decataldo D., Shen S., Mayer L., Baumschlager B., Madau P.: Origin of cold gas in the CGM Fig. 6.Random selection of 30 SPH particles chosen from those having the origin "cooling CGM, outflowing," "cooling CGM, inflowing,", "from disk," and "from satellites."The backtracking starts from z ≃ 0.3 Left column: Trajectory of the selected particles in a 2D plane, following them from z ≃ 0.5 to z ≃ 0.3.The dots become smaller and fainter going back in time, and the color is red and blue when the temperature is respectively higher and lower than 3 × 10 4 .Right column: Radial distance from the center of the halo of each particle as a function of redshift.mentum when gas ejected from the disk interacts with the hot CGM medium.In Fig. 9 we show rotation curves for the Eris2k disk (z plane = 0) and parallel planes at different distances from the disk (z plane = 2−6 kpc).Each plane is actually a slab of thickness 0.4 kpc, and the rotational velocity v rot is computed by aver-aging over 100 equal radial bins spanning the range 0 < r < 10 kpc.The gas retains a high rotational velocity up to 6 kpc, with peak velocities of 50-100 km/s.The rotational velocity globally decreases as we take planes that are more distant from the disk plane, with velocity gradients between -10 and -20 km s −1 kpc −1 , Fig. 7.For each snapshot at a given redshift z, the gas particles in the cold CGM are tracked back until their origin is identified (according to the sketch in Fig. 2).Left: Colored points at a given z representing the total mass of particles having a particular origin, as in the legend.The shaded blue and red areas show the total cold and hot inflowing gas mass, respectively.The orange dashed lines represents the total outflowing gass mass.Right: Same as the left panel, but the points represent the fraction of gas mass with respect to the initial mass of the cold CGM that is tracked back at a given redshift.In a similar fashion, the inflowing and outflowing gas are represented as a fraction.
Fig. 8. Normalized mass-weighted PDF of cold gas particle metallicities for gas classified as having origin "from disk," "cooling CGM," and "cold inflow" at redshift z = 1 (upper panel) and z = 0.3 (lower panel).depending on the radius at which we measure v rot .These values are in good agreement with the above-mentioned observations of extra-planar gas and other theoretical studies (Melioli et al. 2009;Marinacci et al. 2011;Armillotta et al. 2016).
A more detailed exploration of the disk-corona interaction is needed, and in principle it can be pursued by tracking SPH gas particles in the two components, following the evolution of their angular momentum.In a similar way, we suggest that the relative contribution of hot and cold mode accretion could lead to a correlation of the angular momentum of accreting gas with the angular momentum of the CGM.We defer to a future work for a detailed analysis of CMG dynamical properties.

Forward tracking of disk gas
Starting from a range of different redshifts, we track forward the cold gas in the disk to understand its fate.The forward tracking is stopped when the particle is heated, when the particle is ejected Fig. 10.For each snapshot at a given redshift z, gas particles that are cold and in the disk are tracked forward until an outcome is identified: the particle is heated within the disk, the particle is launched in an outflow but is still cold, the particle is launched in an outflow and is heated above T 0 , or the particle is converted into stars.If none of the previous outcomes is found earlier than z = 0.3, it means that the particle remains cold in the disk.Left: Total mass of the gas in each outcome.Right: Fraction of the gas mass in each outcome, with respect to the initial cold disk gas at the redshift z at which the forward tracking is started.
in an outflow (i.e., radial velocity v > 50 km/s, distinguishing between hot and cold outflowing gas), or when the gas is converted into a star particle.Particles that have not undergone any of these outcomes at z = 0.3 are labeled as "staying cold" in the disk (with only a negligible percentage going out from the disk with low velocities v < 50 km s −1 ).The result of this analysis is shown in Fig. 10; the total mass in each outcome is in the left panel, and the fraction with respect to the initial cold disk gas (at the redshift we start the tracking) is in the right panel.Starting the tracking from any redshift z > 1.2 gives approximately the same result: about 30% of the gas is outflowing in hot phase, ∼30% is outflowing in cold phase, and ∼25% is heated but stays in the disk; a smaller fraction (∼10%) remains cold and inside the disk, and a similar amount fuels star formation.A transition occurs at z ≃ 1.5, where the amount of hot gas remaining in the disk increases to ∼ 45%, and the outflowing gas decreases (both hot and cold) to 15-20%, while gas converted into stars also slowly decreases to a few percent.The reduction of gas in outflow was already highlighted in Fig. 7 (orange dashed line), and it is evident here how it is associated with a decrease in star formation activity.As a consequence, there is also more heated gas that is not able to escape the galaxy potential well.

Forward tracking of outflowing gas
Finally, we track forward outflowing particles, in the hot and cold phase, separately.In the first case, the tracking is stopped when a particle cools down (T < 3 × 10 4 K), while conversely we stop tracking a cold particle when it is heated up (T > 3 × 10 4 K) or converted to a star particle.Figure 11 shows the evolution of cold gas heated up (red line) and hot gas cooling down (blue line) starting the tracking from three different redshifts (z ≃ 2, 1.2, and 0.8).In all three cases, around 90% of the cold gas is heated up in the outflow, allowing 10% of gas to survive in cold phase.The timescale for cold gas survival (estimated as the time at which 50% of gas is heated) is around 300-400 Myr in all cases.On the other hand, the amount of gas cooling down in Fig. 11.Forward tracking of outflowing gas (defined as gas with outward velocity v > 50 km s −1 ), starting from redshift z ≃ 1.2 (solid line), z ≃ 1.2 (dashed line), and z ≃ 0.8 (dotted line).The red line shows the cumulative fraction of cold outflowing gas that is heated, while the blue line represents the hot outflowing gas that cools down.hot outflows depends strongly on the initial redshift at which the outflow is launched.In particular, only a small fraction (∼ 10%) of outflows launched at redshift z ≃ 2 cool down; this fraction increases to 40% when the outflow is launched at z ≤ 1.2.Checking the average properties of hot outflowing gas, we find that the average distance from the center changes from ⟨r⟩ z≃2 ≃ 150 kpc to ⟨r⟩ z≃0.8 ≃ 75 kpc (with a corresponding drop in average velocity from ⟨v⟩ z≃2 ≃ 100 km s −1 to ⟨v⟩ z≃0.8 ≃ 80 km s −1 ).Because it is closer to the central region, gas that is launched later is more dense, and it cannot reach greater distances due to lower veloc-Fig.12.For each snapshot at a given redshift z, outflowing gas particles are tracked forward in time.For cold gas particles (solid lines), the tracking is stopped when they are heated above T = 3 × 10 4 K ("heated up") or form stars ("into stars"), while particles not undergoing either outcome at z = 0.3 are labeled "staying cold."Similarly, the tracking of hot particles (dashed lines) is stopped if they cool down or form stars, so that the possible outcomes are "cooling down," "into stars," or otherwise "staying hot" (although a hot gas particle needs to become cold before being eligible for star formation).Left: Total mass of the gas in each outcome.Right: Fraction of the gas mass in each outcome, with respect to the initial cold outflow mass for heated up, into stars, and staying cold outcomes, and with respect to the initial hot outflow mass for cooling down and staying hot outcomes.Fig. 13.Hot outflowing gas is tracked forward from z ≃ 0.5 until it cools down.The two panels show the mass-weighted particle distribution according to their initial radius (r init ) and the radius at which they cool down (r cooling ) in the left panel, and according to their initial radial velocity (v init ) and the velocity at which they cool down (v cooling ) in the right panel.The two panels show the distribution of particles, weighted by their mass, according to their initial radius (x-axis), and the radius at which they cool down (y-axis).
ity; as a consequence, the cooling time is smaller, and this gas is favored to cool down.
A cumulative analysis of the fate of the cold and hot phases of outflowing gas is shown in Fig. 12.The solid lines show the cold gas in outflow heated up, remaining cold, and converted into stars.In the right panel we can see that less than 20% of the cold gas can survive in cold phase, in agreement with what was already observed in Fig. 11, while 5-10% of cold gas is forming stars.On the other hand, the dashed lines show the hot gas cooling down and the gas remaining hot.In this case, we see that only ∼10% of gas cools down at high redshift, while this fraction increases to almost 50% for z < 0.75.As noted in the previous paragraph, this is due to the outflowing gas being launched from a region closer to the center, which is denser and hence favors gas cooling.Nevertheless, we note that radiation feedback from stellar sources is not included in this simulation, and some studies have shown the importance of this mechanism in heating the cold gas in outflows (Decataldo et al. 2017).Hence, our study could overestimate the amount of the cold gas surviving in outflow, and the amount of gas cooling down from the hot phase.
To inspect further the dynamics of cooling gas in outflow, we analyzed the position and velocity of the tracked gas at the moment of cooling (Fig. 13).The left panel shows the distribution of hot outflowing SPH particles according to the radius at the beginning of the tracking (i.e., initial distance from the center) and the radius at which it cools down.Analogously, the right plot Fig. 14.Cold outflowing gas is tracked forward from z ≃ 0.5 until it is heated up.The plot shows the distribution of particles, weighted by their mass, according to their initial radius at z ≃ 0.5 (x-axis) and the radius at which they are heated (y-axis).
shows the initial velocity and the velocity at the moment of cooling.The particles represented here are tracked forward starting from z = 0.5, so r init and v init are not necessarily the position and velocity at the moment of the launch of the outflow; the particle could already in outflow since a previous snapshot at z > 0.5.r init spans a range from 3 kpc to ∼100 kpc, and most particles cool down immediately or after travelling at most ∼ 30 kpc.For some particles we have r cooling < r init , meaning that they have reached the turnover point and they are falling back toward the center of the halo.This is also clear from the velocity distribution: most particles cool down at at v cooling ≃ 0 (i.e., at the turning point) and the distribution is almost symmetrical between particles cooling down while still outflowing or falling back to the center.This configuration points to a galactic fountain, which has already been shown to be one of the main contributors to the total CGM cold gas budget in previous works (Suresh et al. 2019).We also note that particles with higher r init travel only a short distance before cooling because their velocity is in general lower, either because they are far from the central region where stellar feedback is more effective, or they are already close to the turnover point of an outflow launched at an earlier time.
Similarly, in Fig. 14 we plot the distribution of cold outflowing particles relative to the radius at the beginning of the tracking and to the radius at which they are heated.In this case most particles start their journey at a distance r < 30 kpc, where most of cold CGM gas is located (see also Fig. 3).These particles are then heated after the turnover at radius r heated < 10 kpc (i.e., within the disk), suggesting stellar feedback as the main heating mechanism.

The cooling of gas in the CGM
In our analysis of the origin of cold gas in the CGM (Sect.3.2), we found that an important contribution, especially at z < 2 is from gas cooling down from higher temperatures.This means that detachment from the hot phase and the cooling down is a process that can happen at CGM densities and on timescales comparable to the dynamic timescales of the galaxy.In order to asses this point, we compute the temperature as a function of time for gas with different densities and initial temperatures, adopting the same cooling function used in the Eris2k simulation and assuming constant density during the cooling.The results are shown in Fig. 15, for gas at solar metallicity (upper panel) and 10 −3 Z ⊙ .The different colors stand for the different gas densities, from n = 10 −3 to 10 3 cm −3 , while the line styles represent the different initial temperatures T = 10 5 , 10 6 , and 10 7 K. Basically, all gas with density n ≥ 10 −2 cm −3 (which is the average density for gas at a distance r = 10 − 30 kpc at z ≃ 0.3) can cool down in less than 50 Myr at Z = Z ⊙ , and in less than 500 Myr for the low-Z case, which allows the gas in the CGM to cool down spontaneously on a short timescale compared to the galaxy evolution.In addition, gas in outflows typically has higher densities (1 − 100 cm −3 in the simulation), and consequently with cooling timescales on the order of 10 Myr, even at the highest initial temperatures (T ∼ 10 7 K).
The condensation of cold clumps in the warm-hot intergalactic medium (WHIM) was studied analytically by Binney et al. (2009) (who had doubts about the efficiency of this process); they showed via linear perturbation analysis that thermal instabilities in a hot stratified corona are suppressed by buoyancy and thermal conduction.Nevertheless, this mechanism regained favor thanks to studies taking into account the rotation of the corona (Nipoti 2010;Sormani & Sobacchi 2019;Sobacchi & Sormani 2019) and studies that explored the nonlinear nature of this process (McCourt et al. 2012;Sharma et al. 2012).In the latter, the possibility of gas cooling and fragmentation within the CGM appears to be related to the ratio t cool /t ff , where t ff is the free-fall time into the potential well of the galaxy.In particular, simulations in these studies find a value t cool /t ff < 1 − 10 to be a necessary condition for condensation.The precipitationregulated feedback model (Sharma et al. 2012;Voit et al. 2017) manages to explain why galaxies self-regulate to this value of t cool /t ff : outflows promote condensation of adiabatically uplifted gas, which eventually falls back and fuels star formation (and eventually stellar feedback), until the outflows become stronger and the gas is heated too much to allow condensation.However, we mention that Saeedzadeh et al. (2023) highlighted how the t cool /t ff criterion does not always capture the susceptibility of gas to condensation: high-t cool /t ff gas still cools down in the presence of large-scale perturbations (Choudhury et al. 2019) and low-t cool /t ff gas is unable to condense when the ratio decreases only for a short period of time (e.g., when t ff increases suddenly in the wake of a satellite).
We explore the thermal instability of CGM gas in Eris2k by plotting histograms of t cool /t ff for hot (T > 3×10 4 K) gas located at different radial distances from the halo center at z = 0.3 and z = 1 (left panel of Fig. 16).Most gas in the CGM (hence r > 10 kpc, which is the typical disk radius assumed in this work) up to 50 kpc has t cool /t ff ≪ 1, making the precipitation scenario a feasible channel for the origin of cold gas.This region close to the disk is where most of cold gas clumps are located (see Fig. 3), and where the distribution of cooling gas in outflow peaks (see red lines in the left panel of Fig. 13).Gas at greater distances (r > 50 kpc) show a t cool /t ff distribution peaking around ∼ 5, but still with a large fraction of gas with t cool /t ff < 1, especially at higher z.In general, we observe that the amount of unstable gas (t cool /t ff < 1) at high radii increases with redshift.To better show the physical conditions of unstable gas, in the right panel of Fig. 16, we show the 2D n-T phase diagram of the hot CGM (r > 10 kpc) at z = 0.3, where the gas in each bin has been color-coded according to the average value of t cool /t ff .As expected already from the analysis in Fig. 15, high-density gas is unstable because of the short cooling timescale.CGM gas at lower densities (n < Fig. 15.Temperature evolution of gas with different densities (from n = 10 −2 cm −3 to n = 10 3 cm −3 ), starting from different initial temperatures (T = 10 5 , 10 6 , and 10 7 K).The gas has metallicity Z = Z ⊙ in the upper panel, and Z = 10 −3 Z ⊙ in the lower panel.The gray shaded area marks the temperature T = 3 × 10 4 K, which in our work is defined as the threshold temperature between cold and hot gas.Fig. 16.Properties of the hot (T > 3 × 10 4 K) CGM gas according to its t cool /t ff ratio.Left: Distribution of t cool /t ff at two different redshifts, z = 0.3 (solid line) and z = 1 (dashed line).Different colors indicate gas belonging to different radial distance bins, while the shaded regions highlight the t cool /t ff < 1 criterion for precipitation-regulated feedback, and the less stringent t cool /t ff < 10 condition.Right: Distribution of hot CGM gas at z = 0.3 in the n-T phase diagram, color-coded by the average value of t cool /t ff .The contour lines show the t cool /t ff = 1 and 10 boundaries.
10 −3 cm −3 ) might still be unstable gas when the temperature is not too high (T < 10 6 K).
We note that the value of the ratio t cool /t ff for the gas depends on the maximum resolution reached in the simulation since the cooling time is shorter when high-density gas is resolved (t cool ∼ n −1 ).Hence, at higher resolution it is possible that the distribution of t cool /t ff in Fig. 16 would move farther to the left, until convergence, implying even more thermally unstable gas.

Conclusions
In the present work we analyzed the history of the gas in the main halo of the simulation Eris2k, by tracking SPH particles properties through the evolution of the halo.In particular, we focused on cold gas (T < 3 × 10 4 K) in the CGM (10 kpc < r < 1.5 R vir ), making up around 15% of the total CGM gas mass at redshisft z ≃ 0.3, concentrated mainly in small clumps in the region within 50 kpc.
We tracked back particles from cold gas in the CGM, starting from different evolutionary stages of the halo.At high redshift (z > 2.5), the CGM is dominated by gas that flows into the halo via cold accretion streams (cold mode accretion), so that the majority of cold CGM gas is accreted directly from the IGM.At lower redshifts (z < 2), the halo is mainly accreating via the hot mode, and the CGM is dominated by shock-heated and recycled gas.In this regime, we show that most cold CGM gas has cooled down from the hot phase (∼ 40% of the cold CGM gas is approximately static at the moment of cooling, up to 30% is in inflow).Outflows also contribute significantly to the cold gas budget, especially around z ≃ 1.5, which is the peak of outflow activity; cold CGM gas is indeed enriched by gas cooling within the outflow (10-20%), and cold gas lifted from the disk (up to 30% at z ≃ 1.5).We also find that gas can undergo multiple cooling and heating phases, especially particles that are recycled by falling back into the disk after being in outflow.
To inspect further the thermal evolution of the gas in the galaxy, we tracked separately particles in the different components of the galaxy: -Disk: For z < 1.2 most of gas is expelled in outflows both in cold and hot form (∼ 30% each), while later more than 60% of the gas remains in the disk, even when heated, due to high star formation and a decreased outflow activity; -Cold outflows: Starting the tracking from different z, we always find that only ∼ 10% of the gas is able to survive in cold form until z ∼ 0.3; -Hot outflow: a higher fraction of gas launched at lower z tends to cool down (∼ 40% from z ≃ 0.8) since the central region is more compact and the outflowing gas is initially denser.
We note that the amount of cold gas in outflows could be overestimated since the simulation does not include radiative transfer on the fly, and cold clumps are subject to the intense radiation field of stars in the disk.
We obtained an insight into the scenario of gas cooling directly in the CGM by computing the cooling time of gas at different initial temperatures and density values.Gas with density n ≥ 10 −2 cm −3 is able to cool down in less than 50 Myr, which is compatible with the spontaneous cooling that we observe in the simulation.In addition, gas at a distance r < 50 kpc (where most of the cold CGM gas is located) satisfies the condition t cool /t ff < 1−10, which allows for gas fragmentation and cooling in the CGM (precipitation-regulated feedback model).
Our analysis was performed on a zoom-in simulation, which allowed for the high spatial resolution (maximum resolution of ∼ 5 pc in the CGM) needed to resolve the cold phase, which would otherwise be mixed with the hot phase.The drawback of this approach is that it allowed us to study only one sample, specifically a Milky Way-mass galaxy.We expect some differences in systems with different masses, as pointed out by similar analysis on large cosmological boxes (Hafen et al. 2019;Nelson et al. 2020); for example, the contribution of winds to the CGM gas budget is seen to decrease with halo mass, while the amount of gas stripped from satellites increases (Hafen et al. 2019).Hence, to improve the current analysis on the origin of cold gas, the following steps are needed: First, galaxies with different properties should be analyzed (mass, star formation rate), and second, the effect of resolution on the amount of cold gas should be investigated by comparison with simulations on smaller scales.

Fig. 1 .
Fig. 1.Gas resolution (defined as the smoothing length of SPH particles) in the simulation Eris2k as a function of distance from the main halo center.The blue curve represents the median value in the radial bin, while the shaded area includes the variation from the 10th to the 90th percentile.The other two lines show the median resolution for a halo in Illustris TNG50 (M halo ≃ 10 12.5 M ⊙ ) and the high-resolution simulation in Suresh et al. (2019) (with 10th to 90th percentile variation).

Fig. 3 .
Fig. 3. Overview of Eris2k gas temperature at z ≃ 0.3.Upper left: Temperature map on a slice of gas that passes through the halo center.Small green patches of gas with temperatures lower than 10 5 K are visible, ubiquitous, and present especially in the central region within 100 kpc.Upper right: Only gas with temperature T < 3 × 10 4 K is shown from the same slice.Most of the cold gas is located in the central area, and large bubbles are also present at greater distances.Lower left: Distribution of the gas in the whole box in the density-temperature diagram; the red line marks the threshold temperature T = 3 × 10 4 K. Lower right: Distribution in the radial distance-temperature diagram.

Fig. 4 .
Fig. 4. Column densities of HI (left) and MgII (right) at redshift z ≃ 0.3, for a box of size 600 3 kpc 3 around the main halo in Eris2k.The white line shows the virial radius of the halo.

Fig. 5 .
Fig. 5. Radial profiles of HI (left) and MgII (right) column densities, averaging over annuli of thickness 0.5 kpc.The solid line represents the median value, and the shaded region includes the variation from the 10th to the 90th percentile.The observational data from the COS-Halos survey are overplotted, distinguishing between star-forming (blue squares) and non-star-forming (red diamonds) galaxies.The arrows indicate upper and lower limits.

Fig. 9 .
Fig. 9. Rotation curves in the galaxy disk and in parallel planes at different distances z plane .The rotational velocity is computed by defining 100 equal radial bins and averaging particle rotational velocities inside the corresponding annuli.