A Low-mass Cold and Quiescent Core Population in a Massive Star Protocluster

Stars form in the high-density parts of molecular clouds, i.e., dense cores. A dense core prior to the protostellar phase is called a starless core, and a gravitationally bound/unstable starless core is referred to as a pre-stellar core (Ward-Thompson et al. 1994). The pre-stellar phase is considered to be the starting point in the star formation process (Bergin & Tafalla 2007; Caselli & Ceccarelli 2012). Pre-stellar cores can form single or multiple stellar systems under the combined effect of gravity, magnetic fields, and turbulence (Pineda et al. 2015). Understanding the properties of pre-stellar cores is critical for characterizing the initial conditions of cluster formation. The pre-stellar cores in low-mass star-forming regions near nearby molecular clouds, e.g., Perseus, Ophiuchus, Chamaeleon, Serpens, and Taurus (Lee et al. 1999, 2001; Enoch et al. 2008; Schnee et al. 2010; Dunham et al. 2016; Kirk et al. 2017; Tokuda et al. 2020) have been intensively studied. However, the studies of pre-stellar cores in massive cluster-forming regions are still limited by the low number of statistics available (Nony et al. 2018; Li et al. 2019; Sanhueza et al. 2019; Sahu et al. 2021).

The massive infrared dark cloud NGC 6334S (also known as IRDC G350.56+0.44) is located at the southwestern end of the NGC 6334 molecular cloud complex, which is a nearby (1.3 kpc; Chibueze et al. 2014) young and massive "mini-starburst" star-forming region (Willis et al. 2013). With a mass of 1.3 × 10³ M_⊙ (Li et al. 2020b), comparable to the clumps with embedded massive protostars and protoclusters in the complex, NGC 6334S has the potential to form a cluster encompassing both low- (<2 M_⊙) and high-mass (>8 M_⊙) stars. NGC 6334S provides an ideal laboratory to search for and study the early stages (e.g., pre-stellar phase) of star formation in the high-mass regime.

In order to investigate massive star and cluster formation, we performed Atacama Large Millimeter/submillimeter Array (ALMA) and Karl G. Jansky Very Large Array (JVLA) observations of NGC 6334S. A study of 49 continuum dense cores (hereafter continuum cores) revealed by the 3 mm wavelength continuum image was presented in Li et al. (2020b), in which we reveal that the nonthermal motions are predominantly subsonic and transonic in both clump scale and continuum core scale. Here we use molecular lines to study the starless/pre-stellar cores.

2.1. ALMA Observations

2.2. JVLA Observations

3.1. Cold and Quiescent Cores

NH₂D (1₁₁–1₀₁, critical density n_cr ∼ 10⁵ cm⁻³) is a good tracer of cold and dense molecular gas (Crapsi et al. 2007; Sanhueza et al. 2013), which can survive in the gas phase in the dense interior region of pre-stellar cores (Harju et al. 2017). In our ALMA data, the NH₂D (1₁₁–1₀₁) line emission is in general well correlated with the 3 mm continuum emission, but there are some exceptions: there are 17 bright compact structures in the NH₂D emission that are associated with weak (3σ–11σ) or no continuum emission (<3σ) at all (Figures 1 and 2), and with no young stellar object (YSO) counterparts (e.g., Class I/II; Willis et al. 2013). There is no reliable Class 0 catalog in this region. The details on the identification of NH₂D compact structures are summarized in Appendix A. We refer to these NH₂D compact structures as NH₂D cores, naming those in M1, M2, M3 ... in order of descending NH₂D velocity-integrated intensity. They have diameters ranging from 0.018 to 0.04 pc (see Table 1).

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Table 1. Physical Parameters of the NH₂D Cores

Core ID	R.A.	Decl.	$\mathrm{Maj}\times \mathrm{Min}$	P.A.	R	σobs	vLSR	${ \mathcal M }$	${W}_{{\mathrm{NH}}_{2}{\rm{D}}}^{\mathrm{peak}}$	${W}_{{\mathrm{NH}}_{2}{\rm{D}}}$	${S}_{\nu }^{\mathrm{peak}}$	Sν	Mgas	Mvir	αvir	${N}_{{\mathrm{NH}}_{2}{\rm{D}}}$	${N}_{{{\rm{H}}}_{2}}$	X(NH2D)	${D}_{{\mathrm{NH}}_{3}}$
	(J2000)	(J2000)	(''×'')	(deg)	(pc)	(km s−1)	(km s−1)		(Jy/beam km/s)	(Jy/beam km/s)	(mJy beam−1)	(mJy)	(M⊙)	(M⊙)		(cm−2)	(cm−2)
(1)	(2)	(3)	(4)	(5)	(6)	(7)	(8)	(9)	(10)	(11)	(12)	(13)	(14)	(15)	(16)	(17)	(18)	(19)	(20)
M1	17:19:09.38	−36:03:37.96	7.1 × 2.8	107.6	0.028	0.11	−3.59	<0.68	0.69	140.95	0.174	0.294	0.43	<0.68	<1.57	2.2E+15	1.8E+22	1.1E-07	⋯
M2	17:19:04.28	−36:06:57.17	9.3 × 4.3	126.6	0.040	0.30	−1.36	2.41	0.40	136.13	0.186	0.412	0.61	2.54	4.18	6.1E+14	1.6E+22	3.2E-08	>0.39
M3	17:19:05.61	−36:11:01.29	7.5 × 3.9	75.4	0.034	0.20	−2.95	1.38	0.41	104.66	0.292	0.462	0.68	1.25	1.84	7.0E+14	2.4E+22	2.8E-08	⋯
M4	17:19:11.12	−36:05:41.37	9.4 × 3.2	79.5	0.034	0.17	−5.07	1.10	0.31	82.90	0.214	0.589	0.87	1.00	1.15	4.8E+14	2.8E+22	1.8E-08	0.25
M5	17:19:07.39	−36:10:24.52	7.5 × 2.1	72.1	0.025	0.15	−2.47	0.84	0.44	82.36	0.180	0.258	0.38	0.63	1.65	9.3E+14	1.8E+22	4.5E-08	⋯
M6	17:19:12.04	−36:06:55.76	7.3 × 3.9	57.3	0.034	0.18	−2.80	1.13	0.28	69.97	0.230	0.458	0.68	1.08	1.60	4.0E+14	2.4E+22	1.7E-08	>0.25
M7	17:19:09.91	−36:09:05.62	6.6 × 5.2	21.1	0.037	0.23	−4.17	1.69	0.23	65.62	0.206	0.447	0.66	1.83	2.77	1.8E+14	2.2E+22	8.0E-09	⋯
M8	17:19:10.30	−36:09:09.49	5.8 × 4.4	122.3	0.032	0.19	−4.22	1.32	0.24	56.76	0.237	0.450	0.66	1.27	1.91	3.0E+14	2.5E+22	1.2E-08	⋯
M9	17:19:11.40	−36:08:44.12	4.9 × 3.9	78.0	0.028	0.12	−1.68	0.34	0.24	45.38	0.158	0.201	0.30	0.73	2.47	2.7E+14	1.4E+22	1.4E-08	0.22
M10	17:19:05.12	−36:06:45.00	4.2 × 3.7	9.9	0.025	0.15	−3.47	0.78	0.26	45.19	0.310	0.366	0.54	0.79	1.46	6.9E+14	2.8E+22	2.5E-08	>0.11
M11	17:18:58.66	−36:07:04.76	4.1 × 2.2	48.5	0.019	0.11	−2.01	<0.73	0.33	44.30	⋯	⋯	<0.13	<0.51	3.82	7.5E+14	<1.5E+22	>7.1E-08	>0.32
M12	17:19:06.86	−36:06:22.65	6.8 × 3.6	124.3	0.031	0.13	−3.51	0.39	0.18	41.11	0.244	0.538	0.79	0.74	0.93	2.5E+14	3.0E+22	8.9E-09	>0.16
M13	17:19:07.39	−36:05:55.13	6.9 × 3.7	116.2	0.032	0.10	−5.36	<0.57	0.15	36.21	⋯	⋯	<0.13	<0.80	6.03	1.3E+14	<1.5E+22	>7.8E-09	>0.14
M14	17:19:04.55	−36:05:31.47	3.7 × 2.2	29.7	0.018	0.16	−4.20	0.93	0.26	33.59	⋯	⋯	<0.13	0.54	>4.04	2.8E+14	<1.5E+22	>1.2E-08	⋯
M15	17:19:08.88	−36:08:01.55	5.3 × 2.7	78.2	0.024	0.11	−4.26	<0.70	0.19	29.68	⋯	⋯	<0.13	<0.62	4.64	2.0E+14	<1.5E+22	>2.0E-08	>0.17
M16	17:19:07.21	−36:08:16.60	6.4 × 2.3	125.7	0.024	0.15	−4.18	0.81	0.16	29.59	0.160	0.168	0.25	0.60	2.44	1.6E+14	1.2E+22	1.1E-08	>0.14
M17	17:19:07.78	−36:10:07.54	4.8 × 1.8	57.5	0.019	0.13	−2.52	0.49	0.13	18.00	0.158	0.170	0.25	0.43	1.70	1.6E+14	1.6E+22	7.2E-09	⋯
Mean					0.028	0.16		0.96	0.29	62.50	0.211	0.370	0.45	0.94	2.60	5.1E+14	2.0E+22	2.6E-08	0.22
Median					0.028	0.15		0.81	0.26	45.38	0.206	0.412	0.43	0.74	1.91	3.0E+14	1.8E+22	1.7E-08	0.20
Minimum					0.018	0.10		0.34	0.13	18.00	0.158	0.168	0.13	0.43	0.93	1.3E+14	1.2E+22	7.2E-09	0.11
Maximum					0.040	0.30		2.41	0.69	140.95	0.310	0.589	0.87	2.54	6.03	2.2E+15	3.0E+22	1.1E-07	0.39

Notes. For M1, M11, M13, and M15, ${\sigma }_{\mathrm{nt},{\mathrm{NH}}_{2}{\rm{D}}}$ = σ_obs was used as the ${\sigma }_{\mathrm{obs},\mathrm{int}}\leqslant {\sigma }_{\mathrm{nt},{\mathrm{NH}}_{2}{\rm{D}}}(10{\rm{K}})$. (4)–(5) Beam-deconvolved size. (6) Beam-deconvolved effective radius. (7)–(8) σ_obs and v_LSR are derived by the core-averaged spectrum. (9) Mach number ${ \mathcal M }$ = $\sqrt{3}{\sigma }_{\mathrm{nt},{\mathrm{NH}}_{2}{\rm{D}}}/{c}_{{\rm{s}}}$. (10) The peak value of NH₂D velocity-integrated image. (11) The integrated flux density of the NH₂D line emission. (12)–(13) The peak and integrated intensity of continuum emission. (14) Gas mass estimated by continuum emission. (15) Virial mass. (16) Virial parameter. (17) Mean NH₂D column density. (18) Mean H₂ column density. (19) Mean NH₂D abundance fraction. (20) The NH₃ deuterium fractionation ${D}_{{\mathrm{NH}}_{3}}$.

Download table as:  ASCIITypeset image

Among the NH₂D cores, M1 is the most prominent pre-stellar candidate with high S/N and a relatively isolated environment that eliminates the contamination from the outflows driven by the continuum cores. Therefore, the following analysis and discussion will use M1 as a showcase, while the physical parameters of the remaining NH₂D cores will only be included in statistical analyses.

Using the HCO⁺, CS, HCN, SO, HNCO, and H₂O maser lines (see Figure 1), we have searched for signatures of protostellar activity, such as bipolar/monopolar/multipolar outflow activity and H₂O maser emission. All the NH₂D cores show no signs of protostellar activity. Although we currently cannot rule out the existence of unresolved or weak outflows due to the coarse spectral resolution and the lack of better low-mass outflow tracers such as CO, these NH₂D cores appear to be starless with the available evidence. In addition, the NH₂D cores do not show any appreciable (<3σ) continuum emission in the infrared wavelengths from 3.6 up to 70 μm, and therefore they are completely dark in these infrared images. The upper limit on the internal luminosity (L_int) of NH₂D cores is estimated to be ∼1.26 L_⊙ using Herschel/PACS 70 μm data (Tigé et al. 2017; see Appendix C for detailed derivation of internal luminosity). These results suggest that these NH₂D cores are cold and quiescent in terms of star-forming activity, though some undetected very faint embedded low-mass protostars cannot be fully ruled out (Dunham et al. 2008; Evans et al. 2009; Offner & McKee 2011; Fischer et al. 2017). An example of the comparison between NH₂D core (M1) and continuum core (#2) is shown in Figure 1, which includes all of the observed spectral windows. Continuum core #2 is a protostellar core as evidenced by the clear molecular outflows in the CS, HCN, and HCO⁺ line emission, and it is one of the most chemically rich among the continuum cores. Nonetheless, only a few lines are detected toward continuum core #2. This indicates that the continuum cores #2 is at a stage prior to the hot cores/corinos. We note that M1 shows even less chemical complexity than continuum core #2, with fewer lines, molecular species, and narrower line profiles. This suggests that the evolutionary stages of the NH₂D cores are likely earlier than the continuum cores. On the other hand, we cannot completely rule out the possibility that some of continuum cores are relatively evolved pre-stellar cores, because of the lack of high sensitivity and high spectral resolution outflow tracers (e.g., CO, SiO) to distinguish the evolutionary stages (Li et al. 2020a). In this study, we will focus only on the NH₂D cores.

3.2. Small Velocity Dispersions

The NH₂D line cube is fitted pixel by pixel, following the fitting processes described in Li et al. (2020b). In the hyperfine fits, the uncertainties of the centroid velocity and linewidths are small, while those of the excitation temperature and optical depth are relatively large (Figure 3). Figure 4 shows the observed velocity dispersions (σ_obs) distribution for each NH₂D core. All the NH₂D cores show small velocity dispersion (〈σ_obs〉 ∼ 0.16 km s⁻¹) and the majority of them are smaller than the sound speed, c_s(10 K) = 0.19 km s⁻¹. The observed velocity dispersions may be regarded as upper limits, given the limited spectral resolution, the blending of velocity components along the line of sight, and opacity-broadening effects. For the NH₂D cores, the linewidths of detected lines are smaller than those of the majority of continuum cores (see Figures 1 and 4).

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Among all the NH₂D cores, M1 shows the strongest NH₂D emission and the smallest velocity dispersion. For M1, pixel by pixel fit to the NH₂D line cube yields a small mean velocity dispersion 〈σ_obs〉 ≈ 0.107 km s⁻¹ (Figure 4), corresponding to an intrinsic observed velocity dispersion of σ_obs,int = 0.06 km s⁻¹ after removing the smearing effect due to the channel width using ${\sigma }_{\mathrm{obs},\mathrm{int}}=\sqrt{{\sigma }_{\mathrm{obs}}^{2}-{\sigma }_{\mathrm{ch}}^{2}}$, where ${\sigma }_{\mathrm{ch}}={\bigtriangleup }_{\mathrm{ch}}/{\left(2\sqrt{2\mathrm{ln}2}\right)}^{2}$ = 0.089 km s⁻¹ is the channel width. For the extremely narrow line emission, higher spectral resolution data is required to accurately characterize the actual intrinsic velocity dispersion. Considering the 〈σ_obs〉 is 1.2 times larger than the σ_ch, the aforementioned derived σ_obs,int can be regarded as a reasonable approximation for the following analysis (see Appendix C in Li et al. 2020b). Assuming that the observed linewidths of the NH₂D line are caused only by thermal broadening, the gas temperature is found to be about 7 K, ${\sigma }_{\mathrm{th},{\mathrm{NH}}_{2}{\rm{D}}}$(7 K) ²² =σ_obs,int = 0.06 km s⁻¹. This indicates that the gas kinematics is dominated by thermal motions with a kinetic temperature that is 7 K at most. In addition to M1, there are four NH₂D cores (M9, M12, M13, and M15) that also show extremely small velocity dispersions comparable to the NH₂D thermal velocity dispersion in the temperature range of 7–10 K. The velocity dispersions of remaining cores are relatively higher than the NH₂D thermal velocity dispersion at a temperature of 10 K (${\sigma }_{\mathrm{th},{\mathrm{NH}}_{2}{\rm{D}}}\,\sim$ 0.07 km s⁻¹), while they are much smaller than the sound speed (Figure 4). This indicates that the gas kinematics are dominated by thermal motions toward these NH₂D cores.

Figure 3 shows the annularly averaged observed velocity dispersion as a function of the radial distance (R_dist) from the center of M1. The observed velocity dispersion appears to increase with increasing R_dist (see Figure 3), which prevails in pre-stellar and starless cores (Crapsi et al. 2005; Auddy et al. 2019). The increase of σ_obs along R_dist is also detected in the other nine cores (M4, M6, M7, M8, M9, M10, M13, M15, M17), although the annularly averaged σ_obs starts to suffer from low S/N toward the outer edges of the cores. This increasing trend of σ_obs along R_dist indicates that the turbulence is significantly dissipated toward the center of the core. The results indicate that the turbulence dissipation is commonly seen toward these NH₂D cores. Alternatively, in the larger-scale collapse scenario, if the pre-stellar cores undergo accretion, the infall speed decreases towards the center, and thus the line width will be smaller for more central regions (Gómez et al. 2021). Although no clear infall signatures are detected in M1 using coarser spectral resolution (3.0–3.3 km s⁻¹) optically thick lines (e.g., HCO⁺ 1–0, HCN 1–0 and CS 2–1), this possibility cannot be ruled out and should be explored using higher spectral resolution data and other tracers (e.g., HCN 5–4/4–3/3–2 and HCO⁺ 5–4/3–2; Chira et al. 2014).

Overall, the observed small velocity dispersions and their spatial distributions resemble those observed in the well-known pre-stellar cores in Ophiuchus/H-MM1 and L 1544 (Crapsi et al. 2007; Harju et al. 2017; Punanova et al. 2018). This suggests that the identified NH₂D cores may have not yet been affected by the star formation activities.

3.3. Gas Masses and Dynamical States

3.4. NH₂D Abundances and NH₃ Deuterium Fractionation

As shown in Figure 4, the derived NH₂D column densities and fractional abundances are higher than the values of Ophiuchus/H-MM1 (1.1 × 10¹⁴ cm⁻² and 2 × 10⁻⁹; Harju et al. 2017) and L1544 (1.8 × 10¹⁴ cm⁻² and 2 × 10⁻¹⁰; Crapsi et al. 2007). The high abundances of NH₂D could benefit from the significant depletion of CO, which occurs faster at cold and dense environments (e.g., starless or pre-stellar cores; the detail chemical processes related to the NH₂D formation can be found in Sipilä et al. 2015; Harju et al. 2017).

The small velocity dispersion, small Mach number (Table 1), and positive σ_obs–R_dist relation indicate that the turbulence is likely dissipated toward these NH₂D cores. In addition, it implies that star formation has not yet started; outflows motions would widen the lines as can be seen in the continuum cores (Figure 1). Therefore, these results, in conjunction with the high NH₃ deuterium fractionation and the lack of embedded YSOs and 70 μm emission, indicate that the identified NH₂D cores are still in a starless phase. Nine out of 17 starless cores are gravitationally bound with α_vir < 2; therefore, they are identified as excellent pre-stellar core candidates. The dust temperature may drop to ∼5.5 K toward the core center (e.g., L1544; Crapsi et al. 2007). The derived gas mass can increase by a factor of 2.25 if a temperature of 5.5 K is assumed, so that the resulting virial parameter will decrease by a factor of 2.25. The levels of the missing flux are unknown due to the lack of single dish data. However, the estimated masses and radius of cores may not be significantly affected by the missing flux, because the observations are carried out with the most compact 12 m array configuration (C40-1) and the MRS is ∼25'', which is much larger than the cores sizes. Therefore, the estimated virial parameters are not seriously affected by the missing flux. In addition, the ambient pressure from the parental clumps and filaments could also provide additional confinement to make these NH₂D cores bound (e.g., Chen et al. 2019; Li et al. 2020b). Therefore, we cannot completely rule out the possibility that the other eight starless cores could also be pre-stellar core candidates.

The estimated masses of NH₂D cores range from 0.13 to 0.87 M_⊙, suggesting that they are in the low-mass core regime. The majority (15) of NH₂D cores are associated with filamentary structures revealed by the H¹³CO⁺ line emission (S. Li et al. 2021, in preparation). Therefore, the NH₂D cores have the potential to accrete a significant amount of material from the parental clump via filamentary accretion. The masses of continuum cores are found to range from low- to high-mass regime (0.17–14.03 M_⊙; Li et al. 2020b). These results suggest that cores with a large range of masses are simultaneously forming in this cluster-forming cloud. The massive cores are found toward the center of cloud (see Li et al. 2020b), while the low-mass cores are wide spread throughout the region. The evolutionary stages vary significantly over the sample of NH₂D cores, continuum cores, and YSOs (Class I and II). This indicates that there is a temporally continuous formation of low- to high-mass cores.

In summary, NGC 6334S is forming a stellar protocluster, its embedded cores show a wide diversity in masses (0.13–14.03 M_⊙) and evolutionary stages (pre-stellar–Class II phases), and the embedded cores are expected to grow in mass by gas accretion from parental clump via filaments. This seems to be consistent with the competitive accretion massive star formation scenario (Bonnell & Bate 2006), in which the dense cores located near the center of the gravitational potential continue accreting material to form massive stars.

We present ALMA and JVLA high spatial resolution observations toward the massive IRDC NGC 6334S. We have identified 17 low-mass starless core candidates that show bright NH₂D emission, but no YSO counterparts and no signs of protostellar activity, although we cannot completely rule out the possibility of very weak (and thus undetected) outflows or very faint embedded low-mass protostars. These candidates show almost-thermal velocity dispersions (down to σ_tot ∼ 0.06 km s⁻¹), high NH₂D abundances (up to ∼10⁻⁷), high NH₃ deuterium fractionations (up to >0.39), and are completely dark in the infrared wavelengths from 3.6 up to 70 μm. In addition, the σ_obs appear to decrease toward the center of cores in 10 of them. These results suggest that turbulence has significantly dissipated and the NH₂D abundance has significantly enhanced toward these cores. The gas kinematics resemble the well-known pre-stellar cores, e.g., Ophiuchus/H-MM1 and L 1544, but with one to two orders of magnitude greater NH₂D abundance (Crapsi et al. 2007; Harju et al. 2017). Nine out of the 17 NH₂D cores are gravitationally bound, and therefore are identified as pre-stellar core candidates. Therefore, the NH₂D line could be a powerful tracer to reveal the starless and pre-stellar cores that do not show significant dust continuum emission. This is the first detection of a cluster of low-mass starless and pre-stellar core candidates in a massive star cluster-forming cloud. Low- to high-mass cores are simultaneously forming in this cluster-forming cloud.

We thank the anonymous referee for constructive comments that helped improve this Letter. This work was partially supported by National Natural Science Foundation of China (NSFC, grant Nos. 11988101, 11911530226). C.W.L. is supported by the Basic Science Research Program through the National Research Foundation of Korea (NRF) funded by the Ministry of Education, Science and Technology (NRF-2019R1A2C1010851). H.B. acknowledges support from the European Research Council under the European Community's Horizon 2020 framework program (2014-2020) via the ERC Consolidator Grant "From Cloud to Star Formation (CSF)" (project number 648505). H.B. also acknowledges funding from the Deutsche Forschungsgemeinschaft (DFG) via the Collaborative Research Center (SFB 881) "The Milky Way System" (subproject B1). I.J.-S. has received partial support from the Spanish State Research Agency (AEI; project number PID2019-105552RB-C41). K.Q. is partially supported by National Key R&D Program of China No. 2017YFA0402600, and acknowledges the National Natural Science Foundation of China (NSFC) grant U1731237. A.P. acknowledges financial support from the UNAM-PAPIIT IN111421 grant, the Sistema Nacional de Investigadores of CONACyT, and from the CONACyT project number 86372 of the "Ciencia de Frontera 2019" program, entitled "Citlalcóatl: A multiscale study at the new frontier of the formation and early evolution of stars and planetary systems," México. J.M.G. acknowledges the support of the grant AYA2017-84390-C2-R (AEI/FEDER, EU). K.W. acknowledges support by the National Key Research and Development Program of China (2017YFA0402702, 2019YFA0405100), the National Science Foundation of China (11973013, 11721303), NSFC grant 11629302, and the starting grant at the Kavli Institute for Astronomy and Astrophysics, Peking University (7101502287). H.B.L. is supported by the Ministry of Science and Technology (MoST) of Taiwan (grant Nos. 108-2112-M-001-002-MY3). This Letter makes use of the following ALMA data: ADS/JAO.ALMA#2016.1.00951.S. ALMA is a partnership of ESO (representing its member states), NSF (USA) and NINS (Japan), together with NRC (Canada), MOST and ASIAA (Taiwan), and KASI (Republic of Korea), in cooperation with the Republic of Chile. The Joint ALMA Observatory is operated by ESO, AUI/NRAO and NAOJ.

Facilities: ALMA - Atacama Large Millimeter Array, JVLA - , Herschel. -

Software: CASA (McMullin et al. 2007), APLpy (Robitaille & Bressert 2012), Astropy (Astropy Collaborationet al. 2013), Matplotlib (Hunter 2007), PySpecKit (Ginsburg & Mirocha 2011).

The astrodendro ²³ algorithm was used to preselect the compact NH₂D emission structures (i.e., the leaves in the terminology of astrodendro; Rosolowsky et al. 2008) from the NH₂D velocity-integrated intensity image (integrated between −11.5 and 5.6 km s⁻¹). The following parameters are used in the computation of astrodendro: the minimum pixel value is 3σ noise level (1σ ∼ 0.02 Jy beam⁻¹ km s⁻¹ for ${W}_{{\mathrm{NH}}_{2}{\rm{D}}}$); the peak value is ≥6σ; the minimum difference in the peak intensity between neighboring compact structures is 1σ; and the minimum number of pixels required for a structure to be considered an independent entity is 40 (approximately the synthesized beam area). From the identified leaves, we further search for the compact structures that are not associated with continuum cores, YSOs, and show a single peak and no extreme elongations (the major axis being smaller than 3 times the minor axis). The physical parameters (size, flux density, and peak intensity) of these selected compact structures are extracted from the NH₂D velocity-integrated intensity image using the CASA-imfit task (see Table 1). We stress that this identification is likely to be incomplete. Potential NH₂D compact structures could have been missed if they cannot be distinguished from the clumpy or filamentary structures. Four NH₂D cores (M11, M13, M14, M15) are not associated with continuum emission at all (<3σ, 1σ ∼ 30 μJy beam⁻¹ for continuum image), and the rest of the NH₂D cores are partially associated with weak continuum emission (3σ–7σ), except for M10 that has about 10% region at the edge of the core associated with 8σ–11σ continuum emission (Figure 2).

With the observed 3 mm continuum fluxes, we estimated the gas mass (M_gas) of identified cores assuming optically thin modified black body emission in the Rayleigh–Jeans limit, following

$$\begin{eqnarray}&&{M}_{\mathrm{gas}}=\eta \displaystyle \frac{{S}_{\nu }{d}^{2}}{{B}_{\nu }(T){\kappa }_{\nu }},\end{eqnarray} \tag{ B1 }$$

where η = 100 is the assumed gas-to-dust mass ratio, d is the source distance, S_ν is the continuum flux at frequency 98.5 GHz, B_ν (T) is the Planck function at temperature T, and κ_ν is the dust opacity at frequency ν. We adopt κ_{98.5 GHz} = 0.235 cm²g⁻¹, assuming κ_ν = 10(ν/1.2 THz)^β cm² g⁻¹ and β = 1.5 for all of dense cores (Hildebrand 1983).

The column density (${N}_{{{\rm{H}}}_{2}}$) is computed with

$$\begin{eqnarray}&&{N}_{{{\rm{H}}}_{2}}=\eta \displaystyle \frac{{S}_{\nu }^{\mathrm{beam}}}{{B}_{\nu }(T){\kappa }_{\nu }\,{\rm{\Omega }}\,\mu \,{m}_{{\rm{H}}}},\end{eqnarray} \tag{ B2 }$$

where ${S}_{\nu }^{\mathrm{beam}}$ is the continuum flux density, Ω is the beam solid angle, μ = 2.8 is the mean molecular weight of the interstellar medium (Kauffmann et al. 2008), and m_H is the proton mass.

The virial mass is estimated using

$$\begin{eqnarray}&&{M}_{\mathrm{vir}}=\displaystyle \frac{5}{a\beta }\displaystyle \frac{{\sigma }_{\mathrm{tot}}^{2}R}{G},\end{eqnarray} \tag{ B3 }$$

where ${\sigma }_{\mathrm{tot}}\ =\ \sqrt{{\sigma }_{\mathrm{nt},{\mathrm{NH}}_{2}{\rm{D}}}^{2}+{\sigma }_{\mathrm{th},\langle {\rm{m}}\rangle }^{2}}\ =\ \sqrt{({\sigma }_{\mathrm{obs},\mathrm{int}}^{2}-{\sigma }_{\mathrm{th},{\mathrm{NH}}_{2}{\rm{D}}}^{2})+{\sigma }_{\mathrm{th},\langle {\rm{m}}\rangle }^{2}}$ is the total velocity dispersion, R is the effective radius, G is the gravitational constant, parameter a equals to (1 − b/3)/(1 − 2b/5) for a power-law density profile ρ ∝ r^−b , and $\beta =(\arcsin e)/e$ is the geometry factor (see also Li et al. 2013, 2020b). Here 〈m〉 is 2.37m_H (Kauffmann et al. 2008), σ_obs,int is derived by the core-averaged spectrum, and b is assumed to be 2. For M1, M11, M13, and M15, σ_tot = σ_th,〈m〉(10 K) was used as the ${\sigma }_{\mathrm{obs},\mathrm{int}}\leqslant {\sigma }_{\mathrm{th},{\mathrm{NH}}_{2}{\rm{D}}}(10\,{\rm{K}})$.

The internal luminosity of dense object can be estimated by using the empirical relationship between internal luminosity (L_int) and observed 70 μm flux (Dunham et al. 2008):

$$\begin{eqnarray}&&{L}_{\mathrm{int}}({L}_{\odot })=3.3\times {10}^{8}{\left({F}_{70}{\left(\displaystyle \frac{d}{140\mathrm{pc}}\right)}^{2}\right)}^{0.94},\end{eqnarray} \tag{ C1 }$$

where F₇₀ is the flux at 70 μm in cgs units (erg s⁻¹ cm⁻² Hz⁻¹) and d is the source distance in pc. A 3σ sensitivity of ∼0.3 Jy from Herschel/PACS 70 μm observation is used to derived an upper limit (∼1.26 L_⊙) of internal luminosity for the identified NH₂D cores (Tigé et al. 2017), with the exception of M3, M5, and M17, which are illuminated by the western YSOs rather than internal heating. The measured 70 μm fluxes of M3, M5, and M7 are about 5–6σ.