© ESO 2021

1. Introduction

There are still many gaps in our understanding of the condensation of dense structures out of the diffuse interstellar medium (ISM). There are reasons to think that the influence of the magnetic field is preponderant. In recent years, remarkable observational data have been obtained by Planck Collaboration Int. XXXV (2016) that allow a quantitative analysis of the relative orientation of the magnetic field within a set of nearby (d < 450 pc) molecular clouds. These observations have helped to establish the significance of magnetic fields in the formation of dense structures on physical scales ranging from approximately 1 to 10 pc. They show a clear correlation in the direction of elongation of high-density regions (number column density N_H ≳ 10²² cm⁻²), which appears to be perpendicular to the ambient magnetic field direction.

Molecular clouds are known to contain hierarchical nested density structures with, for example, clumps, filaments, and cores (see André et al. 2014; Dobbs et al. 2014; Heyer & Dame 2015). Pokhrel et al. (2018) studied hierarchical structure over five different scales (ranging from ≳10 pc to 10 AU) in the Perseus molecular cloud using new observations from the Herschel, the James Clerk Maxwell Telescope (JCMT), the Submillimeter Array (SMA), the Very Large Array (VLA). They compared the number of fragments with the number of Jeans masses at each scale to calculate the Jeans efficiency, which is the ratio of observed to expected number of fragments. Könyves et al. (2015) used the results of the Herschel Gould Belt survey (HGBS- André et al. 2010) in the Aquila molecular cloud complex, and compared the numbers of prestellar cores in various density bins to the number of young stellar objects (YSOs). They estimated that the lifetime of prestellar cores is ∼1 Myr, which is about four times longer than the core free-fall time, and that the lifetime decreases as the average core density increases. While current observations cannot determine whether ambipolar diffusion (neutral-ion slip) is occurring during the initiation of gravitational collapse, nonideal magnetohydrodynamic (MHD) simulations suggest it plays an important role in establishing mildly supercritical mass-to-flux ratios in prestellar cores, whether starting from small-amplitude perturbations (see Kudoh et al. 2007; Basu et al. 2009a) or from large-scale turbulent or converging flows (see Nakamura & Li 2005; Basu et al. 2009b; Chen & Ostriker 2014). Furthermore, in mildly supercritical regions, the hybrid modes driven by gravity and neutral-ion slip result in preferred length scales and growth times that can significantly exceed the Jeans scale and free-fall time, respectively (see Basu & Ciolek 2004; Ciolek & Basu 2006; Bailey & Basu 2012).

In this Letter, we probe the variation of lifetime and fragmentation scales of dense structures in molecular clouds as a consequence of ambipolar diffusion. We apply the results of a linear analysis of ambipolar-diffusion-driven fragmentation in planar, isothermal, weakly ionized, self-gravitating sheetlike magnetic clouds. The calculated shortest growth timescale and preferred fragmentation mass for collapse are used to explain the observationally estimated lifetime of prestellar cores and the number of enclosed cores in a parent clump. In Sect. 2 we describe the analytic model, in Sect. 3 we present a comparison with the observational findings, and finally, in Sect. 4 we summarize and draw conclusions from our results.

2. Analytic model

We model interstellar molecular clouds as self-gravitating, partially ionized, isothermal, magnetic, planar, thin sheets with infinite extent in the x- and y-directions and a local vertical half-thickness Z(x, y, t). Although sheets are an idealized geometry, a structured background state like a sheet or filament does capture the essential feature of a preferred scale of instability that is related to the local density scale length; this would not appear if assuming a uniform background. We note that the critical length scale and timescale can differ typically by a factor of approximately two as compared to using a spherical (uniform) geometry. Our static initial state also does not include the effect of any large-scale motion that modifies the evolutionary timescale. The nonaxisymmetric equations and formulations of the model have been described in detail in several papers (Ciolek & Basu 2006; Bailey & Basu 2012). The evolution equations include the nonideal MHD effect of ambipolar diffusion, the process by which neutrals are partially coupled to magnetic field through collisions with ions. This effect is quantified by the neutral-ion collision (momentum-exchange) timescale (e.g., Basu & Mouschovias 1994):

$$\begin{array}{c}\hfill {\tau }_{\mathrm{ni}}\equiv 1.4\frac{m_i+m_n}{m_i}\frac{1}{n_i{⟨\sigma w⟩}_{i{\mathrm{H}}_2}},\end{array}$$(1)

where ⟨σw⟩_iH2 is the average collision rate between ions of mass m_i (singly ionized Na, Mg, and HCO, adopted to have a mass of 25 amu) and neutrals of mass m_n(= 2.33 amu). We adopt a neutral–ion collision rate between H₂ and HCO⁺ of 1.69 × 10⁻⁹ cm³ s⁻¹ (McDaniel & Mason 1973). For the ion number density n_i, there is an assumed constant power-law approximation of the form $n_i\propto n_{n,0}^{1/2}$ (Ciolek & Basu 2006; Ciolek & Mouschovias 1998), where n_n, 0 is the initial uniform number density of neutrals. The typical observed ionization fraction in molecular clouds (primarily due to cosmic ray ionization) (Elmegreen 1979; Tielens 2005) is

$$\begin{array}{c}\hfill {\chi }_i={10}^{-7}(\frac{n_{n,0}}{{10}^4{\mathrm{cm}}^{-3}}{)}^{-1/2}.\end{array}$$(2)

Such low ionization means that ambipolar diffusion is unavoidable in molecular clouds, but there is still enough coupling with the charged species for the neutrals to be affected by the magnetic field. This is because of the high polarizability of the neutrals, particularly H₂ molecules (Osterbrock 1961). The threshold for whether a region of a molecular cloud is magnetically dominated or gravitationally dominated is given by the normalized mass-to-flux ratio,

$$\begin{array}{c}\hfill {\mathrm{\mu }}_0\equiv 2\pi G^{1/2}\frac{{\sigma }_{n,0}}{B_{\mathrm{ref}}},\end{array}$$(3)

where σ_n, 0 is the initial uniform column density of the sheet, B_ref is the magnetic field strength of the background reference state, G is the universal gravitational constant, and (2πG^1/2)⁻¹ is the critical mass-to-flux ratio for gravitational collapse (Nakano & Nakamura 1978). Regions with μ₀ < 1 are defined as subcritical, regions with μ₀ > 1 are defined as supercritical, and regions with μ₀ ≈ 1 are defined as transcritical. For small amplitude perturbations, the governing equations can be combined to yield the following dispersion relation:

$$\begin{array}{cc}\hfill (\omega +i\theta )[{\omega }^2-C_{\mathrm{eff},0}^2{k}^2+& 2\pi G{\sigma }_{n,0}k]=\hfill \\ \hfill & \omega [2\pi G{\sigma }_{n,0}k{\mathrm{\mu }}_0^{-2}+k^2{V}_{A,0}^2],\hfill \end{array}$$(4)

where

$$\begin{array}{c}\hfill \theta ={\tau }_{\mathrm{ni},0}(2\pi G{\sigma }_{n,0}{\mathrm{\mu }}_0^{-2}k+k^2{V}_{A,0}^2)={\eta }_{\mathrm{AD},0}(k{Z_0}^{-1}+k^2),\end{array}$$(5)

ω is an angular frequency, and $k^2\equiv k_z^2=k_x^2+k_y^2$, where k_x, k_y, k_z (or k) are the wavenumbers in the x-, y-, and z- directions, respectively. To obtain the dispersion relation as shown in Eq. (4), the linearized perturbed MHD equations are used (see Appendix A). Here, C_eff, 0 and V_A, 0 are the local effective sound speed and the Alfvén speed, respectively. The term C_eff, 0 includes the effects of a restoring force due to an external pressure P_ext. The dimensionless external pressure ${\stackrel{\sim }{P}}_{\mathrm{ext}}$ ($\equiv 2P_{\text{ext}}/(\pi G{\sigma }_{n,0}^2)$) and temperature (T) are kept fixed at 0.1 and 10 K, respectively.

Figure 1 presents the normalized shortest growth timescale, τ_g, m/t_ff, and normalized preferred fragmentation mass, M_g, m/M_T, m, corresponding to this minimum timescale as a function of μ₀ for the case of normalized neutral–ion collision time ${\tau }_{\text{ni},0}/t_0={\stackrel{\sim }{\tau }}_{\text{ni},0}=0.2$ which itself corresponds to Eq. (2). Here, t_ff is the dynamical, i.e., free-fall time ( = Z₀/c_s), M_T, m is the preferred thermal mass based on our model, and $M_{T,m}={(4\pi C_{\text{eff},0}^2/c_s^2)}^2{M}_0$, where M₀ is effectively the Jeans mass. See Appendix A for definitions and typical values of the units of time (t₀), length (L₀), and mass (M₀), and other parameters. For all objects shown in Table B.2, we calculate M_g, m in units of M_c = πσ_n, 0(L₀/2)² = πM₀/4, as the perturbation is taken to be circular with radius L₀/2. As the mass-to-flux ratio goes to the subcritical regime where ambipolar diffusion drives the evolution, the curve of shortest growth timescale approaches a plateau. It is noteworthy that the peak preferred fragmentation mass for collapse exceeds the Jeans mass by a factor of up to ten. Furthermore, for ${\stackrel{\sim }{\tau }}_{\text{ni},0}=0.2,$ the timescale for collapse of a subcritical region is around 10 − 12 times longer than that of a supercritical region; this is the origin of the often-quoted result that the ambipolar diffusion time is approximately ten times the free-fall time (Mouschovias 1991; Ciolek & Basu 2006; Bailey & Basu 2012; Das & Basu 2021). Hereafter, for a better representation, we use n_n, σ_n, and μ instead of n_n, 0, σ_n, 0, μ₀, respectively. In this study, we are interested in the regime where the normalized mass-to-flux ratio μ remains in the range 1 ≲ μ ≲ 2. See Kunz & Mouschovias (2009) for an application of the linear theory in the subcritical regime to model core masses.

Fig. 1. Growth timescale τg, m (in the units of dynamical or free-fall time tff) and the preferred fragmentation mass Mg, m (in the units of preferred thermal fragmentation mass MT, m) of the most unstable mode as a function of the normalized mass-to-flux ratio (μ0). The model utilizes the standard ionization fraction corresponding to a normalized neutral–ion collision time ${\tau }_{\mathrm{ni},0}/t_0={\stackrel{\sim }{\tau }}_{\text{ni},0}=0.2$. The shaded and unshaded zones represent the subcritical (μ0 < 1) and supercritical (μ0 > 1) regimes, respectively. The present paper focuses on the (mildly) supercritical regime.

3. Observational correspondence to prestellar cores

In this section we discuss the relevance of our theoretical results to observational findings, focusing on the mildly supercritical regime.

3.1. Lifetime of prestellar cores

The technique for finding the timescale of the core-formation process was introduced by Beichman et al. (1986) in the context of IRAS sources. These authors studied the embedded YSOs within the core sample of Myers & Benson (1983) and Myers et al. (1983), and found that 35 cores had IRAS sources meeting the color-selection criteria of embedded YSOs and 43 had no embedded IRAS sources. Beichman et al. (1986) calculated the percentage of cores with embedded sources to estimate the lifetime of a core without an embedded YSO by comparing it with the lifetime of the embedded YSO phase. Using an estimated lifetime of cores with embedded class II sources of 1 − 2 Myr (as discussed in Ward-Thompson et al. 1994; Evans et al. 2009), and assuming that the prestellar cores go on to form protostars, the prestellar core lifetime is estimated to be

$$\begin{array}{c}\hfill \tau =\frac{\#\mathrm{of}\mathrm{cores}\mathrm{without}\mathrm{embedded}\mathrm{sources}}{\#\mathrm{of}\mathrm{cores}\mathrm{with}\mathrm{embedded}\mathrm{sources}}\times [1-2]\mathrm{Myr}.\end{array}$$(6)

This formula was used by Jessop & Ward-Thompson (2000) on a catalog of molecular cloud cores from the all-sky IRAS Sky Survey Atlas (ISSA), and their Fig. 6 (often called a “JWT” plot) shows the estimated lifetime versus mean density. In a similar way, Könyves et al. (2015) estimated the lifetime of prestellar cores (see their Fig. 9) by comparing the number of prestellar cores found with Herschel to the number of Class II YSOs detected by Spitzer in the Aquila cloud. This study had the advantage of considering a homogeneous sample of cores from a single cloud (Aquila), measured using a single telescope, tracer (dust), and analysis technique that separated prestellar (gravitationally bound) cores from starless (unbound) cores. Some of the underlying assumptions in these studies are: (1) that prestellar cores will evolve into YSOs in the future; and (2) that star formation proceeds at a roughly constant rate, at least when averaged over an entire cloud.

Figure 2 presents the core lifetime (values in green filled squares) estimated from our model as the instability growth time at a particular density n_n, when adopting a specific model for B_ref as a function of n_n (see Fig. 4 and Eq. (7) below). The corresponding estimated numbers based on the observations of candidate prestellar cores identified with Herschel (Könyves et al. 2015) and the literature data (Ward-Thompson et al. 2007) are also shown. We achieve the good correspondence by varying only one relationship, namely that between the normalized mass-to-flux-ratio μ and the number density n_n, with μ ranging from about 1.1 to 1.5 for n_n in the range of 10⁴ − 10⁶ cm⁻³. To calculate μ we evaluate σ_n using Eq. (A.5). These values of μ are in the range of mildly supercritical values generally obtained from Zeeman detections and use of the Davis-Chandrasekhar-Fermi (DCF) method (Crutcher 2012; Pattle & Fissel 2019). The ionization level is set by the value ${\stackrel{\sim }{\tau }}_{\text{ni},0}=0.2$ corresponding to the standard value set by Eq. (2). Our use of the instability growth time as a proxy for the evolutionary time is similar to the commonly used comparison of the free-fall time at a particular density with the evolutionary time of a core at that density. Figure 2 shows that the typical lifetime of prestellar cores decreases from ∼1.37 Myr for cores with a number density of ∼10⁴ cm⁻³ to ∼0.35 Myr at ∼10⁵ cm⁻³, and to ∼0.05 Myr at ∼10⁶ cm⁻³ (see Table B.1). The timescale for collapse of a core with a number density of ∼10⁴ cm⁻³ is around 10 − 50 times longer than that of a highly dense core of ∼10⁵ − 10⁶ cm⁻³. The highest density cores are significantly supercritical and evolve essentially on a free-fall timescale (t_ff = Z₀/c_s) as shown in Fig. 1. The lower black dashed line presents the dynamical timescale, or free-fall time, as a reference line.

Fig. 2. Estimated lifetime as a function of number density (nn). The green filled squares are obtained from our model for a normalized neutral–ion collision time ${\stackrel{\sim }{\tau }}_{\text{ni},0}=0.2$. The blue filled triangles show the corresponding data points for the population of 446 candidate prestellar cores identified with Herschel in the Aquila cloud (Könyves et al. 2015). Earlier data from Ward-Thompson et al. (2007) are shown by the red circles for better comparison. The black dashed line shows the dynamical time (Z0/cs), i.e., the free-fall time, of our model.

3.2. Number of enclosed cores

Sadavoy et al. (2010) used point-source photometry to explore the dense cores in the Perseus star-forming complex as found in surveys with SCUBA (85 μm) and Spitzer Space Telescope (3.6 − 70 μm). Mercimek et al. (2017) characterized the distribution of these cores inside the clumps. Pokhrel et al. (2018) analyzed the submillimeter starless or protostellar cores in the Herschel column density maps of Mercimek et al. (2017) and used the estimated mass and areas (see Table B.2) to determine the average density of each clump for A_V > 7 mag. These authors used the dust temperatures from Sadavoy et al. (2014) to estimate the thermal support. To calculate the Jeans number (N_J, CLUMP ≡ M/M_J, th, i.e., the number of contained thermal Jeans masses M_J, th) of the clumps, Pokhrel et al. (2018) used the line-of-sight averaged temperatures and mass derived in Sadavoy et al. (2014).

We fit the observations of the number of enclosed cores (Num_CORE) in each clump as a function of the corresponding Jeans number of the clumps, as seen in Fig. 3. We calculate Num_CORE in the context of our model by dividing the total clump mass M by the preferred fragmentation mass M_g, m, adopting ${\stackrel{\sim }{\tau }}_{\text{ni},0}=0.2$. We note that M_g, m significantly exceeds M_J, th for mildly supercritical objects. The Jeans number N_J, CLUMP is the expected number of cores in the context of thermally regulated fragmentation. See Table B.2 for detailed specifications of all the clumps. The clumps are arranged in an increasing order of Num_CORE (i.e., number of cores). Figure 3 shows that the number of cores (blue filled circles) increases with the Jeans number of the clumps, as shown in Fig. 6 of Pokhrel et al. (2018). The filled green crosses represent our model values. The black dashed line shows the efficiency of unity (ϵ^th = 1, i.e., Num_CORE = N_J, CLUMP) corresponding to purely thermally regulated fragmentation.

Fig. 3. Number of enclosed cores (log[NumCORE]) as a function of Jeans number (log[NJ, CLUMP]) of clumps. The blue filled circles represent the observed number of enclosed cores as found by Pokhrel et al. (2018). The green filled crosses are obtained from our model for a normalized neutral–ion collision time ${\stackrel{\sim }{\tau }}_{\text{ni},0}=0.2$. The black dashed line corresponds to an efficiency of unity.

The observations show fewer cores than that predicted with only thermal pressure. This hints at a larger threshold for fragmentation mass that includes effects beyond that set by thermal pressure alone. In our model, M_g, m serves such a purpose as a magnetic-field-dependent instability threshold in contrast to a Jeans mass.

It is worth noting that the number of fragmented cores is also dependent on the clump mass. The relatively massive clumps are able to generate more cores and Num_CORE comes closer to (but stays below) the value of N_J, CLUMP. Therefore, the number density might not be the only key parameter in this context. The nonthermal motions in these massive clouds could also play a role. However, using the nonthermal dispersion in the calculation of the Jeans mass would not fit the observations for these clouds either, predicting many fewer cores than are observed (see the discussion in Pokhrel et al. 2018).

The inability of the observed nonthermal dispersion to be used as the source of an internal pressure when estimating the number of fragmented cores could possibly be attributed to the following reasons. First, the nonthermal motions may arise at least in part due to large-scale velocity gradients, in which case they cannot act as an effective pressure. Simulations in turbulent boxes that have global stability show that small-scale collapse still occurs unless the driving scale and/or power spectrum is peaked at extremely short scales (e.g., Vazquez-Semadeni et al. 1996; Klessen et al. 2000). Second, the nonthermal dispersion may be dominated by motions in the lowest densities that are traced, while a dense layer that actually undergoes fragmentation may have a lower dispersion that is closer to the thermal value. The latter explanation is supported by simulations (Kudoh & Basu 2003, 2006; Folini et al. 2004) that show that the velocity dispersion peaks in the low-density regions of a stratified molecular cloud. In observations as well, the velocity dispersion peaks in low-density regions and starless dense cores correspond to minima in velocity dispersion maps (see, e.g., Fig. 13 of Friesen et al. 2017).

3.3. Estimation of magnetic field

For gravitationally contracting fragments (particularly cores) in magnetic interstellar clouds, the relation between magnetic-field strength B_ref and gas number density n_n is of considerable interest. In Fig. 4, we present log(B_ref/μG) (which is obtained from our model) as a function of log(n_n/cm⁻³) for the two different density regimes shown in Figs. 2 and 3 based on our model for ${\stackrel{\sim }{\tau }}_{\text{ni},0}=0.2$. We use the least-squares method to find a best fit to the data:

$$\begin{array}{c}\hfill \frac{B_{\mathrm{ref}}}{10\mathrm{\mu }\mathrm{G}}=A{\left(\frac{n_n}{{10}^4{\mathrm{cm}}^{-3}}\right)}^{\kappa },\end{array}$$(7)

Fig. 4. Magnetic field log(Bref/μG) vs. number density $log(n_n/{\mathrm{cm}}^{-3})$. The filled squares and crosses represent the number density region corresponding to the fitting of the lifetime of prestellar cores and NumCORE as shown in Figs. 2 and 3 (see Tables B.1 and B.2). The dashed line is the least-squares fit. The vertical error bars correspond to magnetic field variations that yield a total factor of two (greater or lesser by a factor $\sqrt{2}$) values of the growth time (lifetime).

where κ = 0.427 and A = 1.766, and shown by the black dashed line in Fig. 4.

Mestel (1965) argued that for a spherically and isotropically collapsing magnetic cloud, the scaling between the magnetic field strength and the density is $B\propto n_n^{2/3}$. This argument is true for a quasi-spherical collapse where both the mass M and the magnetic-flux Φ are being conserved, and the magnetic field energy is insignificant compared to the gravitational energy. Later, Mouschovias (1976a,b) argued that the plasma β ($=8\pi \rho c_s^2/B^2$) remains constant during self-contraction of a cloud with a dynamically important magnetic field, therefore $B\propto n_n^{1/2}$, although for individual points within the cloud the exponent was in the range 1/3 − 1/2. These theories assume an evolutionary sequence and that the mass-to-flux ratio and internal thermal (or turbulent) velocity dispersion stays fixed. We should keep in mind that the full virial relation is

$$\begin{array}{c}\hfill B\propto \sigma n_n^{1/2}/\mathrm{\mu },\end{array}$$(8)

where σ is the velocity dispersion. With the inclusion of ambipolar diffusion and a systematically increasing value of μ as n_n increases, the slope in the B_ref − n_n relation is expected to be less than in the flux-freezing models. On the other hand, a systematic dependence of σ on n_n can also have an effect. Crutcher et al. (2010) fitted a slope of ≈2/3 to Zeeman magnetic field data for an ensemble of dense molecular gas clouds. The B − n relation shown in Fig. 4 is for prestellar cores while the B − n relation shown by Crutcher et al. (2010) includes a number of massive protostellar cores or clumps. Those objects definitely do not represent an evolutionary sequence, with the higher density objects representing much more massive clouds that also have an increased velocity dispersion σ. A better fit to the ensemble of different clouds measured by the Zeeman effect is B ∝ σn^1/2 (Basu 2000) and was verified by Li et al. (2015) who used the updated data in Crutcher (2012) and pointed out that the ensemble of objects in the data set also have the correlation σ ∝ n^1/6, thereby leading to an apparent B ∝ n^2/3. A similar result was obtained for B determined using the DCF technique (Myers & Basu 2021).

4. Conclusions

We used a semi-analytic model of ambipolar-diffusion-driven gravitational fragmentation in isothermal self-gravitating interstellar molecular clouds. The only requirements in our model are that prestellar cores are transcritical (mildly supercritical), with 1 ≲ μ ≲ 2, and an evolution toward the higher values of μ as the density increases. With this assumption we show that there is a significant and systematic variation of lifetime and fragmentation scale in molecular clouds. Such systematic variations do not exist in standard thermal-pressure-dominated (Jeans) fragmentation theory and are difficult to reproduce in turbulence-regulated models. Our best-fit model for prestellar cores suggests B ∝ n^0.43 (see Fig. 4 and Eq. (7)), which attains a shallower slope than the flux-frozen case due to the effects of ambipolar diffusion. The estimated lifetime of prestellar cores and the possible number of cores within a parent clump/cloud based on the model for ${\stackrel{\sim }{\tau }}_{\text{ni},0}=0.2$ agree well with that of observations presented by Könyves et al. (2015) and Pokhrel et al. (2018), respectively. For the lower end (∼10⁴ cm⁻³) of the density regime, the timescale for collapse of prestellar cores is approximately six times longer than the free-fall timescale. On the other hand, for the case of higher number density (∼10⁶ cm⁻³), the timescale is nearly the same as the free-fall timescale. The mass scales of fragment formation are also significantly greater than the Jeans mass in this mildly supercritical regime.

We adopted the dimensionless neutral-ion collision time ${\stackrel{\sim }{\tau }}_{\text{ni},0}=0.2$ because of its correspondence to the typical ionization fraction (∼10⁻⁷) at a neutral number density (∼10⁴ cm⁻³). In future studies, the role of varying ${\stackrel{\sim }{\tau }}_{\mathrm{ni},0}$ could be explored in order to compare with measured ionization fractions in cores (see Caselli et al. 2002), and perhaps constrain the cosmic ray ionization rate (canonical value ζ_CR = 10⁻¹⁷ s⁻¹).

Our model provides a means to indirectly infer the effect of ambipolar diffusion on mildly supercritical dense regions (prestellar cores) of molecular clouds. The importance of ambipolar diffusion in dense supercritical molecular cloud gas has not been discussed widely and is independent of its possible effects in low-density molecular cloud envelopes.

Acknowledgments

We thank the anonymous referee for comments that improved the manuscript. SB is supported by a Discovery Grant from NSERC. PhA acknowledges support from ‘Ile de France’ regional funding (DIM-ACAV+ Program) and from the French national programs of CNRS/INSU on stellar and ISM physics (PNPS and PCMI).

References

Appendix A: System of equations

We formulate model clouds as rotating, self-gravitating, partially ionized, isothermal, magnetic, planar thin sheets with infinite extent in the x- and y- directions and a local vertical half-thickness Z(x, y, t) (see Fig. 1 from Das & Basu 2021). In our model, we adopt a velocity unit of c_s, the isothermal sound speed, and a column density unit of σ_n, 0, the initial uniform column density. The length unit is $L_0=c_s^2/(2\pi G{\sigma }_{n,0})$, leading to a time unit t₀ = c_s/(2πGσ_n, 0), where G is the universal gravitational constant. The mass unit is $M_0=c_{\text{s}}^4/(4{\pi }^2{G}^2{\sigma }_{n,0})$ and the magnetic field strength unit is B₀ = 2πG^1/2σ_n, 0. The free-fall time is

$$\begin{array}{c}\hfill t_{\mathrm{ff}}=Z_0/c_{\mathrm{s}}=({\stackrel{\sim }{Z}}_0{L}_0)/c_{\mathrm{s}}={\stackrel{\sim }{Z}}_0{t}_0,\end{array}$$(A.1)

where

$$\begin{array}{c}\hfill Z_0={\sigma }_{n,0}/(2{\rho }_{n,0})\end{array}$$(A.2)

is the initial uniform local vertical half-thickness, with dimensionless form

$$\begin{array}{c}\hfill {\stackrel{\sim }{Z}}_0=\frac{2}{(1+{\stackrel{\sim }{P}}_{\mathrm{ext}})}.\end{array}$$(A.3)

Here ${\stackrel{\sim }{P}}_{\mathrm{ext}}=2P_{\mathrm{ext}}/(\pi G{\sigma }_{n,0}^2)$ is the dimensionless external pressure. The initial local effective sound speed, C_eff, 0, comes from

$$\begin{array}{c}\hfill C_{\mathrm{eff},0}^2=\frac{\pi }{2}G{\sigma }_{n,0}^2\frac{[3P_{\mathrm{ext}}+\frac{\pi }{2}G{\sigma }_{n,0}^2]}{{[P_{\mathrm{ext}}+\frac{\pi }{2}G{\sigma }_{n,0}^2]}^2}{c}_s^2,\end{array}$$(A.4)

and reduces to the isothermal sound speed, c_s, in the limit of very low P_ext. It can be obtained from the linearized condition of vertical hydrostatic equilibrium:

$$\begin{array}{c}\hfill {\rho }_{n,0}{c}_{\mathrm{s}}^2=\frac{\pi }{2}G{\sigma }_{n,0}^2+P_{\mathrm{ext}}=\frac{\pi }{2}G{\sigma }_{n,0}^2(1+{\stackrel{\sim }{P}}_{\mathrm{ext}}).\end{array}$$(A.5)

The initial uniform Alfvén speed V_A, 0 is related to the mass-to-magnetic-flux ratio (μ₀) via

$$\begin{array}{c}\hfill V_{A,0}^2\equiv \frac{B_{\mathrm{ref}}^2}{4\pi {\rho }_{n,0}}=2\pi G{\sigma }_{n,0}{\mathrm{\mu }}_0^{-2}{Z}_0,\end{array}$$(A.6)

where ρ_n, 0 is the initial uniform volume density. The initial uniform ambipolar diffusivity can be expressed as (see Das & Basu 2021)

$$\begin{array}{c}\hfill {\eta }_{\mathrm{AD},0}=V_{A,0}^2{\tau }_{\mathrm{ni},0}=2\pi G{\sigma }_{n,0}{\mathrm{\mu }}_0^{-2}{Z}_0{\tau }_{\mathrm{ni},0}.\end{array}$$(A.7)

Now, from the system of governing equations (see Eqs. (11)–(23) of Das & Basu 2021, which are the fundamental MHD equations obtained by integrating over the z-direction), we obtain a simplified form in x- and y- components:

$$\begin{array}{c}\hfill \frac{\partial {\sigma }_n}{\partial t}+\frac{\partial }{\partial x}({\sigma }_n{v}_{n,x})+\frac{\partial }{\partial y}({\sigma }_n{v}_{n,y})=0,\end{array}$$(A.8)

$$\begin{array}{cc}\hfill \frac{\partial }{\partial t}({\sigma }_n{v}_{n,x})& +\frac{\partial }{\partial x}({\sigma }_n{v}_{n,x}^2)+\frac{\partial }{\partial x}({\sigma }_n{v}_{n,x}{v}_{n,y})\hfill \\ \hfill & ={\sigma }_n{g}_x-C_{\mathrm{eff}}^2\frac{\partial {\sigma }_n}{\partial x}+F_{\mathrm{Mag},x}+2{\sigma }_n\mathrm{\Omega }{v}_{n,y},\hfill \end{array}$$(A.9)

$$\begin{array}{cc}\hfill \frac{\partial }{\partial t}({\sigma }_n{v}_{n,y})& +\frac{\partial }{\partial y}({\sigma }_n{v}_{n,x}{v}_{n,y})+\frac{\partial }{\partial y}({\sigma }_n{v}_{n,y}^2)\hfill \\ \hfill & ={\sigma }_n{g}_y-C_{\mathrm{eff}}^2\frac{\partial {\sigma }_n}{\partial y}+F_{\mathrm{Mag},y}-2{\sigma }_n\mathrm{\Omega }{v}_{n,x},\hfill \end{array}$$(A.10)

$$\begin{array}{cc}\hfill F_{\mathrm{Mag},x}& =\frac{B_{z,\mathrm{eq}}}{2\pi }(B_{x,Z}-Z\frac{\partial B_{z,\mathrm{eq}}}{\partial x})\hfill \\ \hfill & +\frac{1}{4\pi }\frac{\partial Z}{\partial x}[B_{x,Z}^2+B_{y,Z}^2+2B_{z,\mathrm{eq}}(B_{x,Z}\frac{\partial Z}{\partial x}+B_{y,Z}\frac{\partial Z}{\partial y})\hfill \\ \hfill & +{(B_{x,Z}\frac{\partial Z}{\partial x}+B_{y,Z}\frac{\partial Z}{\partial y})}^2],\hfill \end{array}$$(A.11)

$$\begin{array}{cc}\hfill F_{\mathrm{Mag},y}& =\frac{B_{z,\mathrm{eq}}}{2\pi }(B_{y,Z}-Z\frac{\partial B_{z,\mathrm{eq}}}{\partial y})\hfill \\ \hfill & +\frac{1}{4\pi }\frac{\partial Z}{\partial y}[B_{x,Z}^2+B_{y,Z}^2+2B_{z,\mathrm{eq}}(B_{x,Z}\frac{\partial Z}{\partial x}+B_{y,Z}\frac{\partial Z}{\partial y})\hfill \\ \hfill & +{(B_{x,Z}\frac{\partial Z}{\partial x}+B_{y,Z}\frac{\partial Z}{\partial y})}^2],\hfill \end{array}$$(A.12)

$$\begin{array}{cc}\hfill \frac{\partial B_{z,\mathrm{eq}}}{\partial t}=-\frac{\partial }{\partial x}& \left(B_{z,\mathrm{eq}}{v}_{i,x}\right)-\frac{\partial }{\partial y}\left(B_{z,\mathrm{eq}}{v}_{i,y}\right)\hfill \\ \hfill & +[\frac{\partial }{\partial x}\left({\eta }_{\mathrm{OD}}\frac{\partial B_{z,\mathrm{eq}}}{\partial x}\right)+\frac{\partial }{\partial y}\left({\eta }_{\mathrm{OD}}\frac{\partial B_{z,\mathrm{eq}}}{\partial y}\right)].\hfill \end{array}$$(A.13)

where

$$\begin{array}{c}\hfill v_{i,x}=v_{n,x}+\frac{{\tau }_{\mathrm{ni},0}}{{\sigma }_n}{\left(\frac{{\rho }_{n,0}}{{\rho }_n}\right)}^{1/2}{F}_{\mathrm{Mag},x},\end{array}$$(A.14)

$$\begin{array}{c}\hfill v_{i,y}=v_{n,y}+\frac{{\tau }_{\mathrm{ni},0}}{{\sigma }_n}{\left(\frac{{\rho }_{n,0}}{{\rho }_n}\right)}^{1/2}{F}_{\mathrm{Mag},y}.\end{array}$$(A.15)

Here, v_i, x, v_i, y, v_n, x, and v_n, y are the x- and y- components of ion and neutral velocities. The planar sheet is rotating with an angular velocity Ω about the z-axis, so that $\mathbf{\Omega }=\mathrm{\Omega }\widehat{z}$. The magnetic field and rotation axis are perpendicular to the sheet. Here, η_OD,0 is the ohmic diffusivity that is considered as a measure of ohmic dissipation. Starting with a static uniform background, any physical quantity of the thin-sheet equations can be expanded by writing it via

$$\begin{array}{c}\hfill f(x,y,t)=f_0+\delta f_a{e}^{i(k_{x}x+k_{y}y-\omega t)},\end{array}$$(A.16)

where f₀ is the unperturbed background state, δf_a is the amplitude of the perturbation, k_x, k_y, and k are the x-, y-, and z- wavenumbers, respectively, and ω is the complex angular frequency. For assumed small-amplitude perturbations such that |δf_a| ≪ f₀, and retaining the linearized form of the perturbed quantities from Eqs. (A.8)–(A.10) and (A.13), the following equations are obtained:

$$\begin{array}{c}\hfill \omega \delta \sigma {\prime }_n=k_x{c}_{\mathrm{s}}\delta v{\prime }_{n,x}+k_y{c}_{\mathrm{s}}\delta v{\prime }_{n,\mathrm{y}},\end{array}$$(A.17)

$$\begin{array}{cc}\hfill \omega c_{\mathrm{s}}\delta v{\prime }_{n,x}=& \frac{k_x}{k}[C_{\mathrm{eff},0}^{2}k-2\pi G{\sigma }_{n,0}]\delta \sigma {\prime }_n+i2\mathrm{\Omega }{c}_{\mathrm{s}}\delta v{\prime }_{n,y}\hfill \\ \hfill & +\frac{k_x}{k}[2\pi G{\sigma }_{n,0}{\mathrm{\mu }}_0^{-1}+k{V}_{A,0}^2{\mathrm{\mu }}_0]\delta B{\prime }_{z,\mathrm{eq}},\hfill \end{array}$$(A.18)

$$\begin{array}{cc}\hfill \omega c_{\mathrm{s}}\delta v{\prime }_{n,y}=& \frac{k_y}{k}[C_{\mathrm{eff},0}^{2}k-2\pi G{\sigma }_{n,0}]\delta \sigma {\prime }_n-i2\mathrm{\Omega }{c}_{\mathrm{s}}\delta v{\prime }_{n,x}\hfill \\ \hfill & +\frac{k_y}{k}[2\pi G{\sigma }_{n,0}{\mathrm{\mu }}_0^{-1}+k{V}_{A,0}^2{\mathrm{\mu }}_0]\delta B{\prime }_{z,\mathrm{eq}},\hfill \end{array}$$(A.19)

$$\begin{array}{cc}\hfill \omega & \delta B{\prime }_{z,\mathrm{eq}}=\frac{k_x}{{\mathrm{\mu }}_0}{c}_{\mathrm{s}}\delta v{\prime }_{n,x}+\frac{k_y}{{\mathrm{\mu }}_0}{c}_{\mathrm{s}}\delta v{\prime }_{n,y}\hfill \\ \hfill & -i[{\eta }_{\mathrm{OD},0}{k}^2+{\tau }_{\mathrm{ni},0}(2\pi G{\sigma }_{n,0}{\mathrm{\mu }}_0^{-2}k+k^2{V}_{A,0}^2)]\delta B{\prime }_{z,\mathrm{eq}},\hfill \end{array}$$(A.20)

where the perturbed eigenfunctions δσ_n, δv_n, x (and δv_n, y), δB_z, eq are normalized by σ_n, 0, c_s, and B₀ (= 2πG^1/2σ_n, 0), respectively, such that $\delta {\sigma }_n^{\prime }=\delta {\sigma }_n/{\sigma }_{n,0}$, $\delta v_{n,x}^{\prime }=\delta v_{n,x}/c_{\text{s}}$, $\delta v_{n,y}^{\prime }=\delta v_{n,\text{y}}/c_{\text{s}}$, and $\delta {B\prime }_{z,\mathrm{eq}}=\delta B_{z,\mathrm{eq}}/B_0$. Now, finding the determinant from the above set of equations one obtains the full dispersion relation

$$\begin{array}{cc}\hfill (\omega +i[\theta +\gamma ])({\omega }^2& -C_{\mathrm{eff},0}^2{k}^2+2\pi G{\sigma }_{n,0}k-4{\mathrm{\Omega }}^2)\hfill \\ \hfill & =\omega [2\pi G{\sigma }_{n,0}k{\mathrm{\mu }}_0^{-2}+k^2{V}_{A,0}^2].\hfill \end{array}$$(A.21)

This is derived in Das & Basu (2021; see Eq. (47)), where

$$\begin{array}{c}\hfill \gamma ={\eta }_{\mathrm{OD},0}{k}^2,\end{array}$$(A.22)

which is same as Eq. (45) of Das & Basu (2021), and θ is described earlier (see Eq. (5) and Eq. (46) of Das & Basu (2021)). Finally, the dispersion relation (see Eq. (4)) is the same as the Eq. (A.21) when setting Ω = 0 and η_OD,0 = 0. In the above, we discuss all the equations in detail for completeness and clarity of our model. See also Eqs. (32)a–d of Ciolek & Basu (2006) for the dimensionless representation of the Eqs. (A.17)–(A.20) for the model with Ω = 0 and η_OD,0 = 0; also see Eqs. (10)–(13) of Bailey & Basu (2012) for the dimensional representation of Eq. (A.17)–(A.20) for the same model.

Based on these parameters, typical values of the units used and other derived quantities are

$$\begin{array}{c}\hfill {\sigma }_{n,0}=\frac{3.63\times {10}^{-3}}{{(1+P_{\mathrm{ext}})}^{\frac{1}{2}}}{\left(\frac{n_{n,0}}{{10}^3{\mathrm{cm}}^{-3}}\right)}^{\frac{1}{2}}{\left(\frac{T}{10\mathrm{K}}\right)}^{\frac{1}{2}}\mathrm{g}{\mathrm{cm}}^{-2},\end{array}$$(A.23)

$$\begin{array}{c}\hfill L_0=7.48\times {10}^{-2}{\left(\frac{T}{10\mathrm{K}}\right)}^{\frac{1}{2}}{\left(\frac{{10}^3{\mathrm{cm}}^{-3}}{n_{n,0}}\right)}^{\frac{1}{2}}{(1+P_{\mathrm{ext}})}^{\frac{1}{2}}\mathrm{pc},\end{array}$$(A.24)

$$\begin{array}{c}\hfill t_0=3.98\times {10}^5{\left(\frac{{10}^3{\mathrm{cm}}^{-3}}{n_{n,0}}\right)}^{\frac{1}{2}}{(1+P_{\mathrm{ext}})}^{\frac{1}{2}}\mathrm{yr},\end{array}$$(A.25)

$$\begin{array}{c}\hfill c_{\mathrm{s}}=0.188{\left(\frac{T}{10\mathrm{K}}\right)}^{\frac{1}{2}}\mathrm{km}{\mathrm{s}}^{-1},\end{array}$$(A.26)

$$\begin{array}{c}\hfill M_0=9.76\times {10}^{-2}{\left(\frac{T}{10\mathrm{K}}\right)}^{3/2}{\left(\frac{{10}^3{\mathrm{cm}}^{-3}}{n_{n,0}}\right)}^{1/2}{(1+P_{\mathrm{ext}})}^{1/2}{M}_{\odot },\end{array}$$(A.27)

$$\begin{array}{c}\hfill {\tau }_{\mathrm{ni},0}=\frac{3.74\times {10}^4}{(1+P_{\mathrm{ext}})}\left(\frac{T}{10\mathrm{K}}\right){\left(\frac{0.01\mathrm{g}{\mathrm{cm}}^{-2}}{{\sigma }_{n,0}}\right)}^2\left(\frac{{10}^{-7}}{{\chi }_{i,0}}\right)\mathrm{yr},\end{array}$$(A.28)

$$\begin{array}{c}\hfill B_{\mathrm{ref}}=\frac{5.89\times {10}^{-6}}{{\mathrm{\mu }}_0}{\left(\frac{n_{n,0}}{{10}^3{\mathrm{cm}}^{-3}}\right)}^{1/2}{\left(\frac{T}{10\mathrm{K}}\right)}^{1/2}{(1+P_{\mathrm{ext}})}^{-1/2}\mathrm{G}.\end{array}$$(A.29)

Fig. A.1. Top: normalized shortest growth time of gravitationally unstable mode (${\tau }_{g,m}^{\prime }={\tau }_{g,m}/t_0$) as a function of normalized mass-to-flux ratio (μ0). Bottom: normalized preferred length scale of the most unstable mode (${\lambda }_{g,m}^{\prime }={\lambda }_{g,m}/L_0$) as a function of the normalized mass-to-flux ratio (μ0). Each panel shows curves for models (as shown in Eq. (4)) with normalized neutral–ion collision time ${\stackrel{\sim }{\tau }}_{\mathrm{ni},0}=0$ (black), 0.001 (red), 0.04 (blue), and 0.2 (green). Here, ${\stackrel{\sim }{\tau }}_{\mathrm{ni},0}=0$ represents the flux-frozen case as a reference. The Toomre-Q rotation parameter and the dimensionless ohmic diffusivity ${\stackrel{\sim }{\eta }}_{\text{OD,0}}$ are set to be zero. As the degree of ambipolar diffusion (i.e., ${\stackrel{\sim }{\tau }}_{\text{ni},0}$) increases, the growth timescale and length scale become shorter (see Ciolek & Basu 2006; Bailey & Basu 2012; Das & Basu 2021) for further details.

Appendix B: Tables

All Tables

All Figures

Fig. 1. Growth timescale τg, m (in the units of dynamical or free-fall time tff) and the preferred fragmentation mass Mg, m (in the units of preferred thermal fragmentation mass MT, m) of the most unstable mode as a function of the normalized mass-to-flux ratio (μ0). The model utilizes the standard ionization fraction corresponding to a normalized neutral–ion collision time ${\tau }_{\mathrm{ni},0}/t_0={\stackrel{\sim }{\tau }}_{\text{ni},0}=0.2$. The shaded and unshaded zones represent the subcritical (μ0 < 1) and supercritical (μ0 > 1) regimes, respectively. The present paper focuses on the (mildly) supercritical regime.

In the text

Fig. 2. Estimated lifetime as a function of number density (nn). The green filled squares are obtained from our model for a normalized neutral–ion collision time ${\stackrel{\sim }{\tau }}_{\text{ni},0}=0.2$. The blue filled triangles show the corresponding data points for the population of 446 candidate prestellar cores identified with Herschel in the Aquila cloud (Könyves et al. 2015). Earlier data from Ward-Thompson et al. (2007) are shown by the red circles for better comparison. The black dashed line shows the dynamical time (Z0/cs), i.e., the free-fall time, of our model.

In the text

Fig. 3. Number of enclosed cores (log[NumCORE]) as a function of Jeans number (log[NJ, CLUMP]) of clumps. The blue filled circles represent the observed number of enclosed cores as found by Pokhrel et al. (2018). The green filled crosses are obtained from our model for a normalized neutral–ion collision time ${\stackrel{\sim }{\tau }}_{\text{ni},0}=0.2$. The black dashed line corresponds to an efficiency of unity.

In the text

Fig. 4. Magnetic field log(Bref/μG) vs. number density $log(n_n/{\mathrm{cm}}^{-3})$. The filled squares and crosses represent the number density region corresponding to the fitting of the lifetime of prestellar cores and NumCORE as shown in Figs. 2 and 3 (see Tables B.1 and B.2). The dashed line is the least-squares fit. The vertical error bars correspond to magnetic field variations that yield a total factor of two (greater or lesser by a factor $\sqrt{2}$) values of the growth time (lifetime).

In the text

Fig. A.1. Top: normalized shortest growth time of gravitationally unstable mode (${\tau }_{g,m}^{\prime }={\tau }_{g,m}/t_0$) as a function of normalized mass-to-flux ratio (μ0). Bottom: normalized preferred length scale of the most unstable mode (${\lambda }_{g,m}^{\prime }={\lambda }_{g,m}/L_0$) as a function of the normalized mass-to-flux ratio (μ0). Each panel shows curves for models (as shown in Eq. (4)) with normalized neutral–ion collision time ${\stackrel{\sim }{\tau }}_{\mathrm{ni},0}=0$ (black), 0.001 (red), 0.04 (blue), and 0.2 (green). Here, ${\stackrel{\sim }{\tau }}_{\mathrm{ni},0}=0$ represents the flux-frozen case as a reference. The Toomre-Q rotation parameter and the dimensionless ohmic diffusivity ${\stackrel{\sim }{\eta }}_{\text{OD,0}}$ are set to be zero. As the degree of ambipolar diffusion (i.e., ${\stackrel{\sim }{\tau }}_{\text{ni},0}$) increases, the growth timescale and length scale become shorter (see Ciolek & Basu 2006; Bailey & Basu 2012; Das & Basu 2021) for further details.

In the text