THE HYDRODYNAMIC STABILITY OF GASEOUS COSMIC FILAMENTS

The thermodynamic state of gas in cosmic-web filaments has important implications for observations and theoretical predictions. It will affect halos that are fed by those filaments, as well as small halos that accrete onto those filaments as part of the cosmic hierarchical growth.

Spherical virialization of gas in halos has been a prediction of galaxy formation models for decades. In particular, it has been shown (Binney 1977; Rees & Ostriker 1977; Silk 1977; White & Frenk 1991) that a comparison between cooling times and ages of galactic halos can predict the transition from galaxies to groups and clusters. Birnboim & Dekel (2003, hereafter BD03) derived a stability criterion against gravitational collapse of the gas in the presence of significant cooling. They find that for halos below M_crit ≃ 10¹²$\,{M}_{\odot }$ a hot gaseous halo is not expected to form, and gas will free-fall until it reaches the disk, at which point it will stop, radiating its kinetic energy abruptly at that point. This has been confirmed in multiple hydrodynamical simulations (e.g., Kereš et al. 2005; Ocvirk et al. 2008; Faucher-Giguère et al. 2011) and successfully reproduces star-forming galaxies at high-z (Dekel & Birnboim 2006; Dekel et al. 2009) and the color–magnitude bi-modality (Cattaneo et al. 2006; Croton et al. 2006; Dekel & Birnboim 2008). Observational indications of this scenario are gradually accumulating (e.g., Dijkstra & Loeb 2009; Kimm et al. 2011; Martin et al. 2015).

In this Letter, we derive a criterion for the stability of virial shocks around filaments and sheets that form the cosmic web, analogous to the BD03 criterion for halos. Following Fillmore & Goldreich (1984, hereafter FG84), we construct self-similar density profiles of filaments. We apply our stability analysis to these profiles to identify filaments around which a stable virial shock is expected to form. This criterion is translated to a more useful form by identifying which halos are expected to be fed by these filaments. We find that filament instability influences a large portion of halos in the universe throughout cosmic age.

The stability of filaments has been addressed before, numerically (Harford & Hamilton 2011) and analytically (Breysse et al. 2014; Freundlich et al. 2014), but without taking into account cooling, and by analyzing the stability of initially static filaments, ignoring the effects of the shock at the filaments' edge.

In Section 2, we derive the stability criterion for the existence of virialized gas in one-, two-, and three-dimensional collapse. In Section 3, we relate our local criterion to cosmic filaments according to the self-similar solutions of FG84. In Section 4, we relate these filaments to typical halo masses that will likely be fed by them. In Section 5, we summarize and conclude.

The analysis in BD03 was performed for spherical accretion in the presence of cooling. In this section, we generalize that derivation for infall onto spherical, cylindrical (onto filaments), and planar (onto sheets or disks) objects.

The ideal gas equation of state (EoS) is

$$\begin{eqnarray}&&P=(\gamma -1)\rho e,\end{eqnarray} \tag{ 1 }$$

with $\gamma ={\left(\tfrac{\partial \mathrm{ln}P}{\partial \mathrm{ln}\rho }\right)}_{s}$ as the adiabatic index, and ρ, e, P as the density, internal energy, and pressure, respectively. The adiabatic index measures the "stiffness" of the EoS, or the adiabatic pressure response to compression. By analogy, we define the effective EoS index, γ_eff, of a parcel of gas undergoing compression as its pressure response along a Lagrangian trajectory:

$$\begin{eqnarray}&&{\gamma }_{\mathrm{eff}}\equiv \displaystyle \frac{d\,\mathrm{ln}P}{d\,\mathrm{ln}\rho }=\displaystyle \frac{\rho }{P}\displaystyle \frac{\dot{P}}{\dot{\rho }},\end{eqnarray} \tag{ 2 }$$

with the upper dot implying a Lagrangian time derivative. Following BD03, this form is transformed into

$$\begin{eqnarray}&&{\gamma }_{\mathrm{eff}}=\gamma -\displaystyle \frac{\rho }{\dot{\rho }}/\displaystyle \frac{e}{q}.\end{eqnarray} \tag{ 3 }$$

with q related to the Sutherland & Dopita (1993) cooling rate according to

$$\begin{eqnarray}&&q=\rho {{\rm{\Lambda }}}_{\mathrm{cool}}(T,Z)=\rho ({N}_{a}^{2}{\chi }^{2}/{\mu }^{2}){{\rm{\Lambda }}}_{\mathrm{mic}}(T),\end{eqnarray} \tag{ 4 }$$

where N_a/μ is the average number of molecules per unit mass and χ is the number of electrons per particle. This parameterization allows for altered ionization state of the gas (see Cantalupo 2010). This equation illustrates how the effective stiffness of the gas is set by two competing timescales: the compression time, $\tfrac{\rho }{\dot{\rho }}$, and the cooling time, $\tfrac{e}{q}$. If radiative losses are significant during a compression time, we get effective softening, γ_eff < γ. In the absence of cooling, we recover γ_eff = γ.

The compression time for post-shock infalling gas depends on the velocity and geometry of the infall through the Lagrangian continuity equation,

$$\begin{eqnarray}&&\dot{\rho }=-\rho {\rm{\nabla }}\cdot u=-\rho \left((n-1)\displaystyle \frac{u}{r}+\displaystyle \frac{\partial u}{\partial r}\right),\end{eqnarray} \tag{ 5 }$$

with u being the velocity of the gas, and n the dimensionality of the infall: n = 1, 2, 3 correspond to planar, filamentary, and spherical collapse, respectively. Assuming that the post-shock flow is homologous,

$$\begin{eqnarray}&&u=\displaystyle \frac{{u}_{s}}{{r}_{s}}r\end{eqnarray} \tag{ 6 }$$

with r_s as the shock radius and u_s as the velocity directly below the shock, we can apply Equation (5) to show that the post-shock gas must contract uniformly, independent of r,

$$\begin{eqnarray}&&\displaystyle \frac{\dot{\rho }}{\rho }=-n\displaystyle \frac{{u}_{s}}{{r}_{s}}=\mathrm{Const}.\end{eqnarray} \tag{ 7 }$$

This homologous behavior is seen in 1D simulations for the spherical case (Birnboim & Dekel 2003), and is a standard assumption in hydrostatic stability calculations. It also roughly matches self-similar solutions of gaseous spherical infall (Bertschinger 1985).

We consider a quasi-static configuration where the shock radius, r_s, and the post-shock velocity profile, Equation (6), are approximately constant. We test the stability of this configuration by assuming that it is initially hydrostatic:

$$\begin{eqnarray}&&\ddot{r}=-\displaystyle \frac{1}{\rho }{\rm{\nabla }}P+{a}_{g}=0\end{eqnarray} \tag{ 8 }$$

(a_g being the gravitational acceleration) but has an inward post-shock velocity, and checking the sign of the auxiliary acceleration or "jerk," $\dddot{r}$, that forms due to this motion. A positive jerk will instigate outward motion implying a stable configuration, while a negative $\dddot{r}$ will lead to collapse, implying instability.²

The gravitational acceleration depends on the dimensionality of the potential well. Since the perturbation is Lagrangian, the mass enclosed below the gas parcel is constant, so

$$\begin{eqnarray}&&{a}_{g}=-{{Ar}}^{1-n},\end{eqnarray} \tag{ 9 }$$

with A as a positive constant.³ We convert the ${\rm{\nabla }}=\tfrac{d}{{dr}}$ operator in Equation (8) to a Lagrangian mass derivative by noting that a mass element is related to a spatial differential according to

$$\begin{eqnarray}&&{dm}=B\rho {r}^{n-1}{dr},\end{eqnarray} \tag{ 10 }$$

with B as a positive constant.⁴ Plugging Equations (9) and (10) into Equation (8), we find

$$\begin{eqnarray}&&\ddot{r}=-{{Br}}^{n-1}\displaystyle \frac{{dP}}{{dm}}-{{Ar}}^{1-n}=0.\end{eqnarray} \tag{ 11 }$$

From here on we shall denote all mass derivatives with '.

The rate of change in the acceleration is the time derivative of Equation (11):

$$\begin{eqnarray}&&\dddot{r}=-B(n-1){r}^{n-2}{uP}^{\prime} -{{Br}}^{n-1}\dot{P^{\prime} }-A(1-n){r}^{-n}u,\end{eqnarray} \tag{ 12 }$$

noting that $u=\dot{r}.$ Eliminating the last term by use of Equation (11), collecting terms, and exchanging the spatial and time derivative of P, one gets

$$\begin{eqnarray}&&\dddot{r}=-{{Br}}^{n-1}\left[2(n-1)\displaystyle \frac{u}{r}P^{\prime} +(\dot{P})^{\prime} \right].\end{eqnarray} \tag{ 13 }$$

The calculation of $(\dot{P})^{\prime}$ is somewhat lengthy. We first derive an expression for $\dot{P}$:

$$\begin{eqnarray}\dot{P} & = & {\gamma }_{\mathrm{eff}}\displaystyle \frac{\dot{\rho }}{\rho }P=\gamma \displaystyle \frac{\dot{\rho }}{\rho }P-\displaystyle \frac{q}{e}P\\ & = & \gamma \displaystyle \frac{\dot{\rho }}{\rho }P-\displaystyle \frac{\rho P{{\rm{\Lambda }}}_{\mathrm{cool}}}{e}=\gamma \displaystyle \frac{\dot{\rho }}{\rho }P-(\gamma -1){\rho }^{2}{{\rm{\Lambda }}}_{\mathrm{cool}},\end{eqnarray} \tag{ 14 }$$

with the first equality due to Equation (2), the second equality to Equation (3), third to Equation (4), and fourth to Equation (1). By Equation (7) we note that the term $\tfrac{\dot{\rho }}{\rho }$ is independent of the spatial derivative. We further assume that the cooling function, Λ_cool, does not change significantly due to the change of temperature of the Lagrangian mass element⁵ , so Λ_cool can also be taken out of the derivative. We assume that nearby mass elements directly below the shock start with the same thermodynamic conditions (i.e., they lie on the same adiabat). Combined with the definition in Equation (2), this indicates that

$$\begin{eqnarray}&&\displaystyle \frac{P^{\prime} }{\rho ^{\prime} }={\gamma }_{\mathrm{eff}}\displaystyle \frac{P}{\rho }.\end{eqnarray} \tag{ 15 }$$

Differentiating Equation (14) then yields

$$\begin{eqnarray}(\dot{P})^{\prime} & = & \gamma \displaystyle \frac{\dot{\rho }}{\rho }P^{\prime} -2(\gamma -1)\rho {{\rm{\Lambda }}}_{\mathrm{cool}}\rho ^{\prime} \\ & = & \gamma \displaystyle \frac{\dot{\rho }}{\rho }P^{\prime} -2(\gamma -1)q\displaystyle \frac{1}{{\gamma }_{\mathrm{eff}}}\displaystyle \frac{\rho }{P}P^{\prime} \\ & = & \gamma \displaystyle \frac{\dot{\rho }}{\rho }P^{\prime} -\displaystyle \frac{2}{{\gamma }_{\mathrm{eff}}}\ \displaystyle \frac{q}{e}P^{\prime} =\gamma \displaystyle \frac{\dot{\rho }}{\rho }P^{\prime} -\displaystyle \frac{2}{{\gamma }_{\mathrm{eff}}}(\gamma -{\gamma }_{\mathrm{eff}})\displaystyle \frac{\dot{\rho }}{\rho }P^{\prime} .\end{eqnarray} \tag{ 16 }$$

The second equality is due to Equation (15), the third by reverse use of Equations (1) and (4), and the fourth by separating $\tfrac{q}{e}$ from Equation (3).

Inserting Equation (16) into Equation (13), and converting $\tfrac{\dot{\rho }}{\rho }$ and $\tfrac{u}{r}$ to $\tfrac{{u}_{s}}{{r}_{s}}$ according to Equations (6) and (7), we finally get

$$\begin{eqnarray}&&\dddot{r}=-{{Br}}^{n-1}\displaystyle \frac{{u}_{s}}{{r}_{s}}P^{\prime} \left[2(n-1)-n\gamma +\displaystyle \frac{2n}{{\gamma }_{\mathrm{eff}}}(\gamma -{\gamma }_{\mathrm{eff}})\right].\end{eqnarray} \tag{ 17 }$$

u_s < 0 because below the standing shock there is inward velocity, and P' < 0 so that the pressure force, −∇P, is positive to balance gravitation in the quasi-hydrostatic halo, so the factor before the square brackets is always negative. A positive jerk, or stability, thus occurs when

$$\begin{eqnarray}&&2(n-1)-n\gamma +\displaystyle \frac{2n}{{\gamma }_{\mathrm{eff}}}(\gamma -{\gamma }_{\mathrm{eff}})\lt 0.\end{eqnarray} \tag{ 18 }$$

We first note that in the absence of cooling γ_eff = γ and the stability condition reduces to

$$\begin{eqnarray}&&\gamma \gt 2-\displaystyle \frac{2}{n}.\end{eqnarray} \tag{ 19 }$$

For the spherical case (n = 3) the stability criterion is $\gamma \gt \tfrac{4}{3},$ recovering a well-known result. For filamentary accretion (n = 2) of adiabatically collapsing gas, stability is gained when γ > 1, and for planar collapse (n = 1) when γ > 0.

When cooling is present the stability criterion is

$$\begin{eqnarray}&&{\gamma }_{\mathrm{eff}}\gt \displaystyle \frac{2n\gamma }{n\gamma +2},\end{eqnarray} \tag{ 20 }$$

which, for monoatomic gas $\left(\gamma =\tfrac{5}{3}\right)$, is $\tfrac{10}{7},\,\tfrac{5}{4},\,\tfrac{10}{11}$ for spherical, cylindrical, and planar collapse, respectively.

We note that for $\gamma =\tfrac{5}{3}$ the critical value of γ_eff in the presence of cooling is always somewhat larger than the adiabatic critical value. Hence, monoatomic gas that is cooling with some local γ_eff (with entropy and energy decreasing due to cooling) is always less stable than an adiabatic gas with a softened EoS γ = γ_eff.

Filaments grow (in girth) by accreting gas from their surroundings. The infall geometry is primarily cylindrical, although most of the gas is channeled along sheets, and becomes spherical in the vicinity of halos. The flow parallel to the filament below and above the shock is continuous and can be factored out locally with a proper shift of the frame of reference. Dark matter (DM) and gas accrete onto filaments together. Separation between gas and DM occurs when gas becomes thermalized and is decelerated due to its pressure. We wish to determine where and how this gas is thermalized, particularly in the presence of cooling, which softens the effective EoS of the gas. To do so, we must connect global properties of filaments: their mass per unit length and their density and velocity profiles, to local cooling and contraction rates that determine the stability of the gas.

An ideal framework for connecting the large-scale properties of filaments to local conditions is the self-similar solutions of FG84. We do not present here the derivation and results and refer the reader to the original paper. Following FG84, we numerically solve for the self-consistent density profile and trajectories of infalling cylindrical DM shells, starting from an initial mass perturbation within an Einstein de-Sitter universe. In this framework, a filament is characterized by its mass per unit length, and by its initial perturbation:

$$\begin{eqnarray}{M}_{\mathrm{fil}}(r) & = & {M}_{\mathrm{fil}}^{0}(r)\left(1+\displaystyle \frac{\delta {M}_{\mathrm{fil}}(r)}{{M}_{\mathrm{fil}}^{0}(r)}\right)\\ {M}_{\mathrm{fil}}^{0}(r) & \equiv & \pi {r}^{2}{\rho }_{u},\end{eqnarray} \tag{ 21 }$$

with ρ_u as the universal density at the initial time. The perturbation is defined as a function of the unperturbed mass:

$$\begin{eqnarray}&&\displaystyle \frac{\delta {M}_{\mathrm{fil}}}{{M}_{\mathrm{fil}}^{0}}={\left(\displaystyle \frac{{M}_{\mathrm{fil}}^{0}}{{M}_{* }}\right)}^{-\epsilon },\end{eqnarray} \tag{ 22 }$$

with M_* as a reference mass. varies between 0 and 1, where 1 corresponds to the most localized perturbation that still grows with ${M}_{\mathrm{fil}}^{0}$, and 0 to a long range perturbation for which the density of the perturbation is still decreasing with ${M}_{\mathrm{fil}}^{0}.$ Figure 1 shows the resulting mass profiles (normalized to the virial radius and the virial mass) for various values of . It is evident that the profile within the virial radius of the filament depends weakly on . The visible discontinuities in gradient correspond to caustics in DM shells as they turn around consecutively. The outermost caustic, or "first shell crossing," is defined as the virial radius. The mass profiles in Figure 1 are normalized to that radius.

Zoom In Zoom Out Reset image size

Figure 2 shows the self-similar trajectory of a DM shell and the filament's mass profile. The trajectory's radius is normalized to the turnaround radius at each time, so before turnaround, at τ < 1, the spatial coordinate λ/Λ > 1. The mass profile as a function of that same spatial variable is also present. The vertical black line corresponds to the event of first crossing (the virial radius) and is roughly where the virial shock will occur for a gaseous shell.

Zoom In Zoom Out Reset image size

From the self-similar solution we extract the infall velocity and density at every radius. The infalling velocity and density are for DM trajectories and correspond to infalling gaseous shells only before they pass through the virial shock, at which point their pressure becomes significant and their trajectories diverge from those of the DM. The stability of the post-shock gas depends on the local compression rate, density, and temperature for the post-shocked gas. These values are approximated from the pre-shocked ones by use of the strong shock approximation, which is valid as long as the pre-shocked velocity is much larger than the pre-shocked speed of sound $({c}_{s}^{0})$. Assuming that the ${c}_{s}^{0}\lesssim 10\,\mathrm{km}\,{{\rm{s}}}^{-1},$ we show later (Figures 3 and 4) that this approximation goes from being marginally satisfied for the smallest filaments to being fully justified for the large filaments. Using the full shock conditions will not change the results significantly and requires knowledge of the thermodynamic state of the cold gas. As the effective EoS becomes softer, the virial shock ceases to expand and starts to collapse. At its critical state, we expect the shock to be at rest. Using the strong shock approximation and assuming the shock is at rest, the post-shock values are

$$\begin{eqnarray}{\rho }_{1} & = & \displaystyle \frac{\gamma +1}{\gamma -1}{\rho }_{0}\\ {u}_{1} & = & \displaystyle \frac{\gamma -1}{\gamma +1}{u}_{0}\\ {e}_{1} & = & \displaystyle \frac{2{u}_{0}^{2}}{{(\gamma +1)}^{2}}\\ {T}_{1} & = & \displaystyle \frac{\mu (\gamma -1)}{{N}_{a}{k}_{{\rm{B}}}}{e}_{1},\end{eqnarray} \tag{ 23 }$$

with subscript 0 and 1 denoting pre-shocked and post-shocked variables, respectively; ρ₀, u₀ given from the numerical solution of the self-similar collapse; and e₁, T₁ are the internal energy and temperature of the post-shock variables. μ = 0.61 is the mean molecular weight for primordial, fully ionized gas, and N_a, k_B are Avogadro's number and Boltzmann's constant. The conversion from DM density to gas density is achieved by multiplying the density by a universal baryonic fraction (f_b = 0.17 throughout this work). This value is reasonable as long as the DM and gas flow together, i.e., for pre-shocked gas.

Zoom In Zoom Out Reset image size

For a filament characterized by M and we calculate the post-shock values for every radius below r_vir and use Equations (7) and (4) inserted into Equation (3) to calculate γ_eff. Figures 3 and 4 show γ_eff as a function of M_fil and $r/{R}_{\mathrm{vir}},$ with R_vir defined as the radius of the first shell crossing of the self-similar solutions (see Figure 2). In regions where γ_eff drops below the threshold for stability, γ_eff < γ_crit = 1.25, the filament cannot sustain a virial shock. Moreover, if a shock were to form in regions where γ_eff < 0, the post-shock pressure would decrease even as it contracts. The mass range for which the filament is unstable at least at some radius for both values of is between 10¹¹ and 10¹³$\,{M}_{\odot }$ Mpc⁻¹. In these plots, radii of density caustics leave horizontal features, and peaks in the cooling curve create features parallel to infall velocities. For each radius, the sharp drop in γ_eff as mass exceeds the lower threshold for instability is due to the post-shock temperature exceeding 10⁴ K, where the cooling rate grows by many orders of magnitude.

Zoom In Zoom Out Reset image size

We now wish to relate the filament masses (per unit length) to the typical masses of halos fed by such filaments. A full analysis of the filament distribution that accrete onto certain halos requires cosmological N-body simulations. For simplicity, we choose an alternative avenue that approximates the relation between halos and their filaments. Danovich et al. (2012) study the filamentary nature of mass accretion onto high-redshift galaxies. They find that typically, halos that originate from high-σ peaks in the initial perturbation accrete most (f_acc ≃ 70%) of their mass in filaments, out of which, 95% originates from the combined flow in the three largest filaments. Although the work analyzes high-redshift galaxies, we expect the same to be true for low-redshift clusters, which are also high-σ peaks. Using these values, we estimate the typical accretion rate through each filament as

$$\begin{eqnarray}&&{\dot{M}}_{\mathrm{acc}}=\displaystyle \frac{{f}_{\mathrm{acc}}}{3}{\dot{M}}_{\mathrm{halo}},\end{eqnarray} \tag{ 24 }$$

with ${\dot{M}}_{\mathrm{acc}},{\dot{M}}_{\mathrm{halo}}$ as the mass flow rate of gas within a filament and halo accretion rate, respectively. For the halo accretion rate we use the fit from Neistein et al. (2006):

$$\begin{eqnarray}&&\displaystyle \frac{{\dot{M}}_{\mathrm{halo}}}{{M}_{\mathrm{halo}}}=\left(\displaystyle \frac{0.03}{\mathrm{Gyr}}\right){(1+z)}^{2.5},\end{eqnarray} \tag{ 25 }$$

with z as the redshift. Finally, the mass of the filament per unit length is related to the flow rate by ${M}_{\mathrm{fil}}={\dot{M}}_{\mathrm{acc}}/{v}_{\mathrm{vir}}({M}_{\mathrm{halo}},z)$, using the halo virial velocity, v_vir, as an estimate for filaments' velocity as they accrete onto halos. A typical filament feeding a high-σ peak halo of mass M_halo will thus have an estimated mass of

$$\begin{eqnarray}&&{M}_{\mathrm{fil}}\simeq \displaystyle \frac{{f}_{\mathrm{acc}}}{3}\left(\displaystyle \frac{0.03}{\mathrm{Gyr}}\right){(1+z)}^{2.5}\displaystyle \frac{{M}_{\mathrm{halo}}}{{v}_{\mathrm{vir}}}.\end{eqnarray} \tag{ 26 }$$

In Figure 5, we use the inverse of this transformation to show the stable fraction of a filament, f_stable, as a function of the halo mass and redshift:

$$\begin{eqnarray}&&{f}_{\mathrm{stable}}\equiv \displaystyle \frac{1}{{R}_{\mathrm{vir}}}{\int }_{0}^{{R}_{\mathrm{vir}}}{\rm{\Theta }}[{\gamma }_{\mathrm{eff}}(r)-1.25]{dr},\end{eqnarray} \tag{ 27 }$$

with R_vir as the virial radius of the filament, and Θ as the Heaviside function. From Figure 5 it is evident that for 10¹⁰ M_⊙ ≲ M_halo ≲ 10¹³ M_⊙ at z = 3, as well as for 10¹² M_⊙ ≲ M_halo ≲ 10¹⁵ M_⊙ at z = 0, halos are expected to be fed by filaments that are not in hydrostatic stability.

Zoom In Zoom Out Reset image size

We have shown that in the presence of significant cooling, the accretion process of gas onto cosmic-web structures will not always proceed according to the standard virialization scenario of the infall-heating-cooling sequence. The analysis shown here can be applied for accretion onto spherical halos, cylindrical filaments, and planar sheets. For filament of 10¹¹–10¹³$\,{M}_{\odot }$ Mpc⁻¹, we show that gas is expected to fall without ever passing a shock, resulting in dense, thin filaments with low entropy. This is in complete analogy to spherical cold accretion onto halos that have been shown in BD03 and demonstrated in observations and simulations.

Using a simplified toy model for the relation between halo mass and redshift to typical filaments that feed it, we show that throughout cosmic history galaxies and clusters are affected by that instability. In particular, high-z star-forming galaxies (M_halo > 10¹⁰$\,{M}_{\odot }$ at z = 3) and low-redshift groups and clusters (10¹²$\,{M}_{\odot }$ ≲ M_halo ≲ 10¹⁵$\,{M}_{\odot }$ at z = 0) will be fed by filaments for which the gas is unstable.

The process that eventually stops the infall is still unclear, and we postulate that it is either angular momentum support from an original helicity of the filament or by reduction of the cooling rate due to self-shielding of the gas. A prediction of this work is thus that filament gas, in the non-stable regime, will be highly rotating and angular momentum supported. Both processes are hard to identify in simulations and have not been examined so far. In their absence, gas in simulations will flow toward the center of the filament until it approaches the numerical gravitational smoothing length, at which point the force will diminish. This indicates that the density and entropy of gas in unstable filaments are a numerical artifact and will not converge to the right values. This problem will be examined in future work.

The lack of virialized gas in filaments is expected to significantly affect the outcome of galaxies falling onto the filament and of halos fed by the filament. Halos falling onto filaments are expected to lose gas through ram pressure stripping and to enrich the filament with metals. Both these processes will be suppressed when galaxies fall into filaments with no stable atmosphere. Penetration of cosmic-web gas directly to galaxies affects the ISM state and the gas available for star formation and AGNs, as well as their feedback efficiencies. Mandelker et al. (2016) analyze the Kelvin–Helmholtz stability of supersonic filaments. They find that filaments lose stability via bulk modes that correspond to standing waves reflecting through edges of the filament. These results do not account for the effects of gravitational attraction toward the center of the filament and to angular momentum support, both expected to stabilize the filament further. These effects will be addressed in future work.

Observationally, the temperature of the filament could affect its detectability through Lyα absorption (Narayanan et al. 2010; Wakker et al. 2015) and emission (Martin et al. 2015). The temperature of the filaments will also affect the soft X-ray background and the total amount of gas in the "warm phase" (Cen & Ostriker 1999; Davé et al. 2001). All these effects are left for analysis in future work.

We thank Oliver Hahn for useful discussions. Y.B. and D.P. have been supported by ISF grant 1059/14.