Which Molecular Cloud Structures Are Bound?

Neal J. Evans II; Mark Heyer; Marc-Antoine Miville-Deschênes; Quang Nguyen-Luong; Manuel Merello

doi:10.3847/1538-4357/ac1425

1. Introduction

The fundamental problem of star formation was first stated by Zuckerman & Palmer (1974) and Zuckerman & Evans (1974). Simply stated, a theoretical estimate of the star formation rate in the Milky Way exceeded by a factor of about 100 the rate averaged over the last Gyr or so, as inferred from observations. Those early estimates were based on preliminary surveys of the Galaxy in CO and estimates of the average star formation rate, but the situation has only worsened after nearly a half century of improved surveys. Heyer & Dame (2015) estimated a mass of (1.0 ± 0.3) × 10⁹ M_⊙ for the total molecular gas in the Galaxy. The freefall time,

$\begin{eqnarray}&&\begin{array}{r}\begin{array}{r}{t}_{{\rm{ff}}}={\left(\displaystyle \frac{3\pi }{32G\rho }\right)}^{0.5}=3.34\times {10}^{7}{n}^{-0.5}{\rm{yr}}\\ =1.66\times {10}^{7}{r}_{{\rm{pc}}}^{1.5}{M}_{{M}_{\odot }}^{-0.5}\,{\rm{yr}},\end{array}\end{array}\end{eqnarray} \tag{ 1 }$

where n is the particle density in cm⁻³, related to ρ by the mean mass per particle of 2.37 (Kauffmann et al. 2008). For a characteristic density of 100 cm⁻³, t_ff = 3.34 × 10⁶ yr, predicting a star formation rate of 300 M_⊙ yr⁻¹. In contrast, observational data indicate much lower rates. Counting young stars from the GLIMPSE (Benjamin et al. 2003; Churchwell et al. 2009) survey, Robitaille & Whitney (2010) constrained the star formation rate to 0.68–1.45 M_⊙ yr⁻¹, averaged over the lifetime of infrared excesses around young stars, roughly the last few Myr. Using a wider range of indicators, Chomiuk & Povich (2011) derived a star formation rate of 1.9 ± 0.4 M_⊙ yr⁻¹, which we adopt. The discrepancy has, if anything, worsened with further research to a factor of 160.

This problem arises because the molecular clouds identified by CO emission are assumed to be gravitationally bound and collapsing on a freefall timescale. This problem is properly called the slowness of star formation (Krumholz et al. 2014; Padoan et al. 2014) and parameterized by the efficiency per freefall time, _ff, which is taken to be 0.01, reflecting the factor of 100 problem identified by the first references, but the latest numbers suggest _ff = 0.006. Another expression of the slowness is that the depletion time, given by the molecular mass divided by the star formation rate, is t_dep = 5 × 10⁸ yr, much larger than t_ff. The depletion time is even longer in other galaxies, t_dep ≈ 2 Gyr (Bigiel et al. 2008).

A somewhat distinct problem is the final efficiency of star formation, _sf, defined by

$\begin{eqnarray}&&{\epsilon }_{\mathrm{sf}}={M}_{* }/{M}_{\mathrm{gas}},\end{eqnarray} \tag{ 2 }$

where M_* is the mass in stars at the end of the event and M_gas is the mass of gas at the beginning. While we cannot measure these for an individual cloud, star clusters that have formed reasonably recently (i.e., leaving aside globular clusters) are much less massive than the mass in clouds. We can compare the mass function of molecular clouds to that of clusters. The mass functions, expressed per logarithmic unit of mass, are expressed as

$\begin{eqnarray}&&\displaystyle \frac{{dN}}{d\mathrm{log}M}\propto {\left(M/{M}_{u,\ x}\right)}^{{{\rm{\Gamma }}}_{x}}\end{eqnarray} \tag{ 3 }$

where M is the mass of entity x, up to a maximum mass M_u,x. For Milky Way OB associations and clusters, Γ_cluster is −1.0 and M_u,cluster is about 6 × 10⁴ M_⊙ (Massey et al. 1995; McKee & Williams 1997; Portegies Zwart et al. 2010). For molecular clouds, the values are Γ_cloud = −0.6 and M_u,cloud = 6 × 10⁶ M_⊙ (Williams & McKee 1997). See also Rosolowsky (2005), who gives results in terms of γ = Γ_cloud + 1, with our definition of Γ_cloud. Heyer et al. (2001) found a slightly more negative value of −0.80 ± 0.03 in the outer Galaxy. The ratio of upper mass limits suggests that _sf = 0.01, but the differences in the slopes suggests that molecular clouds are not the immediate precursors of clusters and associations. Dense clumps, in contrast, have a similar slope and upper mass limit to those of clusters and associations. Shirley et al. (2003) found Γ_clump = −1.0 for dense clumps defined by strong CS J = 5 → 4 emission. Elia et al. (2017) found Γ_clump = −1.05 ± 0.08 in a sample of 10⁴ clumps defined from submillimeter continuum emission. The upper mass limit for clumps, M_u,clump is about 1 × 10⁵ M_⊙ based on various studies (Shirley et al. 2003; Svoboda et al. 2016; Urquhart et al. 2018).

Much theoretical ingenuity has been expended to solve the problem of slow star formation. Magnetic fields may help to "support" clouds, and turbulence can slow collapse. However, numerical simulations generally found that supersonic turbulence should decay rapidly, with or without magnetic fields (Stone et al. 1998; Mac Low 1999). Also, turbulence appears to be driven on scales larger than the cloud, rather than from within the cloud (Larson 1981; Heyer & Brunt 2004; Brunt et al. 2009).

The kind of solution that invokes turbulence to support a cloud has the problem of the inherent instability of a gravitationally bound entity. It is hard to imagine an individual cloud avoiding gravitational collapse for a depletion time that is over 200 times the freefall time. Indeed, Heyer & Dame (2015) placed a strong upper limit on cloud lifetime (t_cloud) of 10⁸ yr and preferred a lifetime of 3 × 10⁷ yr, based on Kawamura et al. (2009), closer to 10 times the canonical value for t_ff. Clearly most of the mass in clouds must disperse after a star-forming event. The low value of _sf suggests that 99.4% is dispersed. Another approach has been to question the use of the current density to estimate the freefall time (Vázquez-Semadeni et al. 2019). A lower initial density would slow the initial rate of collapse, but the conflict with the prediction based on current properties of molecular clouds remains.

Recent theoretical work has concentrated on feedback to remove most of the gas before it forms stars as some star formation either regenerates turbulence or dissipates the cloud (Offner & Liu 2018; Murray 2011). However, the low efficiency and slowness of star formation in clouds forming only low-mass stars are similar to those in clouds forming massive stars, where the feedback is much stronger (Vutisalchavakul et al. 2016). Also, most simulations find that feedback has difficulty accounting for the very low efficiencies (e.g., Kim et al. 2018; Grudić et al. 2019). For example, Kim et al. (2018) compared the effectiveness of radiation pressure and photoionization in disrupting clouds. Even though they started with marginally bound clouds, the final efficiencies, _sf, increased strongly with increasing initial cloud surface density (Σ_M), exceeding 5% by initial Σ_M = 20 M_⊙ pc⁻² (Figure 7 of their paper). Bound clouds overproduce stars even at very low initial Σ_M, even with strong radiative feedback. Another interesting feature of Figure 7 of Kim et al. (2018) is that _sf reaches about 0.5 at Σ_M ∼ 10³ M_⊙ pc⁻², similar to the surface density of the dense clumps studied by Shirley et al. (2003). If these were the bound entities, the ratio of maximum cluster mass to initial gas mass would come out about right. Simulations of unbound clouds look more promising. Kim et al. (2021) have performed radiative magnetohydrodynamical simulations (RMHD) for clouds with starting Σ_M = 80 M_⊙pc⁻², reasonably appropriate to the structures traced by CO. The star formation is limited by feedback from massive stars. They found that both _ff and _sf decrease with increasing initial virial parameter (α_vir,0), reaching values similar to the observed values for α_vir,0 ∼ 5.

In the current situation, it seems prudent to re-examine the assumptions underlying the fundamental problem: that all molecular clouds are gravitationally bound and that they are collapsing at freefall, with the freefall time calculated from the mean cloud density of 100 cm⁻³. The question of whether the cloud is bound is usually addressed by calculating the virial parameter, α_vir (Bertoldi & McKee 1992). The definition of this parameter varies. The most common definition assumes virial equilibrium. Ignoring surface terms and magnetic fields,

$\begin{eqnarray}&&{\alpha }_{\mathrm{vir}}=\displaystyle \frac{2{E}_{\mathrm{kin}}}{| {E}_{g}| },\end{eqnarray} \tag{ 4 }$

where E_g is the gravitational potential energy and E_kin is the kinetic energy. With this definition, the dividing line between bound and unbound is α_vir = 2, where E_kin = ∣E_g∣. Surface pressure from the surrounding medium would increase α_vir and magnetic fields would decrease α_vir (Tan et al. 2013). Because we cannot measure them, we neglect the largely unknown pressure and magnetic terms. Then, a cloud with E_kin > E_g would tend to expand or dissipate, so we take α_vir = 2 as our critical value to define a cloud as bound if α_vir < 2 and unbound if α_vir > 2.

The calculation of E_g depends on cloud shape and structure (Bertoldi & McKee 1992; Kauffmann et al. 2013) but these effects can be considered as correction factors for the expression for a uniform density sphere:

$\begin{eqnarray}&&{E}_{g}=\displaystyle \frac{-3{GM}^{2}}{5R}\end{eqnarray} \tag{ 5 }$

where M and R are the cloud mass and radius. The kinetic energy is estimated from

$\begin{eqnarray}&&{E}_{\mathrm{kin}}=\displaystyle \frac{3{\sigma }_{v}^{2}M}{2}\end{eqnarray} \tag{ 6 }$

where σ_v is the one-dimensional velocity dispersion averaged over the cloud.

Unbound clouds may of course contain bound structures. Measuring α_vir in the complex structures of the interstellar medium is more problematic, as surface terms are likely to be important (Mao et al. 2020). Kim et al. (2021) compared the virial parameters from their model clouds to those that would be derived from observations. They found that the methods used by observers are more likely to overestimate the boundedness of the cloud, but generally by factors of 2 or less. They caution that neglect of tidal effects in crowded regions would make the overestimation worse. This effect was indeed apparent in larger-scale simulations (Mao et al. 2020), which found that many regions that were apparently bound, based on the criteria that α_vir < 2 were in fact unbound. In particular, the boundedness of large massive structures was more likely to be overestimated using the simpler estimations available to observers. If they are correct, our neglect of surface terms for practical reasons is unlikely to result in systematic underestimations of α_vir for clouds. Of course, there are circumstances where surface pressure may play an important role, especially on smaller scales (e.g., Heigl et al. 2018).

Given our simple approach to the virial parameter, ignoring external pressure, magnetic fields, and tides, we use the term "bound" as shorthand for gravitationally bound. A cloud could be unbound in our definition, but contained or "bound" by external pressure. Our focus is on whether gravity or turbulence is dominant. If the results of the simulations discussed above apply, our criteria for boundedness will be more generous than reality. We calculate virial parameters from

$\begin{eqnarray}&&{\alpha }_{\mathrm{vir}}=\displaystyle \frac{2{E}_{\mathrm{kin}}}{| {E}_{g}| }=1160{\sigma }_{v}^{2}{R}_{{pc}}/{M}_{{M}_{\odot }}=209{\rm{\Delta }}{v}^{2}{R}_{{pc}}/{M}_{{M}_{\odot }},\end{eqnarray} \tag{ 7 }$

where ${M}_{{M}_{\odot }}$ is the actual mass of the cloud in solar masses.

Estimating the actual mass of the cloud can be done in a number of ways, but the original source of the fundamental problem relied on observations of the lowest rotational transition of CO. The problem is that this line is usually very optically thick. Many words have been expended explaining why this can work, including an implicit assumption that the clouds are in virial equilibrium (e.g., Bolatto et al. 2013). If this were the only argument, the whole discussion would be circular, but we will take the CO-derived masses as independent of the assumption of virial equilibrium. The correlation of CO emission with dust extinction via the X(CO) factor provides an independent basis for its use in determining column density (Dickman 1978; Heiderman et al. 2010; Pineda et al. 2010). While column density may be underestimated at high column density and overestimated at low column density, both by factors of 2 or more, Bolatto et al. (2013) argued that variations in X(CO) average out, so L(CO) may reflect mass better than the intensity traces column density for a particular line of sight. The relation is M_CO = α_CO L_CO. We follow Bolatto et al. (2013) in taking an α_CO of 4.3 M_⊙ (K km s⁻¹ pc²)⁻¹. Lada & Dame (2020) have recently suggested that α_CO varies strongly with Galactocentric radius and is higher than the canonical value at almost all radii. We do not adopt this approach, but we will discuss the issue of possible variations in α_CO in Section 4.6.

Samples using ¹³CO J = 1 → 0 to get mass usually combine the ¹³CO data with CO data to derive the optical depth and then column density of ¹³CO, N₁₃, assuming that the excitation temperatures of CO and ¹³CO are equal and that the partition function can be approximated by assuming that the excitation temperature applies to all levels (commonly referred to as the "LTE" assumption). The excitation temperature can be found from the CO brightness temperature. Because we use CO data with low spatial resolution, we instead assume an excitation temperature of 8 K for all lines of sight (Roman-Duval et al. 2010). We follow the description by Ripple et al. (2013) but note the error in their Equation (4), in which the exponential factor should be in the denominator. Conversion to the column density of H₂ follows

$\begin{eqnarray}&&N({{\rm{H}}}_{2})={N}_{13}\left[\displaystyle \frac{{{\rm{H}}}_{2}}{\mathrm{CO}}\right]\left[\displaystyle \frac{\mathrm{CO}}{{\mathrm{CO}}^{13}}\right].\end{eqnarray} \tag{ 8 }$

We have used the latest number for H₂/CO of 6000 (Lacy et al.2017). Because the new measurement directly relates the abundance of H₂, not including helium, to that of CO through infrared absorption observations, we use the mean molecular weight of ${\mu }_{{{\rm{H}}}_{2}}{m}_{{\rm{H}}}$ , with ${\mu }_{{{\rm{H}}}_{2}}=2.8$ to convert to mass (Kauffmann et al. 2008). We have corrected data to use these values. Unless otherwise noted, when analyzing the ¹³CO J = 1 → 0 data, we account for a Galactic gradient in CO/¹³CO. We apply the elemental isotope gradient in the Galaxy most recently derived by Jacob et al. (2020),

$\begin{eqnarray}&&\displaystyle \frac{{{\rm{C}}}^{12}}{{{\rm{C}}}^{13}}=5.87{R}_{\mathrm{gal}}+13.25,\end{eqnarray} \tag{ 9 }$

where R_gal is the distance of the cloud from the Galactic Center in kiloparsecs.

2. Samples

The sets of data used in this paper are listed in Table 1. We consider structures defined by CO emission, ¹³CO emission, and dust emission, supplemented by molecular line data. From the catalogs in the papers, we extracted the information on size, mass, mass surface density (Σ_M), and virial parameter (α_vir), or the parameters needed to calculate them. In most cases, some pruning of the catalog entries was necessary. The cuts used are explained below.

Table 1. Sample List

Name	Tracer	Reference	Notes
CO-Rice	CO J = 1 → 0	Rice et al. (2016)	1
CO-MD	CO J = 1 → 0	Miville-Deschênes et al. (2017)	2
CO-CMZ	CO J = 1 → 0	Oka et al. (2001)	3
CO-OGS	CO J = 1 → 0	Heyer et al. (2001)	4
CO-Sun	CO J = 1 → 0	Sun et al. (2020)	5
¹³CO-GRS	¹³CO J = 1 → 0	Roman-Duval et al. 2010	6
¹³CO-EXFC55-100	¹³CO J = 1 → 0	Roman-Duval et al. 2016	7
¹³CO-EXFC135-195	¹³CO J = 1 → 0	Roman-Duval et al. 2016	8
¹³CO-SEDIGISM	¹³CO J = 2 → 1	Duarte-Cabral et al. (2021)	9
Herschel	submillimeter + NH₃	Merello et al. (2019)	10

Notes.

¹Recovers 25% of emission in the survey by Dame et al. (2001).²Recovers 98% of emission in the survey by Dame et al. (2001).³Attempted to restrict to CMZ.⁴Outer Galaxy Survey: $\left\langle {R}_{\mathrm{gal}}\right\rangle =11.8\pm 1.3$ $\left\langle {R}_{\mathrm{gal}}\right\rangle =11.8\pm 1.3$ .⁵Survey of 66 galaxies.⁶Galactic Ring Survey: $\left\langle {R}_{\mathrm{gal}}\right\rangle =5.2\pm 1.0$ $\left\langle {R}_{\mathrm{gal}}\right\rangle =5.2\pm 1.0$ .⁷Exeter-FCRAO Survey for 55 < l < 100: $\left\langle {R}_{\mathrm{gal}}\right\rangle =10.3\pm 1.9$ $\left\langle {R}_{\mathrm{gal}}\right\rangle =10.3\pm 1.9$ .⁸Exeter-FCRAO Survey for 135 < l < 195: $\left\langle {R}_{\mathrm{gal}}\right\rangle =12.9\pm 2.1$ $\left\langle {R}_{\mathrm{gal}}\right\rangle =12.9\pm 2.1$ .⁹Covers R_gal from 0.6 to 15.7 kpc: $\left\langle {R}_{\mathrm{gal}}\right\rangle =5.5\pm 1.8$ $\left\langle {R}_{\mathrm{gal}}\right\rangle =5.5\pm 1.8$ .¹⁰Based on the Herschel Infrared Galactic Plane Survey.

Download table as: ASCII Typeset image

2.1. Structures Defined by CO J = 1 → 0 Emission

The samples in this category rely on masses measured from the luminosity of the J = 1 → 0 transition of CO. The conversion from luminosity to mass (α_CO) differed somewhat among the studies; to standardize them, the value of M_CO in the Milky Way samples was calculated using α_CO = 4.3 M_⊙ (K km s⁻¹ pc²)⁻¹ (Bolatto et al. 2013).

New catalogs of clouds based on decomposition of the CO surveys are now available (Rice et al. 2016; Miville-Deschênes et al. 2017). These are particularly interesting as they used the data from Dame et al. (2001) that was featured by Heyer & Dame (2015) but analyzed it in two different ways. These two studies, labeled CO-Rice and CO-MD in Table 1, take different approaches to breaking the large-scale CO emission in the Galaxy into clouds, and we will discuss the effect of those differences in Section 4.1.

Rice et al. (2016) used a dendrogram analysis to identify 1064 clouds. They recovered 2.5 × 10⁸ M_⊙, about 25% of the total molecular mass of the Galaxy. Table 3 of Rice et al. (2016) provides a radius, σ_v, and M_CO, the mass determined from the CO luminosity, which we use directly. The mass surface density was determined from M_CO and the radius. No mass or size information was given for clouds without clear resolution of the kinematic distance ambiguity, and one cloud was assigned zero size and mass in Table 3 of Rice et al. (2016); after removing these, 1037 clouds remained for our analysis.

Miville-Deschênes et al. (2017) used a Gaussian decomposition and hierarchical cluster analysis to assign essentially all of the CO emission in the Dame survey to 8107 clouds comprising 1.2 × 10⁹ M_⊙. Their online table provides radii and masses for both near and far kinematic distances, along with a code for resolution of the kinematic distance ambiguity, which we use to select the best radius and mass. Miville-Deschênes et al. (2017) calculated a radius from the two projected radii by $r={({r}_{\max }{r}_{\min }{r}_{\min })}^{1/3}$ , arguing that the depth was more likely equal to the smaller size projected on the sky. They also provide σ_v, which we use in calculating α_vir. The masses were based on summing column densities, using X_CO = 2 × 10²⁰ cm⁻² (K km s⁻¹)⁻¹, which corresponds to α_CO = 4.3 M_⊙ (K km s⁻¹ pc²)⁻¹ (Bolatto et al. 2013). The catalog contains a small number of clouds with nominal Galactocentric radii greater than 30 kpc and a large number of low-mass objects with very high α_vir. To select against these, we require that R_gal < 30 kpc, M_CO > 1 M_⊙, and ${\alpha }_{\mathrm{vir}}\lt {\alpha }_{\max }$ , where ${\alpha }_{\max }$ is an arbitrarily chosen maximum value. Entry 2 in Table 2 required α_vir < 100, while entry 3 used α_vir < 20.

Table 2. Statistics

Sample	Number	α_vir			Σ_M (M_⊙ pc⁻²)			f(N)	f(M)	f_c(N)	f_c(M)	Note
		Med	Mean	Std.	Med	Mean	Std.
CO-Rice	1037	2.44	3.45	3.34	21.2	38.3	43.4	0.38	0.73	0.33	0.41	1
CO-MD	7516	7.58	15.20	18.53	21.9	41.1	53.9	0.07	0.19	0.06	0.09	2
CO-MD	5776	5.49	6.91	4.70	31.4	49.4	58.1	0.09	0.19	0.08	0.10	3
CO-CMZ	110	3.26	4.45	3.46	1311.3	1703.2	1271.0	0.29	0.71	0.00	0.00	4
CO-OGS	380	2.27	2.51	1.25	29.3	35.4	20.2	0.41	0.72	0.40	0.56	5
CO-OGS	3714	6.15	7.22	4.79	17.7	21.0	11.7	0.06	0.58	0.06	0.46	6
CO-Sun	102,788	3.47	4.57	4.10	21.6	49.9	144.3	0.21	0.35	0.18	0.15	7
¹³CO(1−0)	737	1.03	1.95	2.76	82.0	85.5	39.8	0.72	0.95	0.58	0.55	8
¹³CO-CO-GRS	289	2.86	3.83	3.02	80.7	103.7	84.0	0.31	0.46	0.17	0.19	9
¹³CO-CO-EXFC55-100	105	0.99	1.92	2.24	25.9	36.1	32.8	0.62	0.73	0.58	0.70	10
¹³CO-CO-EXFC135-195	98	1.37	2.30	3.00	21.1	24.2	13.0	0.63	0.75	0.63	0.75	11
¹³CO(2−1)	6658	1.25	1.94	2.66	74.7	87.7	49.0	0.73	0.79	0.62	0.41	12
Herschel	1067	0.29	0.58	1.08	2497.7	3249.8	3032.2	0.96	1.00	0.00	0.00	13

Notes.

¹Data from Rice et al. (2016). Clouds with ambiguous distance assignments and one cloud with zero mass were eliminated.²Data from Miville-Deschênes et al. (2017) using criteria M_CO > 1 M_⊙, R_gal < 30 kpc, and α_vir < 100.³Data from Miville-Deschênes et al. (2017) using criteria M_CO > 1 M_⊙, R_gal < 30 kpc, and α_vir < 20.⁴Data from Oka et al. (2001), considering only clouds in CMZ.⁵Data from Heyer et al. (2001), using criterion M_CO > 2.5 × 10³ M_⊙.⁶Data from Heyer et al. (2001), using criterion M_CO > 1.0 × 10¹ M_⊙.⁷Data from Sun et al. (2020).⁸GRS Data from Roman-Duval et al. (2010), updated abundance, some selection.⁹Re-analysis of GRS ¹³CO data within CO clouds.¹⁰Re-analysis of EXFC 55-100 ¹³CO data within CO clouds.¹¹Re-analysis of EXFC 135-195 ¹³CO data within CO clouds.¹²Data from Duarte-Cabral et al. (2021).¹³Data from Merello et al. (2019), using criterion M_dust > 1 M_⊙.

Download table as: ASCII Typeset image

The sample labeled CO-CMZ in Table 1 is from Oka et al. (2001) and focuses on clouds in the Central Molecular Zone (CMZ) using data obtained with the Nobeyama Radio Observatory 45 m telescope. Clouds were defined by topologically closed surfaces in (l, b, v) space. Three different thresholds were used and some clouds appear multiple times in the catalog, providing a range of surface densities. The dispersions in all three coordinates were used to define sizes and line widths. A distance of 8.5 kpc was used for all clouds. They accounted for truncation by the threshold by comparing to simulated clouds and provided corrected values for sizes, velocity dispersions, and CO luminosities. These are what we use. To calculate the mass surface density, we take the size as a radius. They also computed virial masses, but we do not use these. We eliminate clouds determined to be outside the CMZ (labeled D in their Table 1) or confused (labeled with C in their table 1). The use of α_CO = 4.3 M_⊙ (K km s⁻¹ pc²)⁻¹ in the CMZ is dubious. Oka et al. (1998) suggested a value 10 times smaller. Elemental abundances certainly are higher in the CMZ and a correction for that might be appropriate, as suggested by Sun et al. (2020) for other galaxies.

The sample labeled CO-OGS in Table 1 is a sample of clouds in the outer Galaxy mapped by the FCRAO 14 m telescope (Heyer et al. 1998) and tabulated by Heyer et al. (2001). The Table provides L_CO, Δv (FWHM), and diameter of the structure. Heyer et al. (2001) used α_CO = 4.1 to convert to mass, but we changed to 4.3 for consistency. We convert diameter to radius for consistency. The linewidth used to calculate M_vir was that of a spectrum averaged over the cloud. These were almost always larger than an alternative estimate using the dispersion of line centroids over the cloud, indicating that the clouds are generally not dominated by systematic motions such as shear. The catalog contains 10,156 structures, but many are local. To select outer Galaxy sources, Heyer et al. (2001) required v_LSR < −25 km s⁻¹. To select against small structures, they also required M_CO > 2.5 × 10³ M_⊙. We give values for this cut, but we relax it to M_CO > 10 M_⊙ for the figures and most analysis because we are more interested in completeness in structures than in defining a mass function.

2.2. Structures Defined by CO J = 2 → 1 Emission in Other Galaxies

Studies of other galaxies are beginning to provide similar information without the line-of-sight confusion of Milky Way studies. For other galaxies, we use the remarkable compilation by Sun et al. (2020) of the PHANGS-ALMA study of the CO J = 2 → 1 transition, labeled CO-Sun in Table 1, which is part of a larger program including Hubble Space Telescope observations (Lee et al. 2021). They collected data on over 10⁵ sightlines toward 66 galaxies, resampled to a uniform spatial resolution of 150 pc, by far the largest database available, albeit with very poor spatial resolution by the standards of studies in our Galaxy. They used a metallicity dependent value for α_CO, and we use their values. Their Table 2 provides Σ_M and α_vir. We compute M_CO from Σ_M and the fixed size of 150 pc and use their α_vir.

2.3. Structures Defined by ¹³CO J = 1 → 0 Emission

Roman-Duval et al. (2010) provided the sample defined by ¹³CO J = 1 → 0 emission from the Galactic Ring Survey (GRS) by Jackson et al. (2006). The ¹³CO data were combined with the CO data from the University of Massachusetts−Stony Brook (UMSB) survey (Sanders et al. 1985; Clemens et al. 1986) to determine optical depth and ¹³CO column density. The clouds in the Roman-Duval et al. (2010) study were actually identified by Rathborne et al. (2009), who presented a detailed description of the method, using CLUMPFIND (Williams et al. 1994) and a number of tests to be able to identify both bright, compact structures and fainter, more diffuse structures. The effective lowest brightness level was 0.2 K, so they should not have selected regions of abnormally high density.

For the original analysis by Roman-Duval et al. (2010) of the GRS survey, we consider these to be structures defined by ¹³CO. Radius, mass, surface density, and α_vir are provided by their Table 1, with 749 entries. Requiring M_CO > 1 M_⊙, Σ_M > 10 M_⊙ pc⁻², and α_vir < 20 pruned the sample to 737. The condition on α_vir removed a few extreme outliers. A total of 170 clouds in the sample were not covered by the UMSB CO survey, so the usual method to determine the column density of ¹³CO was impractical. If we exclude those, the results do not differ substantially. Roman-Duval et al. (2010) used a different definition of α_vir, assuming a correction for the density distribution. For consistency, we recalculated α_vir from Equation (7). These values were larger by a factor of 1.3 on average than those given in their paper when the same assumptions about abundances were made.

For the ¹³CO J = 1 → 0 samples from other parts of the Galaxy, both ¹³CO and CO J = 1 → 0 were observed simultaneously. These are the Exeter-FCRAO surveys, which are divided into EXFC55-100 and EXFC135-195, but our analysis excludes the range of longitude where kinematic distances are unreliable (165 < l < 195), as described by Roman-Duval et al. (2016). For the EXFC surveys, Roman-Duval et al. (2016) attempted to identify all voxels with emission to distinguish diffuse (CO only) from "dense" gas (¹³CO and CO). Structures were not identified and no catalogs were provided. We have re-analyzed all the ¹³CO J = 1 → 0 surveys in a uniform way, using the CO-defined structures from Miville-Deschênes et al. (2017), with the procedure detailed in the Appendix. Because the masses are based on all of the ¹³CO data within the CO-defined structures, we refer to these results as ¹³CO-CO, followed by the survey name, such as GRS.

The GRS sample is concentrated in the inner Galaxy; the subsample that we analyze has median (R_gal) = 5.2 kpc, $\left\langle {R}_{\mathrm{gal}}\right\rangle =5.2\pm 1.0\,\mathrm{kpc}$ . The EXFC55-100 sample is concentrated closer to the solar circle, with 7.2 < R_gal < 14.5, median (R_gal) = 10.7 kpc, $\left\langle {R}_{\mathrm{gal}}\right\rangle =10.3\pm 1.9\,\mathrm{kpc}$ . The EXFC135-195 data are concentrated in the outer Galaxy, with 10.3 < R_gal < 20.4, median (R_gal) = 12.5 kpc, $\left\langle {R}_{\mathrm{gal}}\right\rangle =12.9\pm 2.1\,\mathrm{kpc}$ . They thus provide convenient samples for examining the effects of R_gal.

2.4. Structures Defined by ¹³CO J = 2 → 1 Emission

The SEDIGISM survey (Duarte-Cabral et al. 2021) has provided a catalog of structures defined by emission from the J = 2 → 1 transition of ¹³CO. They used the SCIMES algorithm described by Colombo et al. (2019) to define the structures. They determined column densities by assuming a constant ratio of H₂ column density to integrated intensity of ¹³CO J = 2 → 1 of $X{{(}^{13}\mathrm{CO}\,{\text{}}J=2\to 1)=1\times {10}^{21}\,{\mathrm{cm}}^{-2}({\rm{K}}\,\mathrm{km}\,{{\rm{s}}}^{-1})}^{-1}$ and a mean molecular weight per H₂ of 2.8. Because their choice of X(¹³COJ = 2 → 1) was empirically based (Schuller et al. 2017), we did no re-scaling. The sample was pruned to remove those with ambiguous distances, near edges of maps, or smaller than 114 pixels, roughly re-creating their "science sample." The masses and virial parameters in their tables were used directly, but α_vir was checked and agreed with a calculation from their mass, size, and linewidth information.

2.5. Structures Defined by Herschel and NH₃

Structures identified from the Herschel and ATLASGAL maps of the Milky Way were presented by Merello et al. (2019). They meet a number of requirements, but most relevant to us, they have NH₃ observations, which are used to obtain a linewidth from the inversion transition in the (J, K) = (1, 1) level and a temperature from the ratio of (1, 1) and (2, 2) lines. We use the sample of 1068, called the "Final Sample" by Merello et al. (2019). They provide mass, Δv (FWHM), and radius, from which we compute α_vir and Σ_M. The masses were computed from a fit to the spectral energy distribution, the resulting optical depth at 300 μm, and an opacity of 0.1 cm² g⁻¹ of gas plus dust. The radius was half the FWHM of the emission at 250 μm.

3. Results

For each of the samples, four plots are provided: (1) the mean, median, and standard deviation of α_vir versus mass; (2) the mean, median, and standard deviation of α_vir versus surface density (Σ_M); (3) a histogram of ${\rm{l}}{\rm{o}}{\rm{g}}{{\rm{\Sigma }}}_{M}$ ; (4) a histogram of $\mathrm{log}{\alpha }_{\mathrm{vir}}$ . Statistics of the relevant properties are summarized in Table 2. For each sample, Table 2 has the median, mean, and standard deviation of α_vir and the median, mean, and standard deviation of Σ_M. The table also includes the fraction of the total number of structures that satisfy α_vir < 2. The fraction of the number (f(N)), and the fraction of the total mass (f(M)) satisfying these conditions are both tabulated.

We will compare the criterion for the virial parameter to the surface density criterion for star formation. In local clouds, star formation is strongly concentrated in regions with Σ_M ≥ 120 M_⊙ pc⁻² (Heiderman et al. 2010; Lada et al. 2010, 2012; Evans et al. 2014). For simplicity, we refer to structures with Σ_M <120 M_⊙ pc⁻² as clouds and those with Σ_M ≥ 120 M_⊙ pc⁻² as clumps and the plots of α_vir versus Σ_M include a line demarcating that boundary. We emphasize that this definition of "clump" does not include constraints often used to define clumps, such as size information, and we use the terminology only for simplicity. This usage of "clump" is also distinguished from "dense clumps" as identified in lines of much higher effective density, as in the Herschel sample or samples studied by Wu et al. (2010) and others. We use bins of 0.3 in log Σ_M, on either side of 2.08, corresponding to Σ_M = 120 M_⊙ pc⁻². Similarly, we summarize the boundedness with labels of "bound" or "unbound" on either side of a line at α_vir = 2 in the plots. Finally, the fraction of the number (f_c(N)) and mass (f_c(M)) that satisfy both α_vir < 2 and Σ_M < 120 M_⊙ pc⁻² are tabulated, the subscript (c) indicating that they apply only to structures we have called clouds.

3.1. Structures Defined by CO J = 1 → 0 Emission

We plot α_vir calculated from the data in Rice et al. (2016) in Figure 1. The mean α_vir declines with mass, so that $\left\langle {\alpha }_{\mathrm{vir}}\right\rangle \lt 2$ for log M_CO > 5.7. The mean α_vir also declines with Σ_M, so that $\left\langle {\alpha }_{\mathrm{vir}}\right\rangle \lt 2$ for log Σ_M > 1.7. The histograms indicate that most structures are clouds and most are unbound (f(N) = 0.38), but the fraction of mass in bound structures is 0.73, reflecting the tendency for more massive structures to be bound.

The catalog of Miville-Deschênes et al. (2017), as noted in Section 2, contains some outliers with very high nominal α_vir, so we exclude those greater than a maximum value, ${\alpha }_{\max }$ . Entries in Table 2 are given for two choices, ${\alpha }_{\max }=100$ and ${\alpha }_{\max }=20$ . The mean and standard deviation are plotted versus mass in Figure 2. As we saw with the catalog of Rice et al. (2016), the $\left\langle {\alpha }_{\mathrm{vir}}\right\rangle$ values decline with M_CO. In this case, they reach $\left\langle {\alpha }_{\mathrm{vir}}\right\rangle =2$ only for the most massive bin. The mean and median values of α_vir decline with Σ_M, but never drop below α_vir = 2. The histograms indicate that the vast majority of structures are unbound clouds. The fraction of the mass that is bound is relatively insensitive to the choice of ${\alpha }_{\max }$ , at about 0.19 for all structures. The fact that these results differ substantially from those of Rice et al. (2016) will be discussed in Section 4.1.

**Figure 2.** Upper left: the mean and standard deviation of the virial ratio is plotted vs. the logarithm of the clump mass, for the catalog of Miville-Deschênes et al. (2017) after selecting only those catalog entries with α_vir < 100 and M_CO > 1 M_⊙. Median values are plotted in a magenta dashed line. Upper right: the same quantities are plotted vs. the mass surface density. The horizontal line at α_vir = 2 demarcates nominally unbound clouds above the line from nominally bound clouds below the line. The vertical line in the right panel indicates Σ_M = 120 M_⊙ pc⁻². Lower left: the histogram of values of log α_vir, with a vertical red line at α_vir = 2.0. The number of structures in each category are given in parentheses. Lower right: the histogram of values of log Σ_M, with a vertical red line at Σ_M = 120 M_⊙ pc⁻².
Download figure:
Standard image High-resolution image

Figure 3 shows that the structures in the CMZ in the catalog of Oka et al. (2001) are generally unbound even though they almost all have Σ_M > 120 M_⊙ pc⁻², satisfying the local criterion for star formation. Clearly, the definition of a clump based on the local surface density criterion for active star formation does not work in the CMZ. While most structures are unbound, the tendency for higher mass structures to be bound results in a relatively high f(M) = 0.71, with the caveat that these masses may well be overestimated (Section 2).

The outer Galaxy sample (OGS) from the catalog of Heyer et al. (2001) shows a similar behavior (Figure 4). Most are unbound, but more than half the mass is in bound structures despite having Σ_M < 120 M_⊙ pc⁻², with f(M) = 0.72 (Table 2). To some extent, this result depends on the requirement that M_CO > 2.5 × 10³ M_⊙ imposed by Heyer et al. (2001) to ensure completeness. Since completenesss is less important for our purposes, we relaxed that criterion to M_CO > 10 M_⊙ to see the effect. While f(N) dropped to 0.06, f(M) = 0.58, still substantial because most of the bound mass was in the higher mass structures.

3.2. Structures Defined by CO J = 2 → 1 Emission in Other Galaxies

We were able to use all of the catalog entries with 150 pc resolution from Table B1 of Sun et al. (2020), save one entry with a negative value for α_vir. Sun et al. (2020) themselves found an area-weighted median α_vir = 3.5 with a 1σ spread of 0.6 dex, based on over 10⁵ pixels of size 150 pc in 66 galaxies. The mass-weighted median was smaller at 2.7. All of the caveats about neglect of tidal effects, etc. apply quite strongly to studies with 150 pc resolution. Independent analysis of the data in their Table B1 confirms the median value and finds a mean value of $\left\langle {\alpha }_{\mathrm{vir}}\right\rangle =4.7\pm 4.0$ . The averages of α_vir in Figure 5 all lie above α_vir = 2, but they approach that value for the most massive data points. A similar pattern applies to the plot of α_vir versus Σ_M. While many points have very high α_vir, the fraction with α_vir < 2 is f(N) = 0.21 and the fraction of the total mass in those structures is f(M) = 0.35, reflecting the usual trend for lower α_vir in more massive clouds. The median and mean surface densities are 22 and 50 M_⊙ pc⁻², but with a maximum at Σ_M = 9.7 × 10³ M_⊙ pc⁻². The last value is clearly more characteristic of dense clumps, not clouds. All points above log M_CO = 6.43 M_⊙ correspond to Σ_M > 120 M_⊙ pc⁻² or A_V > 8 mag. Vertical lines on Figure 5 indicate surface densities of 120 M_⊙ pc⁻² and 1 g cm⁻². The second of these corresponds to the criterion for massive star formation (Krumholz & McKee 2008). The fact that α_vir drops to values near two as the surface density approaches this value is interesting. To have such high surface density averaged over a region of 150 pc is quite remarkable. If Σ_M is restricted to be less than 120 M_⊙ pc⁻², the fraction of the total mass with α_vir < 2 becomes 0.15.

**Figure 5.** The mean and standard deviation of the virial ratio is plotted in bins vs. the logarithm of the clump mass, for the catalog of Sun et al. (2020). Upper left: the mean and standard deviation of the virial ratio is plotted vs. the logarithm of the clump mass. Median values are plotted in a magenta dashed line. Upper right: the same quantities are plotted vs. the mass surface density. The horizontal line at α_vir = 2 demarcates nominally unbound clouds above the line from nominally bound clouds below the line. Vertical lines indicate Σ_M = 120 M_⊙ pc⁻² and Σ_M = 1 g cm⁻². Lower left: the histogram of values of log α_vir, with a vertical red line at α_vir = 2.0. The number of structures in each category are given in parentheses. Lower right: the histogram of values of log Σ_M, with a vertical red line at Σ_M =120 M_⊙ pc⁻².
Download figure:
Standard image High-resolution image

3.3. Structures Defined by ¹³CO J = 1 → 0 Emission

To use data from all three surveys and to compare methods for structure identification, we have re-analyzed the three data sets, using the CO clouds defined by Miville-Deschênes et al. (2017) to constrain the ¹³CO voxels included. The detailed method is explained in the Appendix. However, we begin with the original catalog, but with updated abundance assumptions. These results will allow a comparison between different cloud identification methods.

3.3.1. GRS: Original Sample

The tabulated masses and surface densities for the ¹³CO GRS data from Roman-Duval et al. (2010) originally used values of H₂/CO of 12,500 and ${\mu }_{{{\rm{H}}}_{2}}=2.37$ . As discussed in Section 1, we instead use the current best values for H₂/CO and ${\mu }_{{{\rm{H}}}_{2}}$ , resulting in masses and surface densities lower by a factor of 0.567 and virial parameters higher by a factor of 1.76. For a first analysis, (entry 8 in Table 2), we retain the assumption from Roman-Duval et al. (2010) of a constant CO/¹³CO of 45. The usual quantities are plotted in Figure 6 and the results are in Table 2. The majority of the sample would be classified as clouds, with Σ_M < 120 M_⊙ pc⁻², but about 1/6 are denser. The fraction of mass in bound structures, f(M) = 0.95, but the fraction for clouds is much less, f_c(M) = 0.55.

**Figure 6.** Upper left: the mean and standard deviation of the virial ratio is plotted vs. the logarithm of the clump mass, for the catalog of Roman-Duval et al. (2010), after selecting only those catalog entries with Σ_M > 10 and M_CO > 1 M_⊙, and α_vir < 20. Median values are plotted in a magenta dashed line. Upper right: the same quantities are plotted vs. the mass surface density. The horizontal line at α_vir = 2 demarcates nominally unbound clouds above the line from nominally bound clouds below the line. The vertical line in the right panel indicates Σ_M = 120 M_⊙ pc⁻². Lower left: the histogram of values of log α_vir, with a vertical red line at α_vir = 2.0. The number of structures in each category are given in parentheses. Lower right: the histogram of values of log Σ_M, with a vertical red line at Σ_M = 120 M_⊙ pc⁻².
Download figure:
Standard image High-resolution image

3.3.2. GRS: New Analysis

From now on, we use the CO-defined structures and the ¹³CO-derived properties, with the abundance assumptions described at the end of Section 1, including the variation of CO/¹³CO with R_gal. We begin with a re-analysis of the GRS survey in ¹³CO (Roman-Duval et al. 2010). The structures are defined in the catalog of Miville-Deschênes et al. (2017), but only for a subsample. The subsample of CO clouds is somewhat biased toward higher Σ_M and lower α_vir than the full catalog: $\left\langle {{\rm{\Sigma }}}_{M}\right\rangle =95\pm 64$ M_⊙ pc⁻², nearly twice that for the full catalog; $\left\langle {\alpha }_{\mathrm{vir}}\right\rangle =5.62\pm 5.28$ , 0.81 times that for the full catalog.

The analysis of the GRS data with the method described in the Appendix is listed as entry 9 of Table 2 and the usual plots are in Figure 7. The main difference from entry 8 is that the clouds were identified in the CO analysis and the ¹³CO data were taken only within the boundaries of those clouds. The change in cloud definition makes the mean α_vir larger and f(M) smaller by substantial amounts (f(M) = 0.46).

The fraction of the mass of the CO-defined structure contained in the ¹³CO emission (M₁₃/M_CO) is a valuable guide to the fraction of possibly bound structures within unbound CO clouds. For the GRS sample, $\left\langle {M}_{13}/{M}_{\mathrm{CO}}\right\rangle =0.22\pm 0.08$ with a median M₁₃/M_CO = 0.21.

3.3.3. EXFC55-100

For this sample, most structures are identified as clouds, but they are mostly bound, despite having low Σ_M. The usual plots are in Figure 8 and the statistics are given as entry 10 in Table 2. The fraction of the CO mass recovered by the ¹³CO emission, $\left\langle {M}_{13}/{M}_{\mathrm{CO}}\right\rangle =0.79\pm 0.80$ with a median M₁₃/M_CO = 0.61. These values are much higher than those in the GRS.

3.3.4. EXFC135-195

As for EXFC55-100, most structures are identified as clouds, but they are mostly bound, despite having very low Σ_M. The usual plots are in Figure 9 and the statistics are given as entry 11 in Table 2. For the EXFC135-195 sample, $\left\langle {M}_{13}/{M}_{\mathrm{CO}}\right\rangle \,=0.76\pm 0.34$ with a median M₁₃/M_CO = 0.66. These values are much higher than those in the GRS and similar to those of the EXFC55-100 sample.

3.3.5. Summary of Consequences of New Analysis

The main effect of the different method of structure identification is to decrease the fraction of mass in bound structures. In the inner Galaxy, most become unbound, but 46% of the mass is still in bound structures. For larger R_gal, the structures become less dense but more bound. For the outer Galaxy (EXFC135-195) sample, f(M) = 0.75 even though $\left\langle {{\rm{\Sigma }}}_{M}\right\rangle =24$ M_⊙ pc⁻². This trend reflects a strong decrease in velocity dispersion with R_gal. The fraction of the mass recovered by ¹³CO (M₁₃/M_CO) is low in the inner Galaxy, but quite high in the other two samples.

3.4. ¹³CO J = 2 → 1

The situation for the J = 2 → 1 transition of ¹³CO (Duarte-Cabral et al. 2021), shown in Figure 10, is different. Most clouds defined in this sample are nominally bound, even for relatively low mass and Σ_M. However, $\left\langle {\alpha }_{\mathrm{vir}}\right\rangle$ is lower and f(M) is higher than for the J = 1 → 0 data in the GRS sample with the new analysis, but comparable for the other J = 1 → 0 samples.

3.5. Stuctures Defined by Herschel and NH₃

Compact structures identified by their submillimeter emission with Herschel have very low virial parameters, strongly indicative of being gravitationally dominated and f(M) = 1.0. Figure 11 plots the virial parameter versus the logarithm of clump mass from an analysis by Merello et al. (2019). Only the clumps with the lowest mass and Σ_M are unbound. While most of these contain star formation, even the 157 "prestellar" clumps are strongly gravitationally dominated, with $\left\langle {\alpha }_{\mathrm{vir}}\right\rangle =0.57$ , almost identical to the value for the full sample. The values for f_c(N) and f_c(M) of zero reflect the fact that the very few structures in the sample with surface densities small enough to be considered clouds are all unbound.

**Figure 11.** Upper left: the mean and standard deviation of the virial ratio is plotted vs. the logarithm of the clump mass, for the catalog of Merello et al. (2019). Median values are plotted in a magenta dashed line. Upper right: the same quantities are plotted vs. the mass surface density. The horizontal line at α_vir = 2 demarcates nominally unbound clouds above the line from nominally bound clouds below the line. The vertical lines in the right panel indicate Σ_M = 120 M_⊙ pc⁻² and Σ_M = 1 g cm⁻³. The mass and surface density are the total (gas+dust) determined from submillimeter continuum emission. Lower left: the histogram of values of log α_vir, with a vertical red line at α_vir = 2.0. The number of structures in each category are given in parentheses. Lower right: the histogram of values of log Σ_M, with a vertical red line at Σ_M = 120 M_⊙ pc⁻².
Download figure:
Standard image High-resolution image

4. Discussion

In this section, we will explore how the boundedness of structures depends on methods and properties. We include subsections on the method used to identify structures, the tracer used to define the structures and measure the mass, the location in the Galaxy, the mass of the structure, and the surface density of the structure.

4.1. Dependence on Method Used to Identify Structure

The first issue to consider is the fact that the results from two different decompositions of the CO J = 1 → 0 survey by Dame et al. (2001) show quite different results. The difference is clear from the figures, but it can be encapsulated in the values of f(M): the structures identified by Rice et al. (2016) have f(M) = 0.73 and those of Miville-Deschênes et al. (2017) have f(M) = 0.19.

Clearly, the choice of method to identify structures plays a major role in the different results. We can borrow the terminology of taxonomists to distinguish lumpers and splitters. In our context, lumpers would aggregate emission into larger structures, while splitters would break emission into smaller structures. Rice et al. (2016) used a dendrogram analysis to identify 1064 clouds, accounting for about 25% of the molecular mass in the Galaxy. They note that their method avoids merging unrelated clouds into pseudo-clouds. In this sense, they are splitters. Comparing their catalog to that of Dame et al. (1986; based on an early version of the first quadrant survey), they recover all of the same clouds, but some are split and less massive. They use the size−linewidth relation and latitude distribution to break the kinematic distance ambiguity. In contrast, Miville-Deschênes et al. (2017) might reasonably be called lumpers. They performed a Gaussian decomposition of each spectrum and then lumped them together by a clustering procedure. This process led to larger linewidths in the clouds (median σ_v = 3.6 km s⁻¹) than in the individual Gaussian components (median σ_v = 1.65 km s⁻¹). Interestingly, this method allowed them to identify 98% of the CO emission with 8107 clouds. They found that the size−linewidth relation has too much scatter, so they use the relation ${\sigma }_{v}\propto {(R{{\rm{\Sigma }}}_{M})}^{0.5}$ , where R is the cloud size, to resolve the kinematic distance ambiguity. While the mass also increased by this process, the net effect was probably to produce larger α_vir. Were they lumpers or splitters? In the sense that they aggregated velocity components, they were lumpers, but they decomposed four times as much of the CO emission into eight times as many clouds as Rice, making them splitters. Lumping in velocity and splitting in number of clouds is likely to lead to larger α_vir. However, the fact that the scaling relation mentioned above holds indicates that unrelated gas is not being lumped in pseudo-clouds. Despite these differences, the scaling relations for the two catalogs are similar (Lada & Dame 2020).

Another difference between the results of the two methods is apparent in the x-axis values for the plots of α_vir versus M_CO. The plot for Rice et al. (2016) starts at $\mathrm{log}{M}_{\mathrm{CO}}=3.5$ , while that of Miville-Deschênes et al. (2017) starts at $\mathrm{log}{M}_{\mathrm{CO}}=0.5$ . The latter catalog contains many more clouds of lower mass. However, this does not explain the difference in f(M) because the total mass is still dominated by the massive clouds.

Probably the main factor in the difference is the completeness in accounting for the CO emission of the Galaxy. Because Rice et al. (2016) account for only 25% of the CO emission, their catalog is less suitable for addressing the question of slow star formation on the Galactic scale, where the problem was defined by taking the total mass of molecular gas from the CO survey and dividing by a freefall time.

In this context, the data from Sun et al. (2020) are interesting as no method is used to define structures. Instead only the mass within a 150 pc region is tabulated. The fact that f(M) = 0.35 and f_c(M) = 0.15, averaged over many galaxies and galactocentric radii, even with some very high values for Σ_M, suggests that the most likely value for f(M) is relatively small for structures measured by CO, regardless of which of the two transitions is used.

The importance of the identification method is also apparent in the difference between the results from the original structures in the Roman-Duval et al. (2010) GRS sample (entry 8 in Table 2) and the new analysis, where structures were defined by the CO catalog of Miville-Deschênes et al. (2017; entry 9 in Table 2).

These differences suggest that procedures like CLUMPFIND or dendrogram analyses tend to split structures or minimize low-level, extended emission. For addressing the global issue of slow star formation, procedures that account for all of the emission are necessary.

4.2. Dependence on the Tracer Used to Define Structure

We begin with the extremes: the samples of Miville-Deschênes et al. (2017) and Sun et al. (2020) showing low values for f(M) for mass measured from CO emission, versus the sample of Merello et al. (2019) showing f(M) = 1.0 for structures traced by Herschel peaks and NH₃ linewidths. The distinction between clouds, primarily traced by CO, and dense clumps, primarily traced by peaks in submillimeter emission and lines of species rarer than CO, has been based on this difference, but the analysis of α_vir strongly supports this distinction. The much higher mean and median Σ_M for the Herschel-NH₃ data also indicate that these are clearly distinct structures.

What about the structures probed by ¹³CO? They clearly lie between the extremes of clouds and dense clumps. For the GRS sample with our new analysis, f(M) = 0.46, $\left\langle {\alpha }_{\mathrm{vir}}\right\rangle =3.83$ , and $\left\langle {{\rm{\Sigma }}}_{M}\right\rangle =104$ . For the sample defined by ¹³CO J = 2 → 1, from Duarte-Cabral et al. (2021), f(M) = 0.79 and $\left\langle {\alpha }_{\mathrm{vir}}\right\rangle =1.94$ , even though $\left\langle {{\rm{\Sigma }}}_{M}\right\rangle =87.7$ M_⊙ pc⁻², lower than that of ¹³CO J = 1 → 0 in the GRS, but higher than that of ¹³CO J = 1 → 0 in the other two surveys.

The trend with tracer is that α_vir decreases (structures are more likely to be bound) as the tracer changes from CO to ¹³COJ = 1 → 0 to ¹³CO J = 2 → 1 to Herschel-NH₃. This can be related to a trend in the characteristic density of the material being traced. The characteristic density probed by each tracer can be crudely estimated from the effective density of the tracer used to define the structure. The effective density (n_eff), defined by Evans (1999) and developed further by Shirley (2015), is the density of particles needed to produce a 1 K km s⁻¹ observed line for an assumed column density and kinetic temperature. Because effective densities for CO and ¹³CO were not calculated by those references, we calculated them using the online tool RADEX (van der Tak et al. 2007) assuming T_K = 10 K and column density corresponding to A_V = 1 mag, so N(CO) = 1 × 10¹⁷ cm⁻² and N(¹³CO) = 1.7 × 10¹⁵ cm⁻². The resulting effective densities are 15, 60, 800, and 3800 cm⁻³ for CO J = 1 → 0, CO J = 2 → 1, ¹³CO J = 1 → 0, and ¹³CO J = 2 → 1, respectively. The comparable density for the Herschel sample is less clear. The structures were identified from dust continuum emission, but NH₃ observations of both (1, 1) and (2, 2) inversion transitions were required in order to estimate the kinetic temperature, T_K. The effective density for the (2, 2) line at T_K = 15 K is n_eff = 1.6 × 10⁴ cm⁻³ (Shirley 2015), similar to the density of about 2 × 10⁴ cm⁻³ needed to make T_K close to the dust temperature, T_d, as found by Merello et al. (2019). We take n_eff = 1.6 × 10⁴ for the Herschel sample.

The values for f(M) are plotted versus these effective densities in Figure 12. We use only the results for entries 3, 7, 9, 12, and 13, along with the average of 10 and 11 (the EXFC samples), from Table 2 as most representative of structure selection by each tracer. Clearly the fraction of the mass in bound structures increases rapidly with n_eff, reaching values about 0.8 for densities near the effective density of ¹³CO J = 2 → 1 and unity for the Herschel sample. The data roughly follow a line in the semilog plot, but the two points for ¹³CO J = 1 → 0 suggest that location in the Galaxy also matters.

4.3. Dependence on Location in the Galaxy

The contrast between the inner Galaxy (CMZ; Oka et al. 2001) and the outer Galaxy (OGS; Heyer et al. 2001) is perhaps the most striking. The structures identified in the CMZ have very high surface densities ( $\left\langle {{\rm{\Sigma }}}_{M}\right\rangle =1700$ M_⊙ pc⁻²), but f(M) = 0.71. By contrast, the structures in the OGS have very low surface densities ( $\left\langle {{\rm{\Sigma }}}_{M}\right\rangle =17.7$ M_⊙ pc⁻²), but f(M) = 0.58 for entry 6, with the more inclusive lower limit on cloud mass. With the criterion of Σ_M < 120 M_⊙ pc⁻², none of the structures in the CMZ defined by CO are bound (both f_c(N) and f_c(M) are zero), while f_c(M) = 0.46 or 0.56 in the OGS, depending on the mass limits used. Clearly, the local criterion of Σ_M = 120 M_⊙ pc⁻² for boundedness, based on a star formation threshold, is inappropriate for these parts of the Galaxy. In the CMZ, a clump would need to be defined with a much higher Σ_M, higher than 3000 M_⊙ pc⁻². In the outer Galaxy, the opposite applies.

A difficulty with this comparison is that very different criteria were used to define structures in the CMZ compared to the disk or outer Galaxy. A high threshold was needed to isolate structures in the highly confused CMZ. Using the same procedure for all R_gal, Miville-Deschênes et al. (2017) did not find a drop in median of α_vir at large R_gal, but did see a rise inside a few kiloparsecs. The sample of other galaxies (Sun et al. 2020) also uses a common method for all radii. There are 1562 data points in their catalog with R_gal < 0.5 kpc, similar to the definition of the CMZ in the Milky Way. In that restricted sample, the median, mean, and standard deviation of Σ_M are 130, 375, and 670 M_⊙ pc⁻², much higher than for the full sample, but still less than those for the CMZ. The median, mean, and standard deviation of α_vir are 5.45, 8,28, and 9.03, about twice those for the full sample, and larger than those for the CMZ, while the fractions in bound structures are f(N) = 0.11 and f(M) = 0.27, substantially smaller than for the full sample or for the CMZ of the Milky Way. The data from other galaxies supports the idea that a much higher surface density is needed to produce bound structures in the inner regions of galaxies and suggests that most of the molecular gas there is unbound despite very high surface densities.

Comparing the GRS and EXFC samples in ¹³CO J = 1 → 0 reveals a similar pattern. More mass is in bound structures despite the lower Σ_M in the outer parts of the Galaxy.

4.4. Dependence on Mass

A trend common to all of the samples is a decrease in α_vir with increasing mass. At some level, this trend is a result of a trend of decreasing α_vir with increasing Σ_M. This is most obvious for the sample of Sun et al. (2020), where the size of the region considered is constant, so the trend with mass is really a trend with surface density.

Contrary evidence that the most massive clouds may be unbound can be seen in Figure 8 from Dame et al. (1986), who consider the 26 most massive complexes in the first quadrant of the Milky Way. About half have line widths larger than those expected from virial equilibrium. A more recent analysis of the most massive complexes by Nguyen-Luong et al. (2016) found that most were unbound. While some are engaged in what they called "mini-starbursts" that could be unbinding the cloud, those sources do not have a very different distribution of virial parameters. It is these very large complexes that will dominate studies of CO in other galaxies. They also dominate the total mass of molecular gas in our Galaxy because the cloud mass function slope is less than 2, unlike the case for stars. However, the mass function peaks at much lower masses in the outer Galaxy, as shown in Figure 7 of Miville-Deschênes et al. (2017). A similar trend is seen for molecular clouds in M33 (Rosolowsky 2005).

4.5. Dependence on Mass Surface Density

The decrease in median and mean α_vir with increasing Σ_M is clear from the relevant plots for most samples that cover the relevant range of Σ_M. The distinction between clouds and clumps adopted for discussion in this paper is Σ_M = 120 M_⊙ pc⁻². To explore this connection further, we have calculated the mass fraction in bound structures in bins of 0.3 in log Σ_M for samples that cover the transition well. As shown in Figure 13, these all show a strong increase in the fraction of mass that is bound as Σ_M increases.

**Figure 13.** The value of f(M) is plotted vs. log Σ_M for the CO samples of Sun et al. (2020), Miville-Deschênes et al. (2017), Heyer et al. (2001), and the ¹³CO J = 1 → 0 sample for the GRS, with structures identified by CO (entry 9 in Table 2).
Download figure:
Standard image High-resolution image

The dependence on surface density is also reflected in the values of f_c(N) and f_c(M) in Table 2. These are generally smaller than the f(N) and f(M) computed for all structures. Only a small fraction of the gas in structures we have called clouds is bound. This statement applies both to Milky Way structures and to those in other galaxies.

The decline of virial parameter with mass may in part be an artifact of observational biases in the way that the linewidth is measured (Traficante et al. 2018). They point out that the tracers used for determining the velocity dispersion may be weighted more heavily to denser regions, with smaller dispersions that do not reflect the full dispersion of the entire structure. This effect would be most noticeable in samples for which the tracers of mass and velocity dispersion are different. The effect might contribute to the very low values of α_vir and the strong trend with mass in the Herschel sample.

4.6. Caveats and Uncertainties

As discussed extensively in Section 1, there are substantial uncertainties in computing the virial parameter. For any given structure, the uncertainty in α_vir is likely to be a factor of 2 or 3. Analyses of the uncertainties in estimating cloud mass, velocity dispersion, etc. can be found in Beaumont et al. (2013). Masses from CO luminosities have uncertainties, as discussed in the introduction, so that we expect uncertainties of at least a factor of 3, and these may be systematic uncertainties, at least in the case of M_CO. Gong et al. (2018) argue that X(CO) may be overestimated based on "observing" simulations. They get a value of X(CO) about half the standard value of 2 × 10²⁰. If that is correct, masses from CO would be overestimated and virial parameters underestimated. Gong et al. (2020) has further examined this issue and provided conversion factors that are shallower functions of metallicity (Z) than those used by Sun et al. (2020). Finally, we re-emphasize the dependence on the method used to define a cloud.

For ¹³CO, our newer analysis accounts for the latest measurements of the isotopic ratio, CO/¹³CO, as a function of R_gal and the latest value of H₂/CO. However, we have made no correction for changes in the CO abundance with metallicity, which is known to decrease with R_gal. Deharveng et al. (2000) found a gradient in the oxygen abundance over the range of 5–15 kpc: 12 + logO/H = (−3.95 ± 0.49) × 10⁻² R_gal + (8.82 ±0.05). The effect of metallicity on CO emission is not simple, but extragalactic observers commonly use a correction factor: α_CO = 4.35Z^−1.6 M_⊙ (K km s⁻¹ pc²)⁻¹ (Sun et al. 2020 and references therein), where Z is the metallicity relative to solar. If we applied the same correction factor, CO clouds at 3 kpc would have twice higher α_vir, while those at 15 kpc would have about one-third the value of α_vir. Such a correction would strengthen the already strong trend for clouds to be unbound in the inner Galaxy, where most of the molecular mass resides, and to be bound in the outer Galaxy. However, the divergence of opinion on how α_CO depends on Z is currently substantial (see, e.g., Accurso et al. 2017; Gong et al. 2020; Lada & Dame 2020).

We reiterate that methods of structure identification remain one of the most significant sources of uncertainty.

5. Consequences for the Star Formation Rate Problem

The results of this analysis have consequences for the problem of the star formation rate, discussed in the 1. If only a fraction of the mass of molecular gas in the Galaxy is in bound structures, the over-prediction of the star formation rate of 300 M_⊙ yr⁻¹, can be decreased substantially. Based on the catalog of Miville-Deschênes et al. (2017), the only one that accounts for all of the CO emission, only 19% of structures defined by CO J = 1 → 0 are bound, decreasing the theoretical star formation rate to 57 M_⊙ yr⁻¹, using the same assumption about the freefall time as the value found in Section 1.

We can make a more accurate estimate by assessing the star formation rate cloud by cloud, then adding those rates up. If we assume that only clouds with α_vir ≤ 2 form stars in a freefall time based on the mean density of that cloud, the predicted star formation rate for an individual cloud is

$\begin{eqnarray}&&{SFR}_{{\rm{t}}{\rm{h}}}=M/{t}_{{\rm{f}}{\rm{f}}}=6.04\times {10}^{-8}{M}_{{M}_{\odot }}^{1.5}{r}_{{\rm{p}}{\rm{c}}}^{-1.5}.\end{eqnarray} \tag{ 10 }$

Using the catalog of Miville-Deschênes et al. (2017) because it accounts for all of the molecular clouds defined by CO emission, and summing the star formation rate only over bound clouds, we find a total star formation rate for the Milky Way of 46.4 M_⊙ yr⁻¹, a factor of 6.5 less than the value in the 1. This is mostly due to the small fraction of bound clouds, but there is also a contribution from a longer t_ff because the mean density of the clouds identified by Miville-Deschênes et al. (2017) is less than 100 cm⁻³. The resulting depletion time of bound molecular gas is decreased to t_dep = 1.0 × 10⁸ yr, closer to the freefall time, but still considerably larger.

Of course, regions of clouds are likely to be bound even if the entire cloud is not. The fraction of bound mass could be assessed from comparison of masses indicated by different tracers. Ideally, we would do the same analysis for a complete survey of the Galaxy in other tracers, especially ¹³CO. A similar calculation for the ¹³CO GRS survey, using the full sample (entry 8 in Table 2) predicts a star formation rate of 7.4 M_⊙ yr⁻¹. However, the existing ¹³CO surveys do not cover the whole Galaxy, and the new analysis was possible only for a subset of the ¹³CO data, so we cannot do a direct calculation. Instead we calculate the fraction of the mass of the unbound CO-defined structures that is found in the bound ¹³CO-derived properties. This value is 0.22 for the GRS, which applies to the inner Galaxy, which contains most of the molecular mass. The fraction is higher (0.77–0.80) in the outer parts of the Galaxy. So this material would increase the predicted star formation rate. If all of the ¹³CO-derived mass in the unbound CO-defined clouds in the GRS sample form stars at the fiducial t_ff, they would contribute another 55 M_⊙ yr⁻¹ to the Galactic star formation rate. Clearly, a full sky survey of ¹³CO would be very valuable. We have focused on _ff, rather than the final efficiency (_sf), but the simulations do find lower final efficiencies (_sf) for higher α_vir (a factor of 5–10 decrease as α_vir increases from 1 to 5; Kim et al. 2021). Because the most massive clouds seem most likely to be bound (subject to the caveats above), the difference between mass functions of clouds and clusters is, if anything, increased. The resolution to this problem may lie in the inclusion of tidal forces since most of these very massive clouds reside in the inner Galaxy.

6. Conclusions

We list the main conclusions before discussing them:

1.
Clearly, the choice of method to identify clouds plays a major role in the different results.
2.
The tracer (both molecule and transition) that is used to define structures has a strong effect on whether those structures are bound.
3.
At a fixed R_gal, the fraction of the mass in bound structures increases rapidly with the effective density of the tracer, suggesting that regions with larger volume densities contribute more to bound gas.
4.
Many clouds identified by CO are unbound and most of the mass in those clouds is unbound.
5.
Structures defined by dust continuum emission and line widths from NH₃ emission are almost all bound.
6.
Structures defined by ¹³CO emission are more likely to be bound than those defined by CO emission. However, their boundedness also depends on the method used to identify them and on assumptions about isotopic and elemental abundances.
7.
Boundedness correlates strongly with surface density.
8.
Structures within 0.5 kpc of a galaxy center have much higher surface density, but much of the mass is unbound because of higher turbulence. The opposite is true in outer regions.
9.
More generally, structures identified by CO or ¹³CO have lower surface density, but are more likely to be bound, at larger R_gal.
10.
More massive clouds are more likely to be bound, but this trend is partially confused with the trend with surface density.
11.
For the most complete catalog of structures traced by CO (Miville-Deschênes et al. 2017), the fraction of mass in bound structures, f(M) = 0.19, alleviating, but not eliminating, the fundamental problem of slow star formation in the Milky Way.

The most important conclusion is that most of the mass traced by CO emission is in unbound structures. This result was anticipated by theorists, notably Dobbs et al. (2011), whose title asked "Why are most molecular clouds not gravitationally bound?". Their answer involved cloud−cloud collisions and stellar feedback; shredding and merging resulted in clouds not being well-defined entities with the same gas over cloud "lifetimes." Padoan et al. (2017) argued that supernova feedback kept the interstellar medium sufficiently stirred up that many structures were unbound even though their densities were such that they were likely molecular. As mentioned in Section 1, Kim et al. (2021) found that unbound clouds are needed to reach the observational constraints on the star formation rate.

Observers have often found, but seldom emphasized, this result. Roman-Duval et al. (2016) suggested that gas detected in CO, but not ¹³CO might be "diffuse, non-star-forming gas," but did not explicitly suggest that this gas was unbound. Figure 16 of Miville-Deschênes et al. (2017) already showed that almost all of the structures in their catalog are unbound, but this result was not particularly emphasized. Nguyen-Luong et al. (2016) showed that mini-starburst complexes, which have enhanced star formation rate density, are actively forming massive stars but are not necessarily bound, suggesting that being bound is not the major factor regulating their star formation activity. The idea that structures identified by molecular emission are well-defined, bound structures, unlike the more diffuse atomic interstellar medium, is long-standing and creates cognitive dissonance with contrary factual evidence. Uncertainties are often invoked to argue that clouds are bound even when the evidence suggests otherwise (see Dobbs et al. 2011 for discussion of this tendency.). However, there is no good reason to believe that the atomic−molecular transition is identical to the unbound−bound transition. This would be plausible if most of the support were thermal, but it is turbulent (Zuckerman & Evans 1974). The fact that these transitions are close can be coincidental (see page 36 of Elmegreen 1985) and may have led us astray.

We acknowledge stimulating discussions with many colleagues, but especially with E. Ostriker and J.-G. Kim. N.J.E. acknowledges C. Dobbs, A. Burkert, and J. Pringle for their paper that opened his eyes to the possibility that molecular clouds are unbound. The authors acknowledge Paris-Saclay University's Institut Pascal program "The Self-Organized Star Formation Process" and the Interstellar Institute for hosting discussions that nourished the development of the ideas behind this work.

Software: astropy (Astropy Collaboration et al. 2013, 2018), GILDAS (Pety 2005; Gildas Team 2013).

Appendix

We have selected a subset of clouds from the Miville-Deschênes et al. (2017) catalog that fit into the coverage of the ¹³CO J = 1 → 0 surveys from the Five College Radio Astronomy Observatory in order to investigate the boundedness of these clouds using a mostly optically thin tracer of H₂ column density. To connect the ¹³CO emission to the structures defined by the sum of Gaussian components (Miville-Deschênes et al. 2017), we developed the following procedure. For a given cloud, we calculate the fraction of the observed CO emission that is recovered by the summed Gaussian for each position and velocity channel, F(X, Y, V). A typical profile of F is flat-topped with values ranging from 0.9 to 1 in the core of the line (see Figure 14). The fractional values rapidly decrease toward the wings but remain nonzero for the full extent of the spectral axis. To exclude spectral channels beyond the core and wing of a cloud, we set F(X, Y, ∣v − v_cen∣ > dv) = 0 where v_cen and dv are the velocity centroid and FWHM from the catalog of Miville-Deschênes et al. (2017). A ¹³CO data cube, T₁₃(x, y, v), is constructed from the FCRAO surveys at the native 48'' angular resolution and 22'' sampling interval to match the CO coverage of the cloud. The spectral axis of the ¹³CO data is smoothed and resampled to 0.25 km s⁻¹ wide channels. At this stage the ¹³CO data cube may contain emission features not associated with the CO-defined cloud at velocities outside of the core of the CO emission (see e.g., the left panel of Figure 14). To further isolate the ¹³CO emission, we first resample the spectral axis of the fraction cube, F, to be aligned with the ¹³CO spectral axis and then multiply the ¹³CO spectral axis with the corresponding fractional spectrum,

$\begin{eqnarray}&&{T}_{13}(x,y,v)^{\prime} =F(X,Y,v){T}_{13}(x,y,v)\end{eqnarray} \tag{ A1 }$

where X − dX/2 < x < X + dX/2, Y − dY/2 < y < Y + dY/2, and dX, dY are the 0 fdg 125 pixel sizes of the CO data cube. The effect of this multiplication is to suppress any contribution to the integrated ¹³CO emission from signal outside of the CO velocity interval, as shown in Figure 14.

**Figure 14.** The method used to filter unrelated signal from ¹³CO. The spectrum on the left is the original data. The center panel shows the function that isolates the velocity range that is consistent with that of the CO-defined structure. The right panel shows the filtered ¹³CO spectrum that contributes to the ¹³CO-defined mass and other properties.
Download figure:
Standard image High-resolution image

The cloud mass is computed from N(H₂) found from Equation (10) in Section 1 as follows:

$\begin{eqnarray}&&M{(}^{13}\mathrm{CO})={\mu }_{{{\rm{H}}}_{2}}{m}_{{\rm{H}}}{D}^{2}\int d{\rm{\Omega }}N({{\rm{H}}}_{2})(l,b)\end{eqnarray} \tag{ A2 }$

where ${\mu }_{{{\rm{H}}}_{2}}=2.8$ , m_H is the atomic hydrogen mass, and D is the distance to the cloud. For each cloud, the velocity dispersion, σ_v(¹³CO), and cloud radius are calculated following Miville-Deschênes et al. (2017). Both σ_v(¹³CO) and radius are deconvolved from the native spectral and angular resolutions of the data. We also deconvolved the CO velocity dispersions and cloud sizes. Clouds are excluded from any subsequent analyses if either the cloud size or velocity dispersion is not resolved. Finally, the ¹³CO-derived virial mass is derived from Equation (8). Accounting for these selection criteria, there are 289 clouds from the GRS, 105 clouds from the EXFC55-100 survey, and 98 clouds from the EXFC135-195 survey.

Which Molecular Cloud Structures Are Bound?

Article metrics

Author e-mails

Author affiliations

ORCID iDs

Dates

Abstract