The ALMA Survey of 70 {\mu}m Dark High-mass Clumps in Early Stages (ASHES). X: Hot Gas Reveals Deeply Embedded Star Formation

Massive infrared dark clouds (IRDCs) are considered to host the earliest stages of high-mass star formation. In particular, 70 $\mu$m dark IRDCs are the colder and more quiescent clouds. At a scale of about 5000 au using formaldehyde (H2CO) emission, we investigate the kinetic temperature of dense cores in 12 IRDCs obtained from the pilot ALMA Survey of 70 $\mu$m dark High-mass clumps in Early Stages (ASHES). Compared to 1.3 mm dust continuum and other molecular lines, such as C18O and deuterated species, we find that H2CO is mainly sensitive to low-velocity outflow components rather than to quiescent gas expected in the early phases of star formation. The kinetic temperatures of these components range from 26 to 300 K. The Mach number reaches about 15 with an average value of about 4, suggesting that the velocity distribution of gas traced by H2CO is significantly influenced by a supersonic non-thermal component. In addition, we detect warm line emission from HC3N and OCS in 14 protostellar cores, which requires high excitation temperatures (Eu/k ~ 100 K). These results show that some of the embedded cores in the ASHES fields are in an advanced evolutionary stage, previously unexpected for 70 $\mu$m dark IRDCs.


INTRODUCTION
High-mass stars play a major role in the energy budget of galaxies via their radiation, wind, and supernova events, but the picture of their formation remains unclear.Massive infrared dark clouds (IRDCs), dark silhouettes against the Galactic 8 µm mid-infrared background in Galactic plane surveys (e.g., Simon et al. 2006;Peretto & Fuller 2009), are suitable targets for studying the earliest stages of high-mass star formation (Rathborne et al. 2006;Chambers et al. 2009;Sanhueza et al. 2012;Bovino et al. 2019;Sabatini et al. 2020;Redaelli et al. 2021).Among IRDCs, those that are 70 µm dark are especially important because they are the coldest, most quiescent clouds (e.g., Guzmán et al. 2015).They are considered to trace the earliest stages of high-mass star formation (e.g., Zhang et al. 2009;Tan et al. 2013;Sanhueza et al. 2013Sanhueza et al. , 2017;;Contreras et al. 2018;Sabatini et al. 2019).However, previous studies have found that not all cores in 70 µm dark IRDCs are in the prestellar phase.Some evidence of active star formation is revealed by the detection of molecular outflows (e.g., Sanhueza et al. 2010;Wang et al. 2011;Traficante et al. 2015;Feng et al. 2016;Svoboda et al. 2019;Tan et al. 2016;Li et al. 2019;Kong et al. 2019;Morii et al. 2021;Tafoya et al. 2021;Redaelli et al. 2022;Sabatini et al. 2023;Jiao et al. 2023).
We have carried out the pilot Atacama Large Millimeter/submillimeter Array (ALMA) Survey of 70 µm dark High-mass clumps in Early Stages (ASHES) toward 12 IRDC clumps.The sample was selected by combining the ATLASGAL survey (Schuller et al. 2009;Contreras et al. 2013) and a series of studies from the MALT90 survey (Foster et al. 2011(Foster et al. , 2013;;Jackson et al. 2013;Hoq et al. 2013;Guzmán et al. 2015;Rathborne et al. 2016;Contreras et al. 2017;Whitaker et al. 2017).The source selection is described in Sanhueza et al. (2019) and Morii et al. (2023) 1 .The names of the 12 targets are listed in Table 1.From the pilot survey, 294 dense cores were detected in continuum emission using the astropy astrodendro package (Rosolowsky et al. 2008;Robitaille et al. 2013).Among them, 97 cores are thought to be in the protostellar stage (hereafter, protostellar cores) based on the detection of molecular outflows and/or emission from any of the three lines requiring warm excitation temperatures: methanol (CH 3 OH) , and H 2 CO J K A K C = 3 2,1 -2 2,0 (E u /k = 68.11K).The other 197 dense cores are prestellar candidates (hereafter, prestellar cores).In addition, 97 protostellar cores are classified into three evolutionary stages based on searching molecular outflows and/or warm-core line emission.They are named "outflow core", "warm core", and "warm and outflow core".The detailed explanation of the core classification is summarized in Li et al. (2022).Using these data, we have investigatedthe fragmentation process (Sanhueza et al. 2019), outflows (Li et al. 2020;Tafoya et al. 2021), chemical properties (Sakai et al. 2022;Sabatini et al. 2022;Li et al. 2022), active star formation signatures (Morii et al. 2021), and dynamical properties (Li et al. 2023).The summary of the whole ASHES sample has been recently presented in (Morii et al. 2023).
The metastable lines of ammonia (NH 3 ) are frequently used as a standard molecular cloud thermometer within our Galaxy and also in external galaxies (e.g., Ho & Townes 1983;Walmsley & Ungerechts 1993;Danby et al. 1988;Friesen et al. 2009;Mangum et al. 2013b;Keown et al. 2019).However, the NH 3 abundance can vary strongly in different molecular environments (e.g., 10 −5 in dense molecular "hot cores" around newly formed massive stars, Mauersberger et al. 1987; 10 −8 in dark clouds, Benson & Myers 1983;Chira et al. 2013; ∼ 10 −10 in the Orion Bar PDR, Batrla & Wilson 2003) and furthermore, is can be extremely affected by a high UV flux.In contrast to this, H 2 CO has a relatively constant fractional abundance, with variations rarely exceeding one order of magnitude during different stages of star formation (e.g., Mangum et al. 1990;Mangum & Wootten 1993;Caselli et al. 1993;Johnstone et al. 2003;Gerner et al. 2014;Tang et al. 2017aTang et al. ,b, 2018a;;Zhu et al. 2020;Sabatini et al. 2021).For instance, the H 2 CO abundance is the same in the hot core and in the compact ridge of the Orion nebula (Caselli et al. 1993;Mangum & Wootten 1993).The NH 3 (2,2)/(1,1) ratio is sensitive to gas temperatures with temperatures of less than 50 K (Mangum et al. 2013b;Gong et al. 2015).This is similar to the kinetic temperature range that the H 2 CO ratio is most sensitive to (Mangum & Wootten 1993).The NH 3 lines have lower effective excitation densities than the H 2 CO transitions by a few orders of magnitude, n eff (NH 3 (1,1)) ∼ 10 3 cm −1 while n eff (H 2 CO(3 0,3 -2 0,2 )) ∼ 10 5 cm −1 (Shirley 2015).Previous observations reported that the temperatures derived from NH 3 reflect an average temperature of more extended and cooler gas (Henkel et al. 1987;Ginsburg et al. 2016) while H 2 CO traces denser and hotter regions (Ginsburg et al. 2016;Tang et al. 2017aTang et al. , 2018a,b),b).From the above, the three transitions of H 2 CO allow us to probe the kinetic temperature of IRDCs and unveil potential precursors of future massive hot cores.
The paper is organized as follows.The observations are presented in Section 2. Section 3 describes the H 2 CO results and analysis in this IRDC sample.In Section 4, we discuss the kinetic temperatures derived from H 2 CO.Our main conclusions are summarized in Section 5.

OBSERVATIONS
Observations of twelve 70 µm IRDC clumps were performed with ALMA in Band 6 (∼224 GHz; 1.34 mm) using the main 12 m array, the 7 m array, and the total power (TP) array (Project ID: 2015.1.01539.S, PI: P. Sanhueza).The mosaic observations were carried out with the 12 m array and 7 m array to cover a significant portion of the clumps, as defined by single-dish continuum images.The same correlator setup was used for all sources.More details on the observations have been provided in earlier articles of the ASHES series (Sanhueza et al. 2019;Li et al. 2020).
The data calibration was performed using the CASA software package version 4.5.3,4.6, and 4.7, while imaging was carried out using CASA 5.4 (McMullin et al. 2007).The data cubes for lines were produced using the automatic masking procedure YCLEAN (Contreras et al. 2018).Since some sources were observed with different array configurations (Sanhueza et al. 2019), a uv-taper was used for such sources in order to obtain a similar synthesized beam of ∼ 1. ′′ 2 for all sources.We adopted a Briggs robust weighting of 0.5 and 2 for the visibilities of continuum and lines in the imaging process, respectively.This yielded an averaged 1σ root mean square (rms) noise level of ∼ 0.1 mJy beam −1 for continuum images.For H 2 CO lines, the sensitivity is ∼ 3.3 mJy beam −1 per 0.67 km s −1 .

RESULTS
Figure 1 shows an example of an observed spectrum including the H 2 CO transition lines.The peak intensities, line widths, and integrated intensities of the para-H 2 CO emission were derived using Gaussian fitting (red curve in Figure 1).
For H 2 CO(3 2,2 -2 2,1 ) and H 2 CO(3 2,1 -2 2,0 ), the detection rates of protostellar cores are both 46% (45/97).These values are about half of the one for H 2 CO(3 0,3 -2 0,2 ).We note that the detection threshold is 3σ. Figure 3 shows the relation between 1.3 mm flux density (S 1.3mm ) and the detection rate of H 2 CO(3 0,3 -2 0,2 ) for protostellar and prestellar core candidates.The plot indicates that the detection rate of both protostellar and prestellar core candidates grows with increasing S 1.3mm .However, a larger sample is needed for a fully robust statistical assessment.The details on the detection rates are summarized in Table 1.
Figure 5 presents the relation between S 1.3mm and H 2 CO(3 2,2 -2 2,1 ) integrated intensity (I H2CO(32,1) ; left panel) and between S 1.3mm and H 2 CO(3 2,1 -2 2,0 ) integrated intensity (I H2CO(32,1−22,0) ; right panel) for all protostellar cores.Weak correlations are also seen be-3 Spearman's correlation is a nonparametric measure of the strength and direction of association that exists between two variables measured on at least an ordinal scale.The Spearman's coefficient, r, ranges from -1 to 1, with 0 indicating no correlation.The value of 1 implies an exact increasing monotonic relation between two quantities, while -1 implies an exact decreasing monotonic relation.The probability value, p-value, obtained from the calculator is a measure of how likely or probable it is that any observed correlation is due to chance.The p-value ranges between 0 (0 %) to 1 (100 %).The p-value close to 1 suggests no correlation other than due to chance and that the null hypothesis assumption is correct.While the p-value close to 0 suggests that the observed correlation is unlikely to be due to chance, and there is a very high probability that the null hypothesis is wrong.

Velocity distribution
The H 2 CO(3 0,3 -2 0,2 ) intensity-weighted mean velocity (moment 1; left panel), and velocity dispersion (moment 2; right panel) maps for two representative IRDC clumps: G010.991-00.082 and G337.541-00.082 are shown in Figure 6.The velocity range used to generate these moment maps for the 12 IRDC clumps is approximately 5-25 km s −1 .The central velocity values are provided in Table 1 in Sanhueza et al. (2019).We note that these velocity ranges were determined through visual inspection of the spectra for each clump, ensuring that they encompassed all of the line emission.The figure reveals complex structures, including different velocity components and large velocity dispersions.For instance, the difference between the largest and smallest peak velocities estimated from the moment 1 map reaches to more than 5 km s −1 .The largest velocity dispersion estimated from the moment 2 map is more than 2 km s −1 .
The blue-and red-shifted components of H 2 CO(3 0,3 -2 0,2 ) are shown in Figure 9 and appendix (Section A and Figure 20 ).The velocity range of the blue-and red-shifted components, which are selected manually by visually inspecting the spectra, are written in the upperright on Figure 9.In G340.179-00.0242 and G340.222-00.167,blue-and red-shifted components are not detected (emission is less than 3σ).All components seem to be associated with protostellar cores.The direction of the blue-and red-shifted components roughly agrees with these components detected in CO emission that is considered to trace molecular outflows, as investigated by Li et al. (2020) in the ASHES survey.The velocity and spatial range of these components detected in H 2 CO(3 0,3 -2 0,2 ) emission (1.5-20 km s −1 , 0-0.3 pc) is smaller than those detected in CO emission (9-95 km s −1 , 0.04-0.45pc; Li et al. 2020).Therefore, these H 2 CO velocity components are considered to be tracing lowvelocity outflows.This consideration is consistent with the result of a recent study on chemical diagnostics of protostellar sources in the low-mass star-forming regions with ALMA ( Lukasz Tychoniec et al. 2021).Among the 59 outflows identified by Li et al. (2020), 30 outflows are detected in H 2 CO(3 0,3 -2 0,2 ) emission (solid arrows in Figure 9).In G343.489-00.416,two candidate of blueshifted outflows are newly detected by H 2 CO(3 0,3 -2 0,2 ) emission (black arrows in Figure 9).Figure 10 shows examples of the position-velocity map of the outflows.

LTE assumption
We fit the three H 2 CO lines assuming LTE conditions and estimate best-fitting parameters.Under this assumption, the fitting routine constructs line profiles with a range of input parameters, including peak velocities, temperatures, line-widths, and column densities, following Mangum & Shirley (2015).The peak velocities of the lines can be selected manually and are fixed according to their known rest frequencies.The code then minimizes the differences between the constructed and observed lines, using the non-linear least-squares python fitting routine lmfit4 .The best fit is returned, along with the estimated temperature, line-width, column density, and associated errors, which we then take as the optimized fitting results5 .The fitted temperatures, column densities, and line-widths are in the range of 10 -300 K, 1.0 × 10 10 -1.0 × 10 20 cm −2 , and 0.14 -1.37 km s −1 , respectively.When we derive the temperature for each pixel, the peak velocities are automatically selected as the position of H 2 CO(3 0,3 -2 0,2 ) peak intensity and are fixed according to their known rest frequencies.

Radex modeling
We employ the Radex code (van der Tak et al. 2007), by means of the ndRADEX tool6 , to generate a set of different models covering a wide range of conditions.As pointed out by Ao et al. (2013) and appendix (Section B and Figure 22), the H 2 CO(3 2,1 -2 2,0 )/H 2 CO(3 0,3 -2 0,2 ) ratio is more sensitive to densities than the H 2 CO(3 2,2 -2 2,1 )/H 2 CO(3 0,3 -2 0,2 ) ratio, and is therefore used to derive kinetic temperatures in the following.The Radex code needs five input parameters: 1) the background temperature, 2) the kinetic temperature, 3) the H 2 volume density, 4) the H 2 (or H 2 CO) column density, and 5) the line width.For the background temperature, we adopt 2.73 K.The line width is derived from the Gaussian fitting.For the H 2 density and H 2 (or H 2 CO) column density, we adopt the values derived from the 1.3 mm dust continuum (Sanhueza et al. 2019).The fractional abundance of N(H 2 CO)/N(H 2 ) was adopted as 1.0 × 10 −10 from Tang et al. (2017a), which investigated the fractional abundance in IRDCs.When we derive the temperature of each pixel, we use the averaged values of H 2 density, column density, and linewidth over the whole IRDC for simplicity.With the above parameters and the H 2 CO(3 2,1 -2 2,0 )/H 2 CO(3 0,3 -2 0,2 ) integrated intensity ratio derived from the Gaussian fitting, we then compute the temperatures using a model grid for the H 2 CO lines encompassing 580 temperatures ranging from 10 to 300 K.The temperature error is propagated from the optimum of the Gaussian fitting.Radex modeling (T Radex ) for all protostellar cores.Due to the 5σ detection limit, we derived temperatures for 26 among the total of 97 protostellar cores.Figure 12 shows the T LTE and T Radex distribution in two examples of IRDC clumps.For both T LTE and T Radex , the temperature error becomes larger than 20 K for temperatures of more than 100 K.The T LTE and T Radex are roughly similar, but in almost all cases, the T Radex is slightly higher than the T LTE .The difference likely results from (1) the difference between fitting the H 2 CO(3 2,2 -2 2,1 )/H 2 CO(3 0,3 -2 0,2 ) integrated intensity ratio for the Radex modeling and the LTE assumption, (2) the assumption of the H 2 volume (n(H 2 )) and column densities(N(H 2 )) for the Radex modeling, and (3) the LTE assumption.We note that the effect from the assumption of the N(H 2 ) is negligible compared to the assumption of the n(H 2 ) for the Radex modeling (see Section B and Figure 22).These results indicate that the temperature determination can vary significantly depending on the assumed models and n(H 2 ). Figure 13 shows the relation between T LTE and T Radex and the H 2 CO(3 2,1 -2 2,0 )/H 2 CO(3 0,3 -2 0,2 ) integrated intensity ratio.This plot indicates that the temperature determination can vary significantly depending on the assumed models and n(H 2 ).To investigate the uncertainty, we calculate the minimum and maximum temperature (gray filled area in Figure 13) based on the assumption of possible n(H 2 ) and N(H 2 ) ranges: n(H 2 ) = 10 5 -2.0 × 107 cm −3 and N(H 2 ) = 10 22 -10 24 cm −2 .These ranges come from the minimum and maximum values of the cores derived from the 1.3 mm continuum emission (Sanhueza et al. 2019), which include both protostellar cores and prestellar candidates.The uncertainty increases with growing H 2 CO(3 2,1 -2 2,0 )/H 2 CO(3 0,3 -2 0,2 ) integrated intensity ratios.For instance, the uncertainty at the ratio of 0.4 is about 100 K, which is roughly 5 times larger than the one at the ratio of 0.2.In all models and under all assumptions, more than half of the regions traced by H 2 CO show high temperatures of more than 50 K.
In order to investigate the impact of the missing flux on the line emission, we compare the temperatures derived from 12m7m and 12m7mTP data.The ratio of T LTE derived from 12m7m per 12m7mTP data is 0.74-2.06with an average of 1.03±0.06.The ratio of T Radex is 0.19-6.52 with an average of 1.04±0.17.Compared to the uncertainty from the assumption, the impact of the missing flux is likely negligible.The T LTE and T Radex values across all cores range from 36 to 300 K (average: 84±52 K) 7 and from 37 to 292 K (average: 98±57 K), respectively.The T LTE and T Radex temperatures measured in each pixel range from 29 to 300 K (average: 75±50 K) and from 26 to 300 K (average: 91±60 K), respectively.We stress that these values are derived where all three H 2 CO lines are detected above 5σ.Therefore, this also includes regions outside of cores.Previous studies reported that the temperatures derived from H 2 CO lines in IRDCs and outflow objects range between 40 and 60 K and between 31 and 100 K, respectively (Tang et al. 2017a).These previous observations were performed by single-dish telescopes, such as the JCMT, with spatial resolutions around 20-30 ′′ (corresponding to ∼ 0.5 pc; Tang et al. 2017a).This means that their temperatures are considered to be an average temperature of the entire area of the object.Owing to the much higher resolution in the work here, we could resolve the IRDCs into several components, including jets, outflows, and cores ( Lukasz Tychoniec et al. 2021).Therefore, the temperature range of our results is much wider than what is seen in the previous studies with single-dish telescopes.We note that a recent ALMA study targeting the IRDC G10.21-0.31,with a spatial resolution of about 5000 au, reported temperatures derived from H 2 CO lines in two dense cores as 67.7 K and 83.0 K (Jian et al., 2023), consistent with the range of values found in the ASHES sample.

Comparison of temperatures derived from H 2 CO and NH 3
For this comparison, the NH 3 (1,1) and (2,2) lines in IRDCs are obtained from the Complete ATCA (The Australia Telescope Compact Array) Census of High-Mass Clumps survey (CACHMC survey;Allingham et al., 2023 in prep) at an angular resolution of ∼5 ′′ .In order to compare the kinetic temperatures derived from para-H 2 CO lines and NH 3 lines (Li et al. 2022, Allingham et al., 2023 in prep), we smoothed the H 2 CO data cube to a half-power beam width (HPBW) of ∼5 ′′ with a 2D Gaussian function.
Figure 14 shows examples of spatial temperature distributions derived from H 2 CO and NH 3 with a resolution of ∼5 ′′ .The spatial distributions differ significantly.NH 3 shows a much more extended distribution than H 2 CO.The T LTE , T Radex , and NH 3 temperatures across all IRDCs are 24-300 (average: 78±65), 29-300 (average: 89±66), and 7-32 (average: 14±5), respectively.Figure 15 compares temperatures derived from H 2 CO and NH 3 for each IRDC.We only compare the H 2 CO and NH 3 temperatures in those region where both molecules are detected.The H 2 CO temperature of G028.273-00.167 is not derived in the smoothed data because the H 2 CO components in G028.273-00.167detected with a 1. ′′ 2 resolution (see Figure 19 in appendix) are diluted at a 5 ′′ beam.These results from the comparison of H 2 CO and NH 3 are similar to those for outflows and IRDCs from Tang et al. (2017a).However, the H 2 CO temperature range in our results is much larger than in Tang et al. (2017a), and most likely due to our originally higher spatial resolution, despite the factor of smoothing.The NH 3 temperature range is similar to the previous observations with the Effelsberg telescope (∼40 ′′ ; Tang et al. 2017a), despite the resolution being higher by a factor of about 10.These results indicate that in these IRDCs H 2 CO traces components that are different from what is seen in NH 3 .Moreover, even if at the same resolution of 5 ′′ the variation of the temperatures derived from H 2 CO is much larger than the variation of the NH 3 temperatures.This finding might indicate that H 2 CO is more sensitive to star-formation activities, including cores, outflows, and jets, than NH 3 in IRDCs.

Thermal and non-thermal motions
We calculated thermal linewidth (σ T ), non-thermal velocity dispersion (σ N T ), thermal sound speed (a s ), ratio of thermal to non-thermal pressure (R p ), and Mach number (M).The σ T , σ N T , a s , R p , and M are given by (2) where k is the Boltzmann constant, T kin the kinetic temperature of the gas (T LTE or T Radex ), m x the mass of the relevant molecule, ∆V the measured FWHM linewidth, µ = 2.37 the mean molecular weight for molecular clouds, and m H the mass of the hydrogen atom.The FWHM linewidth is measured by two methods: fitting with the LTE assumption (see Section 3.5.1)and the Gaussian fitting.The errors of these quantities are derived by varying the assumed T kin from minimum to maximum values.We note that we only considered the uncertainties (i.e., from minimum to maximum values) of T kin because the uncertainties for the measured FWHM are comparatively small and they can be neglected.

Non-thermal velocity dispersion
The left panel of Figure 16 shows the relation between temperature and σ N T derived from the H 2 CO(3 0,3 -2 0,2 ) line (σ NT(H2CO(30,3−20,2)) ).There is a positive correlation between them.This result is consistent with previous observations of other IRDCs and starforming regions using single-dish telescopes (temperature ∝ σ 0.66−1.26NT(H2CO(30,3−20,2)) ), which suggests that the gas traced by H 2 CO is turbulent and might be subject to turbulent heating (e.g., Immer et al. 2016;Tang et al. 2017aTang et al. , 2018a,b),b).The right panel of Figure 16 shows the same plot as the left panel but excludes the data from G014.492-00.139(hereafter, G014.49).In this case, we find an even stronger positive correlation.This might be due to G014.49having the most complicated velocity structure in this sample (Redaelli et al. 2022).The spectra for many of the cores in G014.19 are different from a single Gaussian (e.g., double-peak structure).Therefore, the fits for the cores in G014.49are likely more uncertain than for other cores.Spearman's correlation coefficients and least-square fitting results are summarized in Table 5.
Figure 17 is the same as Figure 16, but for pixels.We also confirm a positive correlation between T kin and σ NT(H2CO(30,3−20,2)) , although there is large scattering.
The fitting results are also summarized in Table 5.We also find a weak trend in this distribution.The majority of the pixels is located in the area where 1.0 ≤ σ NT(H2CO(30,3−20,2)) ≤ 2.0 km s −1 and 35.0 ≤ T kin ≤ 75.0 K (hereafter, main-area; black dotted box in Figure 17).If the data are derived based on the LTE assumption or based on Radex modeling, 41 % or 39 % of the pixels are located in the main-area, respectively.Discarding G014.49(right panel of Figure 17), it is 56 % or 55 %.In addition to the main-area, some data points are also concentrated in the area where 2.0 ≤ σ NT(H2CO(30,3−20,2)) ≤ 3.0 km s −1 and 60.0 ≤ T kin ≤ 200.0 K (hereafter, sub-area; black dashed box in Figure 17).In the sub-area, the T kin rapidly increases with increasing σ NT(H2CO(30,3−20,2)) .Both when derived under the LTE assumption or with Radex modeling, about 13% of the pixels are located in the sub-area.Without G014.49(right panel of Figure 17), this becomes about 12 % (both LTE and Radex).In space, both data points located in the main-area and sub-area mainly trace outflow structures.However, the data points in the subarea tend to trace the edges of the outflows.This result suggests that the kinematics along the edges of outflows are different from the main part of the outflows.

Mach number
Figure 18 shows the spatial distribution of M for two examples of IRDC clumps.The range of M derived from LTE and Radex temperatures is 1.5-12.3(average: 4.5±2.4)and 1.3-9.0(average: 4.3±2.4),respectively (e.g., gray solid, dashed, and dotted lines in Figures 16  and 17).These average values reveal a supersonic environment, while the largest values point at a hypersonic setting.This finding suggests that the velocity distributions of the gas traced by H 2 CO are significantly influenced by supersonic non-thermal components (e.g., turbulent motions, infall, outflows, shocks, and/or magnetic fields).Furthermore, we find that M is larger at the edges of the outflows in G337.54 and G343.48 (Figure 18).The values reach to more than 5, corresponding to hypersonic conditions.Such high Mach numbers might be caused by shocks between outflows and the surrounding gas.We note that high values for M are also detected in G014.49, which has several outflow components.However, the values for G014.49might be more uncertain due to the complicated velocity structures, as discussed in Section 4.3.1.

Global connection to star formation
In summary, by utilizing the three transitions of H 2 CO, we have successfully derived kinetic temperature distributions for 9 out of 12 70 µm dark IRDCs (26 out of 294 embedded cores in IRDCs), which are regarded as the coldest and most quiescent clouds (e.g., Guzmán et al. 2015).The temperature range spans from 30 K to 300 K, with more than half of the region exhibiting temperatures higher than 50 K.Such high-temperature gas is considered to be associated with a low-velocity outflow component, which is significantly influenced by supersonic non-thermal motions.Furthermore, 14 out of 294 cores show detections in HC 3 N and/or OCS, which require excitation tempertures of more than 100 K These results suggest that some of the embedded cores in 70 µm dark IRDCs have already entered the protostellar phase and may serve as hosts for hot core precursors.

SUMMARY
We report kinetic temperatures of 12 IRDC clumps derived from the three transitions of H 2 CO with a resolution of about 5000 au.The main results are as follows.
The distributions of all three H 2 CO lines is different from that of the 1.3 mm dust continuum emission.This result is not consistent with previous observations of IRDCs from single-dish telescopes with a resolution of about 0.5 pc (e.g., Tang et al. 2017aTang et al. , 2018a).
2. The population of protostellar cores shows a weak positive correlation between the 1.3 mm flux density and the H 2 CO(3 0,3 -2 0,2 ) integrated intensity.
In contrast, the results obtained for prestellar cores show no clear correlation.Therefore, earlier results reporting a positive correlation based on single-dish telescopes (e.g., Tang et al. 2017aTang et al. , 2018a) ) suggest that their targets are dominated by evolved regions which have already entered the protostellar phase.
3. The velocity dispersion of H 2 CO is much larger than that of other molecular lines, such as C 18 O, DCO + , and N 2 D + (Li et al. 2022).We confirmed that H 2 CO mainly traces the low-velocity outflow components rather than the quiescent gas in the early phase of high-mass star formation.This is consistent with results for low-mass star-forming regions (e.g., Lukasz Tychoniec et al. 2021).
4. The kinetic temperatures across all IRDCs range from 26 to 300 K.This is much larger than what is reported in previous studies using single-dish telescopes (e.g., Tang et al. 2017a).This temperature range is also larger than the temperatures derived from NH 3 with the same resolution.This result suggests that H 2 CO is more sensitive to star-formation activity than NH 3 .
5. The Mach numbers derived from H 2 CO reach about 15 with an average of about 4. This indicates that the velocity distribution of gas traced by H 2 CO is significantly influenced by supersonic non-thermal components.High Mach numbers of more than 5 are mainly detected at the edges of outflows.These high numbers might be caused by shocks between outflows and surrounding material.
6. 14 protostellar cores show emission from warmer lines, such as HC 3 N and OCS.This suggests that some cores embedded in 70 µm dark IRDCs may be developing hot core precursors given than these lines are known hot-core tracers (e.g., Chapman et al. 2009;Lukasz Tychoniec et al. 2021).
b Pixels with T kin = 300 K are discarded because these values are considered to be lower limits.

Figure 1 .
Figure 1.Example of detected spectrum in IRDC clumps (clump: G337.541-00.082,core: ALMA1).The black curve is the observed spectrum.The red curve indicates the results of the Gaussian fitting.

Figure 8 .
Figure 8. Relation between the velocity dispersion of C 18 O(2-1) (σ C 18 O ) and H2CO (σ H 2 CO(3 0,3 −2 0,2 ) ) for all IRDC clumps ((a): all cores, (b) only protostellar cores, and (c) pixels).The gray circles and red squares in the left panel indicate the prestellar candidates and protostellar cores, respectively.The gray squares, blue triangles, and yellow circles in the right panel indicate the three classified protostellar cores: outflow core, warm core, and warm and outflow core, respectively.The pixel scale was re-binned to 0. ′′ 6, approximately half of the synthesized beam, from the original value (0. ′′ 2).

Figure 9 .
Figure 9. Blue-(blue contours) and red-shifted (red contours) components of H2CO(30,3-30,2) for four representative IRDCs.The gray-scale background shows the 1.3 mm dust continuum.The cyan and orange circles indicate the positions of prestellar candidates and protostellar cores, respectively, identified in ASHES (Li et al. 2022).Synthesized beams are displayed at the bottom left of each panel (open ellipses: H2CO; filled ellipses: continuum).The red and blue arrows indicate the outflow directions by CO and SiO emission lines in Li et al. (2020), which investigates the outflow structure with the same ASHES data.(solid arrows: also detected in H2CO(30,3-30,2), dotted arrows: not detected in H2CO(30,3-30,2)).The black arrows indicate new outflow candidates detected in H2CO(30,3-30,2), which are not seen in CO or SiO inLi et al. (2020).The black dotted boxes mark the regions for the position-velocity diagrams in Figure10.The images of the other 5 IRDC clumps are presented in the appendix (Figure20).

Figure 10 .
Figure 10.Position-velocity diagrams of H2CO(30,3-30,2) for regions labeled as (a), (b), and (c) in Figure 9.The vertical and horizontal dashed magenta lines represent the positions and systemic velocities of the dense cores, respectively.

Figure 15 .
Figure 15.Comparison of kinetic temperatures (T kin ) derived from H2CO (red circle: LTE, blue square: Radex) and NH3 (yellow triangle) for each IRDC clump.The error bars indicate the standard deviation of T kin (1 σ).The extended bars indicate the minimum-maximum range of T kin .To derive these values, the pixel scale was re-binned to 2. ′′ 4, approximately half of the angular resolution of NH3 data (∼ 5 ′′ ), from the original value (0. ′′ 2).

Table 1 .
Core detection rate of three H2CO emission lines.

Table 6 .
Properties of warm lines.