High Angular Momentum Halo Gas: A Feedback and Code-independent Prediction of LCDM

In the standard Lambda Cold Dark Matter (LCDM) picture of galaxy formation, gas accreting onto a growing dark matter halo shock-heats to the virial temperature of the halo, giving the gas time to virialize and eventually cool out of the hot gaseous halo and sink into the central galaxy (Rees & Ostriker 1977; Silk 1977; White & Rees 1978; White & Frenk 1991; Maller & Bullock 2004). Under this picture of galaxy growth, it is expected that the resulting angular momentum distribution of galaxies should mimic the spin of their dark matter, resulting in rotationally supported galaxy disks (and presumably hot gaseous halos as well) that are proportional to the spin of the dark matter halo (Fall & Efstathiou 1980; Mo et al. 1998), which has been well-studied in dissipationless N-body simulations and semi-analytic merger trees (e.g., Bullock et al. 2001; Maller et al. 2002; Vitvitska et al. 2002; Avila-Reese et al. 2005; D'Onghia & Navarro 2007; Bett et al. 2010; Muñoz-Cuartas et al. 2011; Ishiyama et al. 2013; Trowland et al. 2013; Kim et al. 2015).

However, recent advances in hydrodynamic simulations and galaxy formation theory have increasingly emphasized the importance of "cold flows"—gas accretion onto galaxy halos via filamentary streams with cooling times shorter than the compression time for establishing a stable shock¹³ , either when the halo is below a critical mass threshold, or even for massive halos at sufficiently high redshift (e.g., Binney 1977; Kereš et al. 2005; Dekel & Birnboim 2006; Ocvirk et al. 2008; Brooks et al. 2009; Dekel et al. 2009; Faucher-Giguère & Kereš 2011; Faucher-Giguère et al. 2011; Stewart et al. 2011a; van de Voort et al. 2011, 2015; Hobbs et al. 2015). In the cold flow paradigm, gas that is accreted in the cold mode tends to have specific angular momentum considerably higher than the dark matter (Chen et al. 2003; Sharma & Steinmetz 2005; Agertz et al. 2009; Kereš et al. 2009; Kereš & Hernquist 2009; Brook et al. 2011; Kimm et al. 2011; Stewart et al. 2011b), inconsistent with the previous picture of galaxy angular momentum buildup. The resulting angular momentum of the stellar disk may be rather different from that of the accreted gas because of feedback effects (Maller & Dekel 2002; Brook et al. 2011).

As a result of this changing paradigm for cosmological gas accretion and galaxy growth (for a recent review, see Stewart 2016), a new scenario of angular momentum acquisition in galaxies and galaxy halos seems to be emerging. In this picture (Kimm et al. 2011; Pichon et al. 2011; Stewart et al. 2011b, 2013; Codis et al. 2012, 2015; Danovich et al. 2012, 2015; Prieto et al. 2015; Tillson et al. 2015), the particularly high angular momentum of cold flow gas is related to its coherent, filamentary origin, coupled with the specific geometry of the cosmic web in the environment of a given galaxy. These filamentary cold flows deliver significant angular momentum to galaxy halos, with the cold gas orbiting for ∼1−2 dynamical times before spiraling into the central galaxy. At any given time, galaxy disks typically have lower spin than halo gas, owing to the fact that the specific angular momentum of infalling material increases with time. Halo gas is "younger" and this correlates with higher spin.

Importantly, this scenario is predictive. The high-spin halo gas is often (but not always) coherent in its spin direction, with inspiraling cold streams often forming a thick planar structure of high angular momentum cool gas that co-rotates with the central disk. It is important to emphasize that while this extended gas tends to rotate, it is not angular-momentum supported. Rather, this gas usually spirals in on ∼2 dynamical times. Though not perfectly aligned with the orientation of the galactic disk, the inspiraling halo gas usually has coherent rotation along a preferred plane.

Encouragingly, an increasing number of observations have begun to demonstrate the abundance of high angular momentum material in galaxy halos, qualitatively consistent with this emerging theoretical picture. In the local universe, some of these observations include detection of high angular momentum extended H i disks and XUV disks (Oosterloo et al. 2007; Christlein & Zaritsky 2008; Sancisi et al. 2008; Lemonias et al. 2011; Holwerda et al. 2012), as well as low metallicity high angular momentum gas (presumably from fresh accretion) in polar ring galaxies (Spavone et al. 2010). There is even indication that local extended H i disks may be environmentally dependent on the galaxy's filamentary environment (Courtois et al. 2015). At moderate redshift (z ∼ 0.5–1.5), numerous absorption line studies of the circumgalactic medium of galaxies have begun to emphasize the bimodal properties of absorbers, where absorption along the galaxy's major axis tends to show high angular momentum (co-rotating) inflow, and absorption along the galaxy's minor axis shows observational signatures of outflow (Kacprzak et al. 2010, 2012a, 2012b; Bouché et al. 2012, 2013; Crighton et al. 2013; Nielsen et al. 2015; Bouché et al. 2016; Diamond-Stanic et al. 2016; Ho et al. 2017). At higher redshift (z ∼ 2–3), kinematic studies of Lyα "blobs" have observed large-scale rotation consistent with high angular momentum cold gas accretion (Martin et al. 2014; Prescott et al. 2015). There are also recent detections of massive protogalactic gaseous disks kinematically linked to gas inflow along cosmic filaments, strikingly similar to the theoretical "cold flow disk" structure (Martin et al. 2015, 2016).

In this context, it is important that we ascertain how robust the predictions of these cosmological simulations are—a difficult task, considering that many properties of simulated galaxies depend sensitively on the implementation of uncertain subgrid physics models such as gas cooling, star formation, radiation pressure, and supernova feedback (e.g., Thacker & Couchman 2000; Kay et al. 2002; Gnedin et al. 2011; Piontek & Steinmetz 2011; Martizzi et al. 2012; Scannapieco et al. 2012; Agertz et al. 2013; Vogelsberger et al. 2013, 2014; Ceverino et al. 2014; Genel et al. 2015; Marasco et al. 2015; Agertz & Kravtsov 2016). In addition, even with identical subgrid implementations, there are inherent numerical advantages and disadvantages between different hydrodynamic code implementations—for example, Lagrangian smoothed particle hydrodynamic (SPH) versus Eulerian grid codes—that result in artificial differences between galaxies simulated with different codes (e.g., Morris 1996; Agertz et al. 2007; Wadsley et al. 2008; Cullen & Dehnen 2010; Hahn et al. 2010; Springel 2010; Hopkins 2015; Richardson et al. 2016).

In order to test the validity of the emerging cold flow picture of angular momentum acquisition, we must ascertain the dependency of these predictions on the use of different numerical techniques and a variety of cutting-edge subgrid physics models. In this paper, we run five hydrodynamic zoom-in simulations of a Milky Way-sized galaxy, each with identical cosmological initial conditions but with different codes: Enzo, Ramses, Art, Arepo, and Gizmo-PSPH, each implemented with recent subgrid physics models. In order to ensure uniform analysis for different hydrodynamic codes, we utilize the analysis software yt (which allows a single analysis routine to be run on different code architectures; Turk et al. 2011) to explore the angular momentum content of halo gas and whether or not the expected "cold flow disk" prediction is robust across these disparate platforms. We introduce the simulations in Section 2, present our main results from the comparison Sections 3–5, finding that the same qualitative picture of high angular momentum halo gas in the form of co-directional inspiraling cold streams (which do occasionally take the form of cold flow disks) is present in all simulations—a seemingly natural consequence of filamentary gas accretion in LCDM. We discuss the implications of these results and the growing observational evidence of their existence in Section 6 and summarize and conclude in Section 7.

2.1. Overview

The simulations used in this paper are all part of the Scylla Multi-Code Comparison Project. This project resimulates a Milky Way halo mass zoom-in simulation (originally performed by Ryan Joung with Enzo in Joung et al. 2012) using other cosmological hydrodynamic codes. While we focus on the redshift range 1 < z < 3 in this work, we note that the resulting disk-type galaxy and its gaseous halo have already been studied in detail at low redshift (e.g., Fernández et al. 2012; Joung et al. 2012). Of particular importance to this work, Fernández et al. (2012) determined that the mass (in H i), covering fraction, and spatial distribution of the cold gas halo at z = 0 are consistent with existing observations of nearby spiral galaxies.

The codes are all run with their recent¹⁴ subgrid models in order to compare state of the art simulations across codes. Thus, the project is much like the Aquilla code comparison (Scannapieco et al. 2012) but with higher resolution. Our resolution is similar to the Agora code comparison project (Kim et al. 2014), but that project is seeking to use uniform physics while we are running each code as it has been used for other science papers. The codes used here are Enzo, Art, Ramses, Arepo, and Gizmo-PSPH. For all runs, the cosmology, dark matter particle mass, and box size are identical—the box is 25 Mpc h⁻¹ across, with a much smaller region simulated at high resolution, using dark matter particles of mass 1.75 × 10⁵ M_⊙. All adaptive mesh refinement (AMR) codes reach the same maximum refinement of 95 h⁻¹ comoving pc, which is identical to the force resolution of the Lagrangian codes, with the exception of gas particles in Gizmo-PSPH, which uses an adaptive gravitational softening with a minimum value of 14 h⁻¹ comoving pc. A flat cosmology consistent with WMAP5 (Komatsu et al. 2009) is used throughout, with Ω_m = 0.279, Ω_Λ = 0.721, Ω_b = 0.046, h = 0.70, σ₈ = 0.82, and n_s = 0.96.

The initial conditions for the original Enzo run (Joung et al. 2012) were generated with the code Grafic ¹⁵ (Bertschinger 2011) with a starting redshift of z = 99. The same code with the same seed was used to generate the initial conditions for the Ramses run. For all other runs, the dark matter particles from the Enzo run were used to determine the initial conditions. That is, the dark matter particles were set identical to those in the Enzo run and baryons were added based on the dark matter distribution (no separate transfer function). We expect these differences to be negligible by the redshift where galaxies are forming.

Based on the cosmological model specified above, all Lagrangian codes set the gas mass resolution (m_gas =3.3 × 10⁴ M_⊙) relative to the dark matter particle mass (m_DM = 1.75 × 10⁵ M_⊙). Table 1 outlines many of the pertinent details for each code, including star formation (SF) density thresholds and efficiency parameters, epoch of reionization, UV background model, and the type of stellar feedback model adopted. Below, we describe the gas cooling and feedback physics of each individual run in more detail and include references to recently published science papers that utilize similar subgrid physics models as those implemented here.

Table 1.  Simulation Code Details

	Enzo	Art	Ramses	Arepo	Gizmo-PSPH
Gravity Solver	FFT in the root grid	Multilevel particle mesh	Multigrid particle mesh	Tree multipole expansion particle mesh	Tree multipole expansion particle mesh
Hydrodynamics Solver	Third-order piecewise parabolic method	Second-order Godunov method	Second-order MUSCL scheme	Second-order MUSCL schemea	Pressure–energy SPH
High Res. mDMb	1.75 × 105 M⊙	1.75 × 105 M⊙	1.75 × 105 M⊙	1.75 × 105 M⊙	1.75 × 105 M⊙
Grav. Softening [h−1 comoving pc]	95 (DM, gas)	95 (DM, gas)	95 (DM, gas)	95 (DM, gas)	95 (DM), 14 (gas)
SF Threshold	0.04 cm−3	1 cm−3	1 cm−3	0.13 cm−3	5 cm−3 + self-grav. + molecular
SF Efficiency	$\epsilon =0.03$	$\epsilon =0.03$	$\epsilon =0.03$	tSFR = 2.2 Gyr	$\epsilon =1$ (in self-grav., molecular gas)
Stellar Feedbackc	Thermal	Thermal & Rad.	Kinetic	Kinetic	Mixed [see text]
Temperature Floor	10 K	300 K	100 K	500 K	10 K
UV Background	HM96 (increased Gaussian width)	HM96	HM96	FG09	FG09
Reionization	z = 6	z = 7	z = 10	z = 10	z = 10

Notes. HM96—Haardt & Madau (1996); FG09—Faucher-Giguère et al. (2009).

^aSubsequent versions of Arepo have switched to a different time integration (Pakmor et al. 2016) using Heun's method. ^bFor Lagrangian codes, high resolution gas particle mass is 3.3 × 10⁴ M_⊙. ^cSee the text for detailed descriptions of feedback models.

Download table as:  ASCIITypeset image

For all analyses that follow, we make the distinction between "cold" and "hot" gas by a temperature cutoff of 250,000 K (commonly used as the distinction between "cold-mode" and "hot-mode" gas accretion; e.g., Kereš et al. 2005, 2009; Stewart et al. 2011b, 2013). Using the Rockstar halo finder (Behroozi et al. 2013b), we calculate the virial radius of each simulation at each output redshift, finding that the halo-finding algorithms produce slightly different virial radii at the same redshift for different simulations (±5% from the mean). For the sake of identical comparison between codes, we therefore utilize a fitting function for R_vir(z) that averages over all simulations to adopt an identical estimated virial radius for all simulations (at a given redshift), which varies from 68 physical kpc at z = 3 to 171 physical kpc at z = 1.

In order to guarantee uniform analysis for the varied code architectures and file formats, all analyses presented here have been performed utilizing the yt ¹⁶ analysis software (Turk et al. 2011; Turk & Smith 2011; Turk 2013), an open source project that has been developed and is continually being maintained and improved by the astrophysical community for the intended purpose of supporting cross-code compatible hydrodynamic analysis routines. In plots showing images of the gas distribution, we use a slightly older version of yt (version yt-3.2.3) because after that version yt updated the way gas particles are deposited into cells. Prior to the update, yt used clouds-in-cells deposition to determine the gas properties of a cell, while after the update the sph smoothing kernel is used instead. Although this update gives more accurate deposition for sph particles, it is very inaccurate for Arepo where the particle size relates to the volume of the cells from the Voronoi tessellation and not a smoothing length. The clouds-in-cells deposition gives adequate results for both methods so we use it for images. For quantitative analysis, we use the gas particles for both Arepo and Gizmo-PSPH so that no deposition into cells is required.

2.2. Enzo

2.3. Art

2.4. Ramses

2.5. Arepo

2.6. Gizmo-PSPH

3.1. Large-scale Structure and Mass Growth

We begin with a visual inspection of the region around the galaxy for each simulation. Figure 1 shows the gas density (number density of H; top) and density-weighted temperature (bottom) projections at z = 3 through a cube of width 272 physical kpc (4R_vir at this redshift). The top panels of Figure 1 show qualitative agreement between the simulations on the general geometry and structure of the forming disk galaxy at this redshift, as well as its placement in a large-scale cosmic filament that is continually delivering an inflow of cold gas into the virial radius of the galaxy. However, the detailed structure of the galaxy—and even that of the filament into which the galaxy is embedded—does appear to vary significantly between simulations. For example, the width of the cosmic filament, the size and structure of the galactic disk, and the peak density of infalling satellite galaxies all vary on a noticeable level.

Zoom In Zoom Out Reset image size

Perhaps more striking is the temperature differences among simulations shown in the bottom panels of Figure 1. All simulations demonstrate the presence of a significant gaseous halo around the galaxy, as well as streams of filamentary gas that penetrate the halo and deposit cold gas in the inner halo, near the galactic region. However, the extent that feedback has enriched the CGM and IGM, the density structure of the gaseous halo, the temperature distribution of hot gas, and the precise structure of the cold flows as they interact with the gaseous halo of the galaxy vary significantly among simulations.

To illustrate some of the similarities and differences among the simulations, the two left panels of Figure 2 show the mass growth of the halo as a function of time, including the total mass (black) and galaxy stellar mass (yellow) on the leftmost panel, as well as the cold (blue) and hot (magenta) gas fractions within the virial radius. Note that the total virial mass (left panel) is quite similar among simulations, despite very different feedback implementations. Comparing the linear scale of the middle panel to the log scale of the left panel, we also note that the gas fractions (middle panel) are relatively similar among the simulations, although there are still noticeable variations. The average total/cold/hot gas fractions during the entire redshift range $z=3\mbox{--}1$ averaged over all the simulations is represented by the horizontal dotted lines in the figure. Not surprisingly (given the mass scale of this halo), all simulations show that the dominant supply of halo gas is in a cold phase, rather than a massive reservoir of hot gas.

Zoom In Zoom Out Reset image size

The total galaxy stellar mass (left panel) shows a much more significant variation among simulations, with the Gizmo-PSPH code in particular forming a much smaller stellar mass than any of the other codes used here, likely as a result of strong feedback implementations. As a comparison, the upper limit of baryonic mass (the virial mass times the cosmic baryon fraction) is shown here as a thin black dotted line, and the shaded cyan region shows the expected galaxy stellar mass range for the given virial mass based on abundance matching (median value ±1σ from Behroozi et al. 2013a). Interestingly, the galaxy stellar mass for Gizmo-PSPH is much closer to observational expectations, and may even be slightly underproducing stars in the simulated galaxy, rather than overproducing them, as in the other simulations. As our goal in this work is to focus on similarities between codes, with emphasis on the galaxy halos and not the galaxies themselves, we defer a more detailed discussion of the numerous differences between the simulations and their implications for galaxy formation as a topic for future study.

3.2. Angular Momentum

In order to place this high angular momentum cold halo gas in the proper cosmological context, Figure 3 shows the large-scale environment around the simulated galaxy—where, for purposes of this work, we define the halo environment by box widths of 8R_vir (544 physical kpc at z = 3). The top panels again show the gas density (H number density), similar to Figure 1 but zoomed out by a factor of two and viewed along an orthogonal orientation. The bottom panels show the line-of-sight velocity of all cold gas above a minimum density threshold in hydrogen (all forms) of n_H > 3 × 10⁻⁴ cm⁻³, which was chosen to select only gas sufficiently dense to be embedded in filamentary (or dark matter halo) structures in these large-scale environments.

Zoom In Zoom Out Reset image size

Because the galaxy is the most massive halo in its environment (i.e., not a member of a group or cluster), the cosmic filaments in its environment are strongly affected by the halo potential, with gas, dark matter, and smaller galaxies all flowing along the filaments toward the galaxy, demonstrated by the clear line-of-sight velocity indications in the bottom panels. For example, the filament to the upper left of the galaxy (situated in front of the galaxy along this line of sight) consistently shows redshifted velocities in all simulations, while the two filaments below (and behind) the galaxy are consistently blueshifted.

This result is perhaps not surprising, as any three-dimensional filamentary structure where matter flows along the cosmic web toward a central overdensity (and is viewed along an arbitrary axis) is unlikely to show multiple filaments all flowing perpendicular to the line of sight. Thus, one should naively expect strong line-of-sight velocities to be apparent when viewing large-scale filamentary gas flows. Although this may not be a surprising result, it is important to keep these large-scale gas flows in mind for future discussion of the kinematics of inspiraling cold streams in the galaxy's halo. We will see in Section 5 that these large-scale filamentary flows have a direct impact on the behavior of the cold gas within the virial radius of the halo.

Note that the line-of-sight velocity structure of the filament flowing in from the right of the galaxy shows considerably more variation between the simulations. This occurs because this filament does happen to be roughly perpendicular to the line of sight. Thus, the velocities along the rightmost filament are more sensitive to the peculiar velocities of galaxies, gas streams, and outflows, which vary more strongly between simulations than the gross large-scale flows toward the central halo.

We begin investigating the morphology and kinematics of cold halo gas in Figure 4, which is analogous to Figure 3, except that the panels now focus only on material within the virial radius (panel widths of 136 physical kpc at z = 3). The bottom panels again show line-of-sight velocity maps of cold dense gas¹⁸ in the halo, except that we have increased the minimum density threshold by a factor of 10 when compared to Figure 3, to a hydrogen density of n_H > 3 × 10⁻³ cm⁻³ (this should correspond to a minimum hydrogen column density of N_{H i} ≳ 10¹⁷ cm⁻²; Schaye 2001; Altay et al. 2011).

Zoom In Zoom Out Reset image size

The exact morphology of gas in the halo of the galaxy varies considerably among the simulations, which is not surprising, given the vastly different feedback mechanisms implemented in each simulation, some of which drive explosive spherical outflows that violently shred the ISM and CGM of the galaxy (e.g., Gizmo-PSPH) and some of which instead drive high-velocity bi-conical outflows out of the plane of the galaxy (e.g., Enzo). However, we also note that some of the morphological differences may also be influenced by the precise timing of galaxy mergers. For example, the Gizmo-PSPH simulation is in the midst of a violent outflow at this epoch, due to a recent merger, which partially explains the significantly more chaotic structure shown in Figure 4 (we will demonstrate in Section 5.2 that the Gizmo-PSPH simulation's line-of-sight velocity structure is much more similar to the other simulations immediately before and after this merger-driven outflow event). While the same general merger and accretion history takes place for each simulation, the exact timing of these mergers at a given epoch may vary, and any coherent velocity structure for cold gas in the galaxy's halo is typically destroyed during a sufficiently strong outflow event.

Despite these varied differences in morphology, the middle panels of Figure 4 show a similar qualitative picture. As was the case with the large-scale environment, the cold gas entering the virial radius from the upper-left filament shows a dramatic redshift in each simulation, while the cold gas entering from the bottom filaments show strong blueshifts. (As before, the line-of-sight velocity of the material in the upper-right quadrant of these panels is less uniform, as it probes a gas accreting along a filament that is roughly perpendicular to the line of sight.) The bottom panels of Figure 4 show that while the angular momentum content of the halo gas is high, the bulk of this inflowing gas does not have enough angular momentum to be fully rotationally supported. That is, most—but certainly not all—of the cold dense halo gas at R > 0.1R_vir has a tangential velocity V_tan/V_circ < 1 (where ${V}_{\mathrm{circ}}=\sqrt{{GM}(\lt R)/R}$ is the circular velocity at a given radius). Thus, despite the clear velocity structure shown in the middle panels of the figure, this high angular momentum gas should not be considered rotationally supported, but rather spiraling in toward the center of the halo, consistent with the short "sinking times" of ∼1–2 halo dynamical times previously reported by Stewart et al. (2011b).

The qualitative result in each case is a clear co-directional inflow signature, with cold dense halo gas easily divided by a single cutting plane into the redshifted versus blueshifted half of the halo, flowing through the halo via a chaotic assortment of high angular momentum inspiraling cold streams that are kinematically linked to inflow from the cosmic web.

Similar structures have been noted a number of times in the literature, but with a variety of terminologies, including the "messy region" (Ceverino et al. 2010), "cold flow disks" (Stewart et al. 2011b, 2013), the "AM sphere" (Danovich et al. 2012), or "extended rings" (Danovich et al. 2015). Indeed, depending on the simulation code utilized, one can easily see in Figure 4 how the kinematics and morphology of the inspiraling streams may or may not be well-described as a "messy region" (e.g., Ramses) or a more orderly disk-like structure (e.g., Art). Thus, while we find that the exact morphology—including size, orientation, clumpiness, thickness—of any structure that results from the inspiraling cold streams may vary significantly among simulations, each code does produces a qualitatively similar picture in which there is a clear line-of-sight velocity structure within the virial radius of the halo that is kinematically linked to that of the large-scale filamentary environment of Figure 3 (with the exception of the Gizmo-PSPH simulation at this epoch; see Section 5.2).

Figure 5 shows line-of-sight velocity maps (along an orthogonal orientation) for dense cold gas at z = 2 (top), z = 1.5 (middle), and z = 1 (bottom) for various subsets of the simulation runs (as labeled). The left panels of this figure look at the large-scale environment (analogous to Figure 3. While the basic filamentary nature of the gaseous inflows becomes less apparent at decreasing redshift (when the filaments are less dense), we can still note the same qualitative behavior of inflowing gas. On environmental scales, filamentary inflow results in the same clear line-of-sight velocity signature as before; across all simulations, gas flowing into the virial radius from the top of the panels is blueshifted, while gas flowing in from the bottom is redshifted, with the only notable exception being the Ramses code at z = 1, which is likely the result of the lack of self-shielding from the UV background, leaving very little cold gas above our minimum density threshold, so almost no cold dense inflow is still visible in the figure. While the detailed structure of the inflowing gas again varies among simulations, it seems apparent that filamentary gas accretion along a three-dimensional cosmic web onto an overdense region (at this mass scale in the redshift range 1 < z < 3) tends to produce the same qualitative picture across all the simulations, regardless of the subgrid physics.

Zoom In Zoom Out Reset image size

The right panels of Figure 5 shows an analogous line-of-sight velocity analysis, but zoomed-in to the halo virial radius for z = 2 (top), z = 1.5 (middle), and z = 1 (bottom), and again increasing the density threshold by a factor of 10 (similar to Figure 4). Again, the precise structure of the inspiraling cold streams varies among the simulations, but most of the simulations produce qualitatively similar pictures; there continues to be a clear large-scale velocity structure within the virial radius of the halo that is kinematically linked to the large-scale filamentary inflow shown in the left panels, again with the exception of Ramses at z = 1, which has evacuated most of its halo of cold dense gas altogether.

5.1. Co-directional Halo Gas

As a means of quantifying this result, we define the "co-directional mass fraction" in the following way, at any given epoch. For three arbitrary orthogonal projections, we define a cutting plane (passing through the center of the halo) that best divides the halo into positive versus negative line-of-sight velocities for all gas within the virial radius. Each particle (or cell, depending on code architecture) can then be defined as co-directional (along this projection) if its line-of-sight velocity was correctly categorized by this cutting plane.¹⁹ We then select the projection with the highest overall co-directional fraction for all gas in the halo (but not the galaxy: 0.1 <R/R_vir < 1.0). This selection typically corresponds to the projection in which the galaxy is seen closest to edge-on, though we note that this may not always be the case, if there is a significant misalignment between the angular momentum direction of the inspiraling gas and that of the galactic disk. The co-directional mass fraction of any given component (dark matter, cold gas, hot gas, or all gas) is the mass fraction that has been categorized as co-directional along this preferred projection.

Using this definition, the left panel of Figure 6 shows the co-directional mass fraction as a function of radius at z = 3 (for all gas versus dark matter). As might be expected, the average co-directional mass fraction for dark matter (taking into account all the simulations) is ∼50%, though with considerable variation depending on the exact orientation of the co-directional cutting plane. However, among all five simulations, the co-directional mass fraction for gas shows remarkably similar behavior with radius—declining smoothly from ∼85% at the galactic region (R = 0.1R_vir) to ∼70% at the virial radius. Indeed, even extending to ∼2.5R_vir, the co-directional mass fraction of gas in the cosmic web remains significantly higher than that of the dark matter, as expected if the line-of-sight velocity structure of the cold gas is kinematically linked to the filamentary gas flows beyond the virial radius of the halo.

Zoom In Zoom Out Reset image size

We explore the co-directional mass fraction over cosmic time in the middle and right panels of Figure 6, which use the same procedure outlined above to compute a single value at each epoch for the total co-directional mass fraction in the halo (but not the galaxy: 0.1 < R/R_vir < 1.0) for dark matter versus gas (middle panel), and further distinguishing between cold gas versus hot gas (right panel). In both panels, the average for dark matter, all gas, cold gas, and hot gas among all simulations and over the entire redshift range from $z=3-1$ are given by the horizontal dotted lines. As shown in the figure, the dark matter co-directional mass fraction varies somewhat sporadically around ∼50%, depending on the cutting plane orientation, while the halo gas shows significantly higher co-directional mass fractions (75% for all halo gas). While there are significant variations in the detailed results among the simulations, with some codes showing stronger co-directional kinematics than others, all simulations also demonstrate a higher tendency for cold gas to show this co-directional velocity structure in the halo over hot gas (79% versus 64%).

Since the halo gas mass is dominated by its cold component (see Figure 2), it is worth exploring whether the above trend in co-directional mass fraction between the cold and hot components might be the result of an offset in angular momentum direction between the hot versus cold gaseous halos. After repeating the above analysis, but using the line-of-sight velocity of the hot gas to define the co-directional cutting plane, we find only small variations in the above results. For example, the co-directional mass fraction in the halo at z = 3 for [Enzo, Art, Ramses, Arepo, Gizmo-PSPH] decreases slightly (if at all) for cold gas from [84%, 82%, 81%, 72%, 69%] to [79%, 81%, 81%, 72%, 69%], respectively, and increases slightly (if at all) for hot gas from [53%, 62%, 65%, 58%, 64%] to [55%, 62%, 68%, 58%, 64%]. Therefore, we conclude that the overall trends shown in Figure 6 are not highly sensitive to the way we define the co-directional cutting plane.

Taken together with our previous, more qualitative results, our findings suggest that across a broad range in hydrodynamic code types and subgrid physics models of galaxy formation, the presence of co-directional inspiraling cold streams in galaxy halos at z > 1 is a natural consequence of high angular momentum filamentary inflow along the cosmic web, and represents a robust prediction of cosmological gas accretion in LCDM.

5.2. The Rapid Destruction and Re-formation of Coherent Inspiraling Gas at z = 3

In the discussion of Figures 3 and 4, we noted that the velocity structure of the galaxy in the Gizmo-PSPH simulation is not nearly as clean and orderly as the other simulations at z = 3, and therefore does not seem to host the same clear co-directional velocity structure. While we thought it important to show all galaxies at precisely the same epoch, we note that in Gizmo-PSPH, the galaxy happens to be in the midst of a post-merger starburst, accompanied by a violent outflow event at this epoch, due to this code's strong feedback physics. Halo–halo mergers, of course, tend to occur at broadly similar times in all codes, but differences in galaxy masses and halo baryonic mass distributions mean that the galaxy–galaxy mergers can and do occur at significantly different times, and with different mass ratios and corresponding consequences for star formation, at the halo center (see, e.g., Stewart et al. 2009; Hopkins et al. 2010). The obvious clumpiness of the outflows may owe, at least partially, to well-known numerical difficulties capturing fluid-mixing instabilities in SPH, even in the improved P-SPH implementation; this is supported by early results from the FIRE-2 simulations, which use a different, mesh-free Godunov-type finite volume method to solve the hydrodynamics (P. Hopkins 2017, private communication).

In Figure 7, we show the structure of this galaxy immediately before and after this violent merger event. The time sequence begins in the left panel at z = 3.3, where inflowing cold gas demonstrates the same line-of-sight velocity structure as in Figure 4, including the presence of co-directional inspiraling cold streams, which have initially taken the apparent form of an extended disk-like structure reminiscent of the "cold flow disks" reported previously in cosmological simulations (e.g., Stewart et al. 2011b). To aid the eye in comparing the images, some of which are quite chaotic during the outflow, an identical ellipse has been overlaid on each image roughly corresponding to this coherent inspiraling gas region. At this epoch, the recent influx of fresh gas onto the central regions of the galaxy results in a spike in star formation, and consequently a violent spherical outflow event from z = 3.2–3.0 that effectively destroys the ISM of the galaxy (leaving a deficit of gas in the center of the galaxy, as seen at z = 3.1) and disrupting the inflowing filamentary gas in the CGM of the halo. However, the filamentary gas continues to flow into the halo, and this inflow continues to contribute substantial angular momentum. As a result, a new co-directional velocity structure becomes apparent almost immediately after the outflow event has subsided, with the co-directional mass fraction of cold dense gas in the halo (defined as described in Section 5.1) changing rapidly from 81% at z = 3.3 (before the outflow) to ∼60% during the outflow, back up to 76% at z = 2.8 (post-merger). By the rightmost panels, the co-directional inspiraling cold streams have once again formed a roughly disk-like structure, along a very similar orientation to the original inspiraling gas structure. We argue that the bursty nature (Muratov et al. 2015) of the subgrid physics as implemented in Gizmo-PSPH coupled with this demonstration of the near-immediate regrowth of the coherent inspiraling gas structure after a massive outflow event only reinforces the robust nature of inspiraling cold streams in the halos of massive galaxies in LCDM.

Zoom In Zoom Out Reset image size

Past studies of galaxy formation simulations have reported the existence of co-rotating structures of cool gas in the outskirts of galaxy halos. Our results indicate that inspiraling halo gas of this kind is robust to different feedback models and hydrodynamic solvers. One implication of this result is that extended, high-angular momentum cold stream configurations offer a testable observational signature of the LCDM galaxy formation paradigm. For example, the fact that a large fraction of the halo gas have velocities that are co-directional means that observations of the gas (e.g., from quasar absorption systems) will show a blueshifted and redshifted side that is usually interpreted as rotation. We emphasize that in this case most of the gas is far from being angular-momentum supported and that "inspiraling" halo gas is a more accurate description of its coherent motion than "rotating."

Encouragingly, there is a growing body of observational evidence that seems to indicate that co-directional halo gas is indeed seen around real galaxies. For example, kinematic studies of some Lyα nebulae suggest rotational velocities and inflow rates consistent with those expected for these inspiraling streams (Martin et al. 2014; Prescott et al. 2015). Similarly, absorption line studies are beginning to emphasize the bimodal distribution of absorption detections, where detections along the galaxy's minor axis tend to show absorption properties consistent with outflowing gas, while detections roughly along the galaxy's major axis demonstrate properties (such as co-rotational inflow) that are consistent with inspiraling cold streams (Kacprzak et al. 2010, 2012a, 2012b; Bouché et al. 2012, 2013; Crighton et al. 2013; Nielsen et al. 2015; Bouché et al. 2016; Ho et al. 2017).

Perhaps the most direct confirmation of the existence of inspiraling halo gas comes from Martin et al. (2015), who performed a spectroscopic analysis on the cosmic filament (illuminated by two nearby QSOs) first detected by Cantalupo et al. (2014) at z ∼ 2. They found that a substantial fraction of the illuminated region was in fact a huge co-rotating gaseous structure. The extremely extended gaseous "disk" (extending to ∼R_vir/2, corresponding to a width of 125 physical kpc) showed smooth rotation kinematics, with one side of the disk kinematically linked to the inflow velocity of the nearby cosmic filament. This very closely resembles what we have presented here for co-directional inspiraling cold streams, though we note that the particular system observed by Martin et al. (2015) was estimated to be a much more massive halo than what we have simulated here (M_vir ∼ 10¹³ M_⊙), and it therefore reported a correspondingly more massive and extended protogalactic disk than found in our simulations, as might be expected for a larger, more massive halo. A similar cold flow protodisk, again fed by a cosmic filament that was first detected in Lyα emission, was also reported in Martin et al. (2016), suggesting that inspiraling disk-like structures may be common phenomena for massive galaxies at high redshift.

While not seen in our particular simulations, we also speculate that polar ring galaxies—which have previously been suggested as evidence of cold flow gas accretion onto galaxies (Macciò et al. 2006; Brook et al. 2008; Spavone et al. 2010)—may be a result of a similar phenomenon. Such galaxies could reasonably occur when strong central torques (e.g., from a major galaxy merger) result in a near perpendicular misalignment between the angular momentum of the central galaxy and that of the inflowing cold-mode gas.

We note that the inspiraling cold streams in our simulations are significantly more massive and extended (relative to the halo virial radius) at high redshift, when cosmic filaments are more narrowly defined and contain higher density gas flows. However, Figure 2 and previous work (e.g., Stewart et al. 2013) both demonstrate that accreting cold gas continues to have high angular momentum, even at later times where the rotational signature of a continuous gaseous structure may be less clear. Additionally, Figure 6 demonstrates that the co-directional mass fraction of cold gas in the halo stayed consistently high over cosmic time (at least until z = 1). We speculate that it may be possible that this high angular momentum accretion helps to explain observations of extended XUV disks (e.g., Thilker et al. 2005, 2007; Lemonias et al. 2011; Holwerda et al. 2012), local extended H i disks (e.g., García-Ruiz et al. 2002; Oosterloo et al. 2007; Christlein & Zaritsky 2008; Sancisi et al. 2008; Walter et al. 2008; Wang et al. 2013; Huang et al. 2014; Courtois et al. 2015), and co-rotating cold halo gas around local Milky Way analogs (e.g., Diamond-Stanic et al. 2016).

Indeed, these growing observations of high angular momentum material in the outskirts of galaxy halos would be quite difficult to explain if one were to assume the canonical picture of galaxy formation whereby baryons in galaxy halos share the same distribution of angular momentum as the dark matter. In contrast, the cold flow paradigm naturally predicts that halo gas (and particularly the cold halo gas) preferentially constitutes recent gas accretion from the cosmic web, with ∼3–5 times the angular momentum of the dark matter, naturally explaining the kinds of high angular momentum phenomena being observed. We caution, however, that we have not focused on the gaseous halos in our simulations at z < 1 here, and leave a more detailed comparison between simulations and low-z observations as a topic of further study.

We have simulated the evolution of a Milky Way-sized galaxy from identical cosmological initial conditions with a variety of simulation codes: Enzo, Art, Ramses, Arepo, and Gizmo-PSPH. Each code has used subgrid physics models drawn from scientific literature common to each simulation type, and we have compared the simulations in an attempt to draw robust conclusions about galaxy formation in LCDM (focusing on z > 1) that are not sensitive to uncertain aspects of galaxy formation simulations. To ensure uniform analysis among the various code types, we have used the hydrodynamic analysis software yt, which enables the same analysis routine to be performed on each code.

While we found many aspects of the simulated galaxies that did vary substantially among the simulations (e.g., morphology, stellar mass, hot gas halo temperature, and mass, to be discussed further in future work), we found the following qualitative features common to all simulations, regardless of which subgrid physics model or hydrodynamic code was used:

• 1.  
  Gas in the galaxy halo has substantially higher specific angular momentum than the dark matter in the halo, with mean values of ${j}_{\mathrm{cold}}\simeq 4{j}_{\mathrm{DM}}$ and ${j}_{\mathrm{hot}}\simeq 2{j}_{\mathrm{DM}}$ (though with considerable scatter), leading to a typical cold halo gas spin parameter of λ_cold ≃ 0.12.
• 2.  
  The large-scale filamentary structure is qualitatively similar in all simulations (with minor variations, for example regarding lower mass streams of secondary importance to the galaxy's growth). The three-dimensional geometry of these filaments, which are all flowing toward the central galaxy (the highest overdensity in its environment) results in a strong line-of-sight velocity structure. Filaments flowing onto the galaxy from opposite directions (along an arbitrary line of sight) tend to show alternating blueshifted and redshifted velocities relative to the galaxy as they flow toward the galaxy center.
• 3.  
  As the filamentary gas accretion enters the virial radius, the large-scale velocity structure of the accreting filaments inevitably results in inspiraling cold streams in the halo of the galaxy that carry significant angular momentum as they spiral in from the virial radius to the galactic region. For the Milky Way-sized halo simulated here, the maximum line-of-sight velocity expected for these inspiraling streams is ∼250 km s⁻¹ (corresponding to roughly 1.5 times the virial velocity of the halo). As a result, the vast majority (∼80%) of cold halo gas follows a clear co-directional velocity structure (with a single cutting plane dividing positive versus negative line-of-sight velocities) as the cold streams spiral toward the center of the halo. In contrast, the co-directional mass fractions in the halo are considerably lower for the hot gas (∼65%) or the dark matter (∼50%).
• 4.  
  Inspiraling cold streams occasionally take the previously reported morphology of cold flow disks: high angular momentum cold gas that is transitioning from the cosmic web, though the halo as a roughly disk-like structure (except that the gas is not angular-momentum supported) and will eventually accrete onto the galactic disk. These coherent inspiraling structures represent continuous and dynamic flows from the cosmic web; even after a violent outflow event disrupts the CGM in one of the simulations, the newly inflowing gas rapidly re-forms a similar inspiraling structure along the same orientation soon after the outflow has subsided.

In this work, we have limited our analysis to the growth of a single Milky Way-sized halo at z > 1 using a variety of different hydrodynamic codes and feedback physics implementations. It is therefore difficult to draw general conclusions about galaxy formation from the simulation of a single halo; however, a number of theoretical works have previously established the high angular momentum nature of filamentary gas accretion, using various hydrodynamic codes, larger cosmological volumes, and/or analysis of multiple zoom-in simulations. For example, Pichon et al. (2011) analyzed ∼15,000 halos at z > 1.5 from a (lower resolution) cosmological-scale simulation using the Ramses code; in a companion work, Kimm et al. (2011) also included ∼900 intermediate-resolution halos and two high-resolution zoom-in simulations to z = 0 using Ramses; Stewart et al. (2011b, 2013) analyzed four zoom-in simulations to z = 0 using the SPH code Gasoline; and Danovich et al. (2015) analyzed 29 zoom-in simulations at z > 1.5 using the Art code. The results presented here demonstrate that the high angular momentum nature of cold gas accretion in LCDM is not likely to change (in the qualitative sense) among a broad range of different physics implementations and hydrodynamic codes, suggesting that the presence of inspiraling cold streams in galaxy-size halos appears to be a robust expectation of LCDM.

However, we note that there are considerable variations in the quantitative nature (morphology, rotational velocity, size, temperature, density, etc.) of the inspiraling cold streams in each of the simulations performed in this work. The detailed properties and prevalence of these inspiraling streams are yet to be fully understood, and cannot be determined from the single high-resolution simulation presented here. The co-directional velocity structure noted here is also likely to depend on the geometry and kinematics of the cosmic web in the galaxy's environment, so we speculate that there are likely to be significant environment effects, even at fixed halo mass. For example, Milky Way-sized galaxies near the outskirts of galaxy clusters would not be expected to dominate the gravitational potential of the cosmic web in their large-scale environment, so we may not expect to find the same clear co-directional velocity signature for filamentary inflow (i.e., Figure 3) for such systems either.

We also take special note that in several of the simulation codes used here (and in previous works), these inspiraling streams result in the formation of transient inspiraling disk-like structures qualitatively similar to the "cold flow disks" of, e.g., Stewart et al. (2013). The exploration of the prevalence of these inspiraling gaseous structures in simulations, for different environments and halo masses, would be a useful topic of further study, especially in light of recent observations of large co-rotating gaseous structures at z ∼ 2 (Martin et al. 2015, 2016) that are strikingly similar to the qualitative results presented here.

Computations described in this work were performed using the publicly available Enzo code (http://enzo-project.org), and the publicly available yt toolkit (http://yt-project.org/), both of which are the products of collaborative efforts of many independent scientists from numerous institutions around the world. Their commitment to open science has helped make this work possible. K.R.S. would especially like to thank Matt Turk and Nathan Goldbaum for their support with yt. K.R.S. also thanks Volker Springel, Lars Hernquist, and Paul Torrey for running the Arepo simulation discussed in this work, for allowing us to include its results as part of this comparison project, and for providing useful discussions and comments.

K.R.S. was supported by HST-GO-14268.026-A. J.S.B. was supported by HST AR-12836. C.A.F.G. was supported by NSF grants AST-1412836 and AST-1517491, by NASA grant NNX15AB22G, and by STScI grants HST-AR-14293.001-A and HST-GO-14268.022-A. J.D. acknowledges support from the Spin(e) grant ANR-13-BS05-0005 of the French Agence Nationale de la Recherche (http://cosmicorigin.org). Support for P.F.H. was provided by an Alfred P. Sloan Research Fellowship, NASA ATP Grant NNX14AH35G, and NSF Collaborative Research Grant #1411920 and CAREER grant #1455342. The Art run was performed at the National Energy Research Scientific Computing Center (NERSC) at Lawrence Berkeley National Laboratory. DC acknowledges support from the European Research Council under the European Community's Seventh Framework Programme (FP7/2007-2013) via the ERC Advanced Grant "STARLIGHT: Formation of the First Stars" (project number 339177). The Ramses simulation was performed on the DiRAC Facility, jointly funded by BIS and STFC. J.D.'s research is partly funded by Adrian Beecroft and the Oxford Martin School. The Gizmo-PSPH run was carried out on the CUNY HPC which is supported, in part, under National Science Foundation Grants CNS-0958379, CNS-0855217, and ACI-1126113, and the City University of New York High Performance Computing Center at the College of Staten Island. Storage and analysis of the simulations was performed on the CTP cluster, which is supported by GRTI grant CT04AGR15001 and the Physics Division of the U.S. Army Research Office grant #64775-PH-REP.