The abundance of dark galaxies

Abstract

We show that gas in a large fraction of low-mass dark matter haloes may form Toomre-stable discs, if angular momentum is conserved when the gas contracts. Such haloes would be stable to star formation and therefore remain dark. This may potentially explain deviation between the predicted and the observed faint-end slope of the luminosity function and the discrepancy between the predicted and observed number of dwarf satellites in the Local Group. The above mechanism does not require a strong variation of the baryon fraction with the virial mass of the dark halo. We show that model fits to rotation curves are also consistent with this hypothesis: none of the observed galaxies lies in the region of parameter space forbidden by the Toomre stability criterion.

galaxies: formation, galaxies: kinematics and dynamics, galaxies: spiral, cosmology: theory

1 Introduction

There is now growing evidence that the observed number of dwarf galaxies surrounding the Milky Way (MW) or M31 does not agree with the predictions of high-resolution N-body simulations for the currently favoured Cold Dark Matter cosmological model (flat, low matter-density universe – LCDM) (Klypin et al. 1999; Moore et al. 1999). Given the successes of the current LCDM paradigm in reproducing many of the large-scale observables (e.g. Jaffe et al. 2001; Peacock et al. 2001; Lahav et al. 2002; Verde et al. 2001), it seems natural to try to make changes to the existing paradigm without modifying its large-scale properties. There are two routes in this direction: one is to prevent the formation of small dark matter haloes, the second is to hide them by suppressing star formation. Some modifications belonging to the first group consist of changing the nature of the dark matter itself to have different physical properties: self-interacting (Spergel & Steinhardt 2000), warm (Colombi, Dodelson & Widrow 1996; Bode, Ostriker & Turok 2001) or other variants (e.g. Goodman 2000; Hu, Barkana & Gruzinov 2000; Cen 2001). It has also been proposed that the solution to the problem may lie within the nature of the initial conditions and that some modification to the shape of the inflation potential may suppress small scale structure (Kamionkowski & Liddle 2000), although this solution faces some problems in reproducing the observed amount of structure in the Lyman-α forest. On the other hand, the solution of suppressing star formation via supernova feedback was suggested a long time ago (e.g. White & Rees 1978; Kauffmann, White & Guiderdoni 1993). This process expels the gas from the halo and therefore suppresses star formation. The accretion and cooling of gas can also be suppressed in the presence of a strong photoionizing background. This route has also been known for a number of years now (Doroshkevich, Zel'Dovich & Novikov 1967) and has been investigated recently in the context of semi-analytic models of galaxy formation (Benson et al. 2001a; Bullock, Kravtsov & Weinberg 2001b; Chiu, Gnedin & Ostriker 2001; Somerville 2002).

A closely-related problem is the fact that the observed faint end of the luminosity function has a shallower slope (Blanton et al. 2001) than predicted from high-resolution LCDM simulations, if one simply assumes a constant mass-to-light ratio (e.g. Jenkins et al. 2001). Stellar feedback is again often used in semi-analytic models to suppress star formation in haloes with shallow potential wells (e.g. Kauffmann et al. 1993; Cole et al. 2000).

All of the above models seek to increase systematically the mass-to-light ratio in small haloes, by depleting the cold-gas fraction available for star formation. In this paper we show that such drastic expulsion of gas is, perhaps, unnecessary: provided that the gas conserves angular momentum during collapse, a large fraction of low-mass haloes will be Toomre-stable (Toomre 1964; Kennicutt 1989), i.e. Q=Σ_c/Σ > 1 (where Σ is the disc surface mass density and Σ_c is the critical surface density to trigger gravitational instability and thus start formation). If the baryon fraction (f_d) in discs is f_d∼ 0.5 Ω_b/Ω_m, all discs in haloes with masses M < 10⁹ M_⊙ will be Toomre-stable and fail to form stars. The ‘missing’ dwarf galaxies predicted to be around the Milky Way and M31 may fall in this category. We also show that, if only a small fraction of baryons initially present in the dark matter halo settle in a disc, as some observational evidence suggests (Fukugita, Hogan & Peebles 1998; Burkert 2000; van den Bosch, Burkert & Swaters 2001; Guzik & Seljak 2001; Jimenez, Verde & Oh 2002), even-more-massive haloes have some probability of remaining dark, provided they have high spin. This could potentially reconcile the observed luminosity function with the mass function predicted in high-resolution LCDM N-body simulations. We show that the predictions of this model are consistent with our previous analysis of the rotation curves of galaxies (Jimenez et al. 2002 hereafter JVO02): none of the galaxies lies in the region of parameter space forbidden by the Toomre instability criterion.

Throughout this paper we assume an LCDM cosmology given by: (Ω_m,Ω_Λ,Ω_b,h,σ₈h⁻¹) = (0.3,0.7,0.039,0.7,1.0).

2 The mass distribution of dark galaxies

We model disc galaxies as exponential discs embedded in a cold dark matter halo. Following Dalcanton, Spergel & Summers (1997); Jimenez et al. (1997); Mo, Mao & White (1998), we assume an NFW (Navarro, Frenk & White 1997) profile ρ∝[(r/r_c)(1+r/r_c)²]⁻¹ where r_c is the break radius and the surface mass density of the disc is given by Σ=Σ₀ exp(−r/R_d), where R_d is the scalelength of the disc.

For a baryonic disc to be locally gravitationally unstable, despite the stabilizing effects of tidal shears and pressure forces, we require the Toomre parameter, Q, to be less than 1, where  

1

Σ is the disc surface mass density, κ= 1.41 (V/r)(1 + d ln V/d ln r^1/2 is the epicyclic frequency and c_s is the gas sound speed.

Regions in local disc galaxies where Q>1 are observationally associated with very little star formation, indicating that the Toomre criterion is obeyed remarkably well, if a gas sound speed of c_s∼ 6 km s⁻¹ is assumed (Kennicutt 1989). Ferguson et al. (1996) noted that for one galaxy (NGC 6946) the agreement with the Toomre criterion was not perfect: there was no truncation of star formation for Q > 1.6 if the observed vertical velocity dispersion was adopted. However, agreement with the Toomre criterion could be obtained if the velocity dispersion was assumed to be constant with radius. Wong & Blitz (2002) have measured the value of Q in the inner part of rotation curves for a sample of eight galaxies. For half of the sample the predicted value from Toomre theory for the abundance of gas is in excellent agreement with observations. In the other cases the agreement is not as remarkable, but still there is no evidence for large deviations from the Toomre stability predictions.

Jimenez et al. (1997)developed a disc model to study the dependence ofQon the mass, radius and spin parameter of the dark matter halo and found that, for a certain range of the above parameters, some discs would be dark. Their disc model also assumed aNavarro et al. (1997)NFW profile for the dark matter and that the baryon component settles into a disc until it is rotationally supported, conserving angular momentum. They then computed the final equilibrium configuration of the baryonic disc, assuming circular motion, and the profile of the final surface density of the gas, which they then compared with the Toomre criterion. An important ingredient in their study was the assumption of a constant baryon fraction, namely they adopted the nucleosynthesis value fromWalker et al. (1991) . Using this constant disc fraction for all galaxies, they concluded that for an Einstein–de Sitter universe the minimum dark matter halo mass for a visible disc galaxy would be about10⁹ M_⊙ . On the other hand, for a low matter-density universe (e.g. LCDM), this minimum mass drops substantially, by about two orders of magnitude, and thus the fraction of dark haloes also falls drastically. The reason for this is that, for a low-density universe with baryon density fixed by nucleosynthesis, the baryons constitute a higher fraction of the halo mass.

It seems unlikely that all baryons initially present in a dark matter halo will be able to cool and form stars. Recent studies from galaxy–galaxy lensing from the Sloan Digital Sky Survey (SDSS) show that for L* galaxies the baryon fraction is about half of the nucleosynthesis value (Guzik & Seljak 2001). On the other hand, theoretical considerations suggest that the baryon fraction for galaxies, either less massive or more massive than L_*, is likely to be much smaller due to feedback and inefficient cooling, respectively (e.g. Benson et al. 2001b; van den Bosch 2002). This implies that a large fraction of discs with masses ≲10¹² M_⊙ should be Toomre-stable. Below, we quantify how this affects the predicted number density of observable galaxies.

In Section 3 we show that even assuming a radically different dark matter halo profile (a truncated isothermal sphere, which we refer to as a pseudoisothermal sphere) very similar regions in the galaxy parameters space are Toomre-stable and thus give rise to dark galaxies. From this consideration we can safely conclude that our results will not be very sensitive to the particular choice of the dark matter profile.

2.1 Luminosity function

CDM models typically predict many more galaxies at the faint end of the luminosity function than are actually observed. This can most easily be seen by comparing the galaxy circular velocity function with theoretical predictions. The galaxy velocity function is obtained by combining the observed luminosity function with empirically determined luminosity–velocity relations (the Tully–Fisher relation for spirals, and the Faber–Jackson relation for ellipticals). Gonzalez et al. (2000) find that treating the entire galaxy population as spirals does not significantly alter the derived velocity function. In Fig. 1 we plot a representative velocity function from their paper, derived from the SSRS2 (Marzke et al. 1998) B-band survey and the Yasuda, Fukugita & Okamura (1997) Tully–Fisher relation. We also construct a new velocity function from the Sloan Digital Sky Survey (SDSS), which has published luminosity functions in five bands (Blanton et al. 2001). We choose to use the R band, as the R-band Tully–Fisher relation is tighter and has less scatter than the Tully–Fisher relation at shorter wavelengths. Blanton et al. (2001) perform the colour transformations from SDSS magnitudes to LCRS (Shectman et al. 1996) R-band magnitudes (see their table 3); we perform the colour transformations and corrections for internal extinction as in Gonzalez et al. (2000) to the R_courteau magnitudes used by Courteau (1997). We then use the Courteau (1997),R_courteau-band Tully–Fisher relation to obtain the velocity function. Note that the Sloan velocity function has a significantly higher normalization than the SSRS2 velocity function. This is due to the higher sensitivity of the SDSS, which allows detection of galaxies of significantly lower surface brightness. The luminosity density of galaxies in the SDSS is therefore significantly higher than other redshift surveys. The SDSS luminosity function is very similar to the luminosity functions derived by other surveys such as the Las Campanas Redshift Survey (LCRS) and Two Degree Field Galaxy Redshift Survey (2dFGRS) if they re-analyse their data using the shallower isophotal limits for galaxy magnitude employed by these surveys; see Blanton et al. (2001) for a discussion.

1

The present-day velocity function. Lines labelled by Sloan and SRSS2 are obtained from the SDSS and the SRSS2 surveys, respectively. The solid line labelledf_d = 0.06is the theoretical velocity function derived from the Sheth & Tormen (1999) mass function using a constant disc–baryon fractionf_d = 0.06. The line labelled NDG (no dark galaxies) is obtained assuming our Ansatz for the dependence off_don mass (see Fig. 2 ). Finally, the thick solid line labelled DG (dark galaxies) is obtained considering that only discs that are Toomre-unstable will be able to form stars and thus be visible (see text for more details).

The theoretical velocity function is given by  

2

where V_c is the circular velocity in the flat part of the rotation curve and f_bright is the fraction of haloes that form stars. For the mass function dN/dM, we use the Sheth & Tormen (1999) mass function, which gives a good fit to the mass function seen in the Hubble Volume simulations (see Jenkins et al. 2001). For evaluating the derivative dM/d log(V_c), we use the fitting formula of Mo et al. (1998), which gives the circular velocity V_c as a function of the virial mass of the dark halo (M₂₀₀), the concentration parameter (c), the spin of the dark halo (λ) and f_d. As in Navarro et al. (1997), for a given M₂₀₀ we evaluate the concentration parameter from the typical collapse redshift. Below, we discuss further our choice of (f_d, λ) and evaluate f_bright using Toomre instability arguments.

Let us start by assuming f_bright= 1, a constant disc mass fraction f_d= 0.06 (i.e. about half the nucleosynthesis value of Fukugita et al. 1998), and the disc spin parameter to be the mean value found in N-body simulations λ_mean= 0.042 (Bullock et al. 2001a). In Fig. 1 we plot the theoretical velocity function (upper solid line labelled by f_d= 0.06). At low velocities, the predicted number density of galaxies is much larger than observed, which is the discrepancy we seek to resolve.

This cannot be resolved by making f_d a function of M₂₀₀, unless the f_d(M₂₀₀) dependence is unrealistically drastic. To illustrate this we assume a different functional form for f_d(M₂₀₀), as shown in Fig. 2. Our Ansatz for the f_d(M₂₀₀) relation models a fairly mild factor of 2 decline in the value of f_d from high-mass to low-mass haloes, which fits the median value of f_d at high masses (Mo et al. 1998) and low masses (van den Bosch & Swaters 2001). If this were the only process at play it would produce only a mild change in the theoretical velocity function, as shown by the thin solid line labelled by NDG (no dark galaxies) in Fig. 1.

2

Our Ansatz for the mass dependence of the disc–baryon fractionf_d. This models a factor of 2 decline in the value off_dfrom high-mass to low-mass haloes, which fits the median value off_dat high masses ( Mo et al. 1998 ) and low masses ( van den Bosch & Swaters 2001 ).

To include the effects of the Toomre instability criterion, we proceed as follows. We define a critical spin parameter λ_crit such that the disc is marginally stable: Q[f_d(M₂₀₀), λ_crit, M₂₀₀, c(M₂₀₀)]= 1.5, where Q is evaluated at one disc scalelength, and we assume a gas sound speed c_s= 6 km s⁻¹ (the gas will in general be photoionized by the meta-galactic ionizing background to T∼ 10⁴ K). The fraction of bright galaxies is given by , where p(λ) is the probability distribution of λ derived from numerical simulations (Bullock et al. 2001a). Note that since haloes with λ > λ_crit are assumed not to host visible galaxies, the probability distribution of λ in observed galaxies is skewed. We evaluate the circular velocity assuming a median value λ_median, defined by the implicit equation . The velocity function we obtain from computing f_bright in this manner is shown in Fig. 1 by the thick solid line labelled by DG (dark galaxies). Thus, once we include the suppression of high-spin galaxies by explicitly computing f_bright, the theoretical predictions match the observations well.

The suppression of accretion after reionization (Benson et al. 2001a; Bullock et al. 2001b; Chiu et al. 2001; Somerville 2002; Tully et al. 2002), is an attractive mechanism for regulating f_d since (unlike, for instance, supernova explosions) it affects dark galaxies as well. We model this effect by using the fitting formula (Gnedin 2000; Somerville 2002):  

3

where we approximate M_C as the mass corresponding to a halo with virial velocity V_c= 50 km s⁻¹ (Somerville 2002) at the redshift of formation z_f. For a given halo mass M₂₀₀, extended Press–Schechter theory (Lacey & Cole 1994) gives the probability distribution of collapse redshifts p(z_f); we can then use equation (3) to obtain the probability distribution p(f_d). Note that f_d=Ω_b/Ω_m for z > z_reion, where z_reion is the redshift of reionization; we assume z_reion∼ 7 (the results depend only very weakly on z_reion).

In the photoionization model, the fraction of dark galaxies can be computed via equation (4), with the probability distribution p(f_d) given above. Significant suppression of the number density of low-circular-velocity galaxies only takes place off the scales plotted in Fig. 1. Photoionization only successfully suppresses gas accretion in very small haloes; it succeeds in alleviating the dwarf satellite problem, but leaves unaffected haloes of V_c > 75 km s⁻¹. By itself it therefore cannot produce the low values of f_d required to reconcile the theoretical and observed velocity function.

We do not attempt to model the f_d(M₂₀₀) relation in more detail (which is the task of hydrodynamic galaxy formation modelling), but merely note that it is consistent with observations and theoretical prejudice that gas expulsion and suppression of accretion are more important for shallow potential wells. Possible mechanisms include suppression of accretion due to entropy injection in the intergalactic medium by supernovae and/or AGN jets, or ram pressure stripping as a smaller halo falls into a larger halo. We note merely that by invoking the Toomre instability criterion we do not require the very strong variation of f_d with M₂₀₀ required in most semi-analytic models of galaxy formation. In particular, we do not require a very sharp drop in f_d at the faint end.

2.2 The fraction of dark galaxies

In Fig. 3 we plot the fraction of dark galaxies f_dark as a function of halo mass.

3

Fraction of dark galaxiesf_darkas a function of halo mass. The dashed line isf_dark = 1 − f_bright, wheref_brighthas been computed in Section 2 , by considering the Toomre instability criterion and a mild mass dependence of the disc–baryon fractionf_dshown in Fig. 2 . The dotted line shows the resultingf_darkif suppression of accretion due to reionization determinesf_d. The solid line is the dark fraction derived using the distribution off_dempirically obtained by fitting the rotation curves.

The dashed line is f_dark= 1−f_bright, where f_bright has been computed in the previous section. The dotted line shows the resulting f_dark if photoionization was the only process at play. As already noted this mechanism only suppresses accretion in haloes of low circular velocity and cannot reconcile the observed and theoretical velocity function. Assuming f_d as a deterministic function of M₂₀₀ might be an oversimplification: in reality f_d is likely to be highly stochastic (e.g. Benson et al. 2001b; van den Bosch 2002). Also, analysis of rotation curves (JVO02) confirms that there is a significant scatter in the relation (see Section 3). We thus compute also the expected fraction of dark galaxies by neglecting any possible mass dependence and using the distribution of f_dp(f_d) empirically obtained by fitting the rotation curves for the NFW profile (see JVO02 and Section 3):  

4

where λ_crit and f_{d crit} are defined by Q(M₂₀₀, f_{d crit}, λ_crit, c) = 1. The result is shown in Fig. 3 as a solid line. The two approaches (the one described in Section 2.1 or the one using the empirically-obtained f_d distribution) yield very similar values for f_dark, indicative of the fact that the main mechanism creating dark galaxies is the Toomre criterion.

2.3 Dwarf satellites

From Fig. 3 it is clear that star formation within dark matter haloes with masses below 10¹⁰ M_⊙ will be almost completely suppressed (not visible because they will be Toomre-stable). Thus this provides a mechanism for hiding these low-mass haloes. To make a more quantitative estimate of the viability of this effect, in Fig. 4 we compare the observed cumulative mass function of dwarf satellites with the predictions from N-body simulations. The solid line is from Moore et al. (1999), where the velocity has been converted in mass.1 The grey area is an estimate of the error from fig. 2 of Moore et al. (1999). The solid and dashed lines are the predicted cumulative number of satellites once the probability that some galaxies will remain dark has be taken into account. In particular we have used the predictions of the solid and dashed lines of Fig. 3, respectively. This prediction agrees nicely with observations (triangles). We have not shown results for masses below 10⁹ M_⊙ because the uncertainties in the procedure for obtaining the corrected line become too important; only a very small fraction (<1 per cent) of dark matter haloes will be visible thus making this fraction very dependent on the exact value of the predicted number of dark galaxies.

4

Cumulative mass function of dwarf satellites of the Milky Way. Small triangles are the observations, the solid line with the grey area is the prediction from N -body simulations with an estimate of the uncertainty. Finally, the solid and dashed lines show the predictions once the the probability that some haloes will remain dark (owing to their spin parameter being above the critical value) has been taken into account. Solid and dashed lines have been computed using thef_darkfraction of Fig. 4 .

3 Rotation Curves

In this section we present some indication of the existence of dark galaxies from analysis of a large sample of rotation curves.

We use the same set of rotation curves as in JVO02, where we determined the best-fit disc parameters within the context of both a disc within an NFW profile halo and within a pseudo-isothermal profile halo (see JVO02 for details).

The two sets of observed rotation curves are as follows.

• The Courteau catalogue (Courteau 1997), which consists of optical (H_α) long-slit rotation curves for 300 Sb–Sc UGC galaxies.

• The 64 spiral galaxies (Sa–Sd) observed in H_α by Palunas & Williams (2000).

In total we use 364 optical rotation curves (to avoid any problems with the smearing of the beam in radio observations, see van den Bosch & Swaters 2001) corrected for inclination and with error bars. For these galaxies the scalelength of the exponential disc (R_d) is also measured. The two models have three free parameters: M₂₀₀, r_c and f_d. In principle there is another disc parameter: the modified spin parameter, λ′. If λ is the spin parameter of the dark halo, then λ′=j_d/f_dλ, where j_d≡J_d/J and J_d,J are the total angular momentum of the disc and halo, respectively. λ′ can be readily obtained using the measured value of R_d and the best-fit parameters (see JVO02). If angular momentum is conserved, then λ=λ′, i.e. j_d/f_d= 1.

In JVO02 we used the sample of 364 galaxies for which R_d is measured, we found the best-fitting parameters (M₂₀₀, f_d and c) for the above two models and recovered λ′. Within both models for the majority of the galaxies, a meaningful set of fitting parameters was recovered. For the remaining galaxies we obtained values for f_d that where either implausibly low (f_d≃ 0.001) or may be too high (f_d≳ 0.2).

3.1 Fit to the disc models

As described in detail in JVO02, for each galaxy rotation curve and for both models, the best-fitting parameters were obtained using a standard χ² minimization and exploring the whole likelihood surface. The sample of galaxies studied is very non-uniform: some curves have much smaller error-bars than others, some curves do not extend to large enough radii to show the flattening of the rotation curve, while others show strong evidence of spiral arms and bars in the rotation pattern, moreover we have not attempted to model any bulge or bar component. Thus it is important to keep in mind that galaxies might be more complicated than our model for galaxy rotation curves. However, for most of these galaxies the high quality of the rotation-curves measurements allows degeneracies among the disc parameters to be lifted (JVO02). Here, we are not interested in a detailed modelling of individual galaxy dynamics but rather the general statistical trends of the recovered disc parameters from the whole sample.

Fig. 5 shows the distribution of the best-fitting disc parameters in thef_d −λ′ andM₂₀₀ − f_d planes, for the NFW profile (top) and pseudoisothermal profile (bottom); diamonds correspond to galaxies from the Courteau (1997) sample and triangles are from Palunas & Williams (2000). In JVO02 we showed that low-surface-brightness galaxies (LSBs) follow the same trend as the rest of the sample, so any observational bias against LSBs does not seem to affect this result: the only difference is that LSBs have systematically higher values of λ′ for a fixed value off_d. In both models there is a correlation between λ′ andf_d, albeit with a scatter, and a weak anti-correlation betweenM₂₀₀ andf_d .

5

Distribution of the best-fitting disc parameters in thef_d −λ′plane and in thef_d − M₂₀₀plane for the NFW (Top) and pseudoisothermal profile (Bottom). We argue that the Toomre instability criterion correctly predicts the zone of avoidance (grey-shaded areas) for combinations of disc parameters that involve lowM₂₀₀and lowf_d. Note that there seems to be another small void region in the bottom-right corner of thef_d −λ′plot. This might arise because for these parameter values discs might be too concentrated and form a spheroid instead of a disc. The dotted line on the upper-left panel corresponds to the stability region under global instabilities (see text). Galaxies above are stable, thus dark galaxies are in the stable region.

For f_d < 0.02, λ′ values are between 0.01 and 0.03, while for f_d > 0.02, λ′ takes all values between 0.02 and 0.4.

The origin of this correlation is controversial: van den Bosch (2002) argues that when feedback by star formation is included in models of galaxy formation, such a correlation between λ′ and f_d arises naturally for low-mass galaxies, where feedback is most efficient; Burkert (2000) showed that it can be an effect of disc growth happening inside-out; or discs might preferentially lose baryons with high angular momentum. On the other hand, Burkert, van den Bosch & Swaters (2002) argue that the correlation may be entirely due to intrinsic degeneracies due to the fact that from rotation curves analysis the total dark halo mass is not well constrained. In JVO02 (fig. 14) we find that this is the case for rotation curves that are not measured at large enough radii, but, for most of the galaxies in the sample (70 per cent), the high-quality rotation curves allow model degeneracies to be lifted. This indicates that the random error in recovered parameters is small. On the other hand, the systematic errors (especially due to the uncertainty in the dark matter profile) may be larger, possibly causing the correlations.

We argue here that the most significant feature of the left panels of Fig. 1is not the correlation itself, but the fact that spirals avoid the upper-left portion of the λ′−f_d plane and that this is in agreement with the Toomre instability criterion. In passing we also note that disc instability might make galaxies with disc parameters in the bottom-left corner of the plot drop out of the sample: for these parameter values, discs might be too concentrated and form a spheroidal instead of a disc, or galaxies might become bulge dominated.

The right panels show an interesting, although weak, correlation between M₂₀₀ and f_d: low-mass spirals (M₂₀₀ < 10¹¹ M_⊙) tend to have f_d > 0.1, since no galaxies are found in the lower-left portion of the plot. There can be several different effects that affect the dependence of f_d on the halo mass. Linear theory evolution from inflation predicts a weak anti-correlation between spin parameter and peak height (Heavens & Peacock 1988), and hence mass of the halo. The expected anti-correlation is, however, very weak. On the other hand, cooling is less efficient in massive haloes (M₂₀₀≳ 10¹² M_⊙), thus f_d is expected to decrease strongly with increasing M₂₀₀ (e.g. van den Bosch 2002) at the high-mass end. However, in small-mass haloes, reheating of cooled gas through feedback and photoionization processes are much more efficient than in massive systems, thus implying that f_d should increase with M₂₀₀ (e.g. Benson et al. 2001b) at the low-mass end. The correlation we find from rotation curves analysis is in rough agreement with the high-mass end trend (see fig. 5 of van den Bosch 2002). For his no-feedback model, the baryon fraction decreases for increasing masses. For the feedback model this is only true for M₂₀₀ bigger than ∼5 × 10¹¹ M_⊙. For lower masses this trend reverses and we do not find evidence for it from our analysis of rotation curves. This may be due to the fact that, for these low masses, only those galaxies with highest f_d are visible.

Here we show that Toomre-instability-criterion predictions are in agreement with the zone of avoidance for combinations of disc parameters that involve low M₂₀₀ and low f_d or low λ.

3.2 Stability of discs

The void areas found in Fig. 5 can be explained by resorting to the Toomre instability criterion: for Q > 1 the gas disc is stable and thus star formation is strongly suppressed. Exressions for Q and Σ for the NFW profile have been presented in Jimenez et al. (1997); here we derive the corresponding expressions for the pseudoisothermal profile.

Recall that in the pseudoisothermal profile the density ρ as a function of the radius r is ρ(r) =ρ₀[1 + (r/r_c)²]⁻¹, where ρ₀ denotes the central finite density. The disc surface density as a function of the radius is Σ(r) =Σ₀ exp(−r/R_d), where Σ₀=f_dM₂₀₀/(2 πR_d²). Appendix I of JVO02 gives a relation between ρ₀ and the concentration parameter c and a fitting formula for R_d as a function of the disc parameters: M₂₀₀, ρ₀, r_c and λ′.

The shaded areas in Fig. 5 represent the void regions predicted by the Toomre criterion assuming c_s= 6 km s⁻¹ for all values of masses and radii. The dashed line shows how the avoidance region would change if c_s= 2 km s⁻¹. Such a low value is based on the findings of Corbelli & Salpeter (1995), who argue that when the quasar UV ionizing background is the only source of ionization of the H _i gas, as should be the case for dark galaxies with no supernova, then the sound speed of the gas cannot be higher than 2 km s⁻¹. Measured values of the sound speed in the Milky Way and in nearby galaxies (see e.g. Kennicutt 1989) range between 3 and 10 km s⁻¹. Note that using a low value for c_s is a conservative assumption, since for small values of the velocity dispersion the formation of dark galaxies is more difficult.

More specifically, for a grid in M₂₀₀, λ′, f_d, we compute the surface density of the exponential disc Σ and whether Q>1. The shaded area in the left panels of Fig. 5 are obtained assuming a fixed mass of 5 × 10¹¹ M_⊙, approximately in the middle of the mass range considered. The shaded area in the right panels is obtained as follows: we assume a log-normal distribution of spin parameters with mean 0.042 and dispersion 0.5 (Bullock et al. 2001a), thus 68 per cent of galaxies are expected to have λ′>λ′₆₈≡ 0.033. We then assume that a combination of f_d and M₂₀₀ parameters corresponds to a dark galaxy when for all λ′>λ′₆₈ the disc is Toomre-stable. The Toomre criterion seems to work better for the NFW profile than for the pseudoisothermal.

However, for the pseudoisothermal profile model there is no theoretical prediction for the expected λ′ distribution and the dependence of the concentration parameter c≡R₂₀₀/r_c on the halo mass. To produce the grey areas in Fig. 5 we assume the same log-normal distribution for λ as predicted from CDM simulations. This assumption should be valid since JVO02 find that the empirically recovered λ′ distribution is well approximated by the CDM prediction. For setting the concentration parameter we consider that, when a galaxy rotation curve is fitted by the pseudoisothermal profile, r_c obtained is about 1/10 of that obtained when the curve is fitted by the NFW profile (fig. 5 of JVO02). Thus, for a given mass, we compute the concentration parameter as expected in the NFW profile, and convert it to the pseudoisothermal one. The scatter around this relation is, however, quite large, thus introducing some uncertainty in the grey areas so obtained.

The Toomre criterion is a local instability criterion; global instabilities can in principle trigger star formation. For example, discs in which the self gravity is dominant are likely to be unstable to the formation of a bar. Here we use the results of Efstathiou, Lake & Negroponte (1982) and Mo et al. (1998) to conclude that the onset of these instabilities does not affect the grey areas. The dotted line in Fig. 5 shows the transition between stability and instability according the above criteria. Discs below the dotted line are unstable under global instabilities, while the ones above are stable. The grey area is entirely in the stable region. The dependence of the line position on concentration parameter is very weak (see Mo et al. 1998).

We conclude that the Toomre criterion, at least in the context of an NFW profile, seems to predict correctly the region of the parameter space that galaxies tend to avoid (i.e. these galaxies should be dark).

4 Discussion and Conclusions

In this paper we propose that gas in a large fraction of low-mass dark matter haloes may form Toomre-stable discs. Such haloes would be stable to star formation and therefore remain dark. This may potentially explain the deviation between the predicted and the observed faint-end slope of the luminosity function as well as the discrepancy between the predicted and observed number of dwarf satellites in the Local Group. We show that model fits to rotation curves are consistent with this hypothesis: none of the observed galaxies lies in the region of parameter space forbidden by the Toomre stability criterion. Such Toomre-stable discs may be the origin of damped Lyα systems seen at high redshift – see, e.g., Kauffmann (1996); Jimenez et al. (1998) and Mo et al. (1998) (especially their section 4).

Most semi-analytic models of galaxy formation achieve a reconciliation between the observed and predicted abundance of low luminosity galaxies by drastically decreasing f_d for faint galaxies. At present, there is no evidence from rotation curve modelling that low-circular-velocity discs are dark matter dominated (which would be the case if f_d were very small). Our model only requires a very gentle decline in f_d toward low masses. It magnifies the effect of previously proposed mechanisms (e.g. supernovae feedback, suppression of accretion, ram pressure stripping, etc.), since Q∝f_d⁻¹ and low-f_d discs are more likely to be Toomre-stable; such mechanisms can therefore be ‘tuned down’ to lower levels. This may be relevant to the problem of trying to fit the luminosity function and the Tully–Fisher relation simultaneously. In our model, the mass-to-light ratio is much more stable as a function of halo mass (except in dark haloes, where it is infinite). Indeed, by construction we fit both the Tully–Fisher relation and the luminosity function, since both were used to obtain the velocity function.

Jimenez et al. (1997)proposed that discs in high-spin haloes may not be low-surface-brightness galaxies but instead be completely dark. This was at variance with the proposal that a large fraction of the ‘missing galaxies’ lies in low-surface-brightness galaxies that formed preferentially in high-spin haloes ( Dalcanton et al. 1997 ). Here, we have extended theJimenez et al. (1997) model to an LCDM scenario and show that it is supported from an analysis of rotation curves. If it is true that discs in such high-spin haloes may not be low surface brightness but instead be completely dark, then there should be a cut-off in the surface brightness distribution of galaxies; as surveys probe increasingly-low-surface-brightness galaxies, the luminosity density will not continue to increase but instead will plateau.

It is important to consider independent indications for the existence of this dark-galaxy population and devise possible observational tests. For example, Augusto & Wilkinson (2001) have investigated the detectability of dark objects.

These dark galaxies have very-low surface densities of H _i (below 10²⁰ cm⁻² at one scalelength of the disc), still high enough to be detected by current surveys (Zwaan 2001) if the H _i is in its neutral phase. Current observations that are sensitive enough to detect such amounts of neutral H _i have been carried out recently by Charlton, Churchill & Rigby (2000), Zwaan & Briggs (2000) and Zwaan (2001). In particular, Zwaan (2001) found no significant detection of H _i, down to a detection limit of 7 × 10⁶ M_⊙, not associated with visible galaxies. Therefore, it seems unlikely that large masses of neutral H _i are harboured in dark galaxies. This should not come as a surprise, since, given the low surface densities, the hydrogen might be ionized by the extragalactic background. At densities below 10¹⁹ cm⁻², H _i is not able to self-shield from the external radiation field and most of it will be in its ionized phase. Another possibility is that, owing to pressure ram stripping from the hot gas of the hosting large galaxy, most of the H _i gas could be removed from the dark galaxy (Quilis, Moore & Bower 2000). Blitz & Robishaw (2000) have shown that even for moderate amounts of hot gas in the haloes of both M31 and MW, the ram pressure is very efficient at swiping out the atomic gas in the ‘satellites’.

On the other hand, Blitz & Robishaw (2000) have shown that 10 of the 21 dwarf spheroidal galaxies in the Local Group are associated with large reservoirs of atomic gas. This is in good agreement with our predictions, where it is not discarded that the small central region may suffer some star formation (see below), but the majority of the disc remains stable to star formation. Mayer et al. (2001) show that small satellites entering a galaxy halo will develop global instabilities. In particular, they show that low-surface-brightness and high-surface-brightness galaxies will develop into dSph's and dE's, respectively. Thus the strengh of the instability seems to depend on how concentrated the galaxies are. The majority of our galaxies will have densities lower than those considered in their paper.

A more promising route to detect the presence of dark galaxies is through their gravitational lensing properties. Recently, Metcalf & Zhao (2002), Dalal & Kochanek (2002), Keeton (2001) and Bradac et al. (2001) have discussed the need of substructure to explain the relative fluxes of multiple-lensed systems. They show that the flux ratios observed in each lens can only be explained by the presence of substructure within a large smooth halo. This population would be similar to our dark galaxies. Dark galaxies with dark-halo masses above 10¹⁰ M_⊙ should be directly visible via strong gravitational lensing. They will produce splitting angles of ≲0.1 arcsec, which is just at the detection limit of current technology.

Using the metallicity distributions of globular clusters, Cote, West & Marzke (2002) conclude that the mass spectrum of proto-galactic fragments for the galaxies in their sample is a power law with index ∼−2, indistinguishable from that predicted from N-body simulations. They argue that the missing satellite population must therefore belong to a class of dark galaxies, similar to the ones we consider here.

Stoehr et al. (2002)have shown that the kinematics of the 11 known satellites of the Milky Way are in excellent agreement with those from high-resolutionN -body LCDM simulations and can be identified with many of the most massive simulation satellites(M≫ 10⁹ M_⊙) . Thus they argue that the rest of the satellites in their simulation must be inefficient in making stars. This agrees well with our findings.

The biggest caveat in our model is the assumption of conservation of angular momentum. In principle, angular momentum loss is possible owing to torquing by substructure in the dark matter halo; indeed, it is invariably seen in SPH simulations of disc galaxy formation. This results in more compact discs that are much more susceptible to gravitational instability. However, both fits to the scalelengths of discs (Mo et al. 1998) and modelling of the rotation curves (van den Bosch et al. 2001; Jimenez et al. 2002) recover spin parameters consistent with the log-normal distribution seen in dissipationless N-body simulations, indicating little or no loss of (the modulus of) angular momentum. Moreover, van den Bosch et al. (2002) show that the spin and angular momentum distribution of dark matter haloes and discs are very similar. The angular momentum vectors are poorly aligned, but the effect on the surface density of the discs (which is what our model depends on) is very small.

In addition, we have ignored the effect of mergers; as haloes merge, their discs are likely to be disrupted and transformed into spheroidal systems, which again are more likely to become self-gravitating. The merger history of haloes is well handled in semi-analytic models of galaxy formation that employ Press–Schechter-based Monte Carlo merger trees. It would be very interesting to incorporate the additional physics of Toomre stability into pre-existing semi-analytic (or hydrodynamic) galaxy formation models to compare its importance against various other proposed schemes for suppressing star formation in low-mass haloes.

Acknowledgments

We thank David Spergel for insightful comments and stimulating discussions. We are also grateful to Fabio Governato (the referee) and Simon White for insights and useful discussions. LV is supported in part by NASA grant NAG5-7154. SPO is supported by NSF grant AST-0096023. LV and RJ thank the TAPIR group at Caltech for hospitality.

Footnotes

1

We have assumed that v_c³∝M_vir where v_c denotes the circular velocity (and normalized the result to the fact that there is only one LMC with mass ∼2 × 10¹⁰ M_⊙ in the local group). This assumption is not strictly valid because dark matter is stripped from satellites during infall. However, the mass function for satellites is expected to be n(M)∝M⁻², thus the cumulative mass function is N(M)∝M⁻¹. If dark matter stripping in the satellites is important (and plausibly more important for low-mass haloes than high-mass ones) then the cumulative number of satellites as a function of v_c should deviate significantly from N(v_c)∝v_c⁻³. The predicted N(v_c) from fig. 2 of Moore et al. (1999) is well approximated by a function proportional to v_c⁻³. Thus, for the level of accuracy of this argument we can neglect tidal stripping effects.