# Probing the physics of star formation (ProPStar): III. No evidence of dissipation of turbulence down to 20 mpc (4000 au) scale

Jaime E. Pineda, Juan D. Soler, Stella Offner, Eric W. Koch, Dominique M.

Segura-Cox, Roberto Neri, Michael Kuffmeier, Alexei V. Ivlev, Maria Teresa

Valdivia-Mena, Olli Sipilä, et al.

<figure><img src="image_1.png" /></figure>

► To cite this version:

Jaime E. Pineda, Juan D. Soler, Stella Offner, Eric W. Koch, Dominique M. Segura-Cox, et al. Probing the physics of star formation (ProPStar): III. No evidence of dissipation of turbulence down to 20 mpc (4000 au) scale. Astronomy & Astrophysics - A&A, 2024, 690, 10.1051/0004-6361/202451208. insu-04836805

HAL Id: insu-04836805

https://insu.hal.science/insu-04836805v1

Submitted on 15 Dec 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>
LETTER TO THE EDITOR

Probing the physics of star formation (ProPStar)

III. No evidence of dissipation of turbulence down to 20 mpc (4000 au) scale*

Jaime E. Pineda¹,**, D, Juan D. Soler², D, Stella Offner³, D, Eric W. Koch⁴, D, Dominique M. Segura-Cox³,¹,**** D,
Roberto Neri⁵, D, Michael Kuffmeier⁶,¹, D, Alexei V. Ivlev¹, D, Maria Teresa Valdivia-Mena¹, D, Olli Sipilä¹, D,
Maria Jose Maureira¹, D, Paola Caselli¹, D, Nichol Cunningham⁷,⁸, D, Anika Schmiedeke⁹, D, Caroline Gieser¹, D,
Michael Chen¹⁰, D, and Silvia Spezzano¹, D

(Affiliations can be found after the references)

Received 21 June 2024 / Accepted 26 August 2024

ABSTRACT

*Context.* Turbulence is a key component of molecular cloud structure. It is usually described by a cascade of energy down to the dissipation scale. The power spectrum for subsonic incompressible turbulence is ∝k⁻⁵/³, while for supersonic turbulence it is ∝k⁻².

**Aims.** We determine the power spectrum in an actively star-forming molecular cloud, from parsec scales down to the expected magnetohydrodynamical (MHD) wave cutoff (dissipation scale).

**Methods.** We analyzed observations of the nearby NGC 1333 star-forming region in three different tracers to cover the different scales from ~10 pc down to 20 mpc. The largest scales are covered with the low-density gas tracer $^{13}$CO (1–0) obtained with a single dish, the intermediate scales are covered with single-dish observations of the C$^{18}$O (3–2) line, while the smallest scales are covered in H$^{13}$CO$^+$ (1–0) and HNC (1–0) with a combination of NOEMA interferometer and IRAM 30m single-dish observations. The complementarity of these observations enables us to generate a combined power spectrum covering more than two orders of magnitude in spatial scale.

**Results.** We derive the power spectrum in an active star-forming region spanning more than 2 decades of spatial scales. The power spectrum of the intensity maps shows a single power-law behavior, with an exponent of 2.9 ± 0.1 and no evidence of dissipation. Moreover, there is evidence that the power spectrum of the ions to have more power at smaller scales than the neutrals, which is opposite to the theoretical expectations.

**Conclusions.** We show new possibilities for studying the dissipation of energy at small scales in star-forming regions provided by interferometric observations.

**Key words.** turbulence – stars: formation – ISM: clouds – ISM: molecules – ISM: individual objects: NGC 1333

**1. Introduction**

Dense cores are the places where stars form (Pineda et al. 2023). These are the places with subsonic levels of turbulence (Goodman et al. 1998; Pineda et al. 2010; Friesen et al. 2017), and this transition to coherence has been discussed as possibly being linked to the dissipation of turbulence.

Turbulence in clouds has been typically studied using different tracers, which probe different scales of the cloud. For example, Larson (1981) found a correlation between linewidth and cloud size, $\sigma_v \propto L^{0.38}$. Further analysis using principal component analysis (PCA) showed that the velocity structure functions as a function of spatial scale, $l$, and is universal in molecular clouds (Heyer & Brunt 2004), $\delta v \propto l^{0.49}$; although further work showed some small variation (Heyer et al. 2009; Roman-Duval et al. 2011), the relation is approximately universal (see also references in McKee & Ostriker 2007; Heyer & Dame 2015).

In particular, the spatial power spectrum is used to study the lower density section of the cold neutral medium.

Miville-Deschênes et al. (2010) measured the spatial power spectrum in the Polaris cloud using Herschel data, while Miville-Deschênes et al. (2016) showed a lack of evidence for the cutoff. These works focused on cirrus clouds, where gravity is not expected to play a major role.

In the case of star-forming molecular clouds, various works predict a turnover in the power spectrum at the dissipation scale (e.g., Elmegreen & Scalo 2004). Moreover, in the case of magnetohydrodynamic (MHD) turbulence, ions and neutrals at higher angular resolution are expected to display a difference in their linewidths at the same characteristic scale $l$ (Li & Houde 2008). This is because there is a minimum MHD wave for wave propagation ($\lambda_{min}$), below which MHD waves will not transmit, which would imply that ions should exhibit a turnover in the power spectrum at larger scales than neutrals.

In this study we present a comprehensive analysis of the power spectrum of spatial frequencies of gas in the NGC 1333 young cluster in the Perseus molecular cloud. We combine three data sets to study the power spectrum from parsec scales down to 20 mpc (4000 au).

**2. Data**

In order to properly probe the power spectrum over a wide spatial dynamical range we used three different data sets with different

<footnote>* Based on observations carried out under project number S21AD with the IRAM NOEMA Interferometer and 090-21 with the IRAM 30m telescope. IRAM is supported by INSU/CNRS (France), MPG (Germany) and IGN (Spain).</footnote>

<footnote>** Corresponding author; jpineda@mpe.mpg.de</footnote>

<footnote>*** NSF Astronomy and Astrophysics Postdoctoral Fellow.</footnote>
tracers. This enabled us to determine the power spectrum over an extended range of physical scales, which combined enabled us to cover more than two decades of physical scales. Figure 1 shows the integrated intensity maps used to analyze the scales covered, as well as the footprint of the higher angular resolution (and density tracer) observations used for the next spatial scale. The median spectra for the different data sets are presented in Appendix B and are shown in Fig. B.1.

For the largest scale, we used the ¹³CO (1–0) map of Perseus obtained with the FCRAO telescope by the COMPLETE survey (Ridge et al. 2006). This data has a resolution of 46" and a noise level of 0.12 K.

The integrated intensity was calculated by adding the emission in the velocity range of −1 and 11.5 km s⁻<fnref n="1" />. The top left panel of Fig. 1 shows the integrated intensity map.

For intermediate scales, we used the C¹⁸O (3–2) observations of NGC 1333 taken with HARPS at JCMT (Curtis et al. 2010). The angular resolution of these observations is 17.7" and it has a noise level of 0.18 K (Curtis et al. 2010).

The integrated intensity maps are calculated between 5.5 and 10 km s⁻<fnref n="1" />, which covers all the emission seen. The top right panel of Fig. 1 shows the integrated intensity map.

To trace small scales, we selected HNC and H¹³CO⁺ (1–0) because they are the brightest pair of molecular lines tracing the cloud, and have little contamination from outflows or shocks. They both cover a large fraction of the map and do not display resolved hyperfine components. The details on the image processing for these data are presented in Appendix A.

# 3. Analysis

## 3.1. Intensity power spectrum

The power spectrum of the integrated intensities is calculated using the Spatial Power Spectrum function, PS₂, in Turbstat (Koch et al. 2019) on the integrated intensity maps of the different transitions used (see Table A.1). We performed a beam correction on the power spectrum available on Turbstat, and in addition we only used spatial frequencies down to 3× the beam size, which allowed us to avoid the regime where the beam affects the power spectrum measurements. We also apodized the images with the CosineBell kernel, which is available in Turbstat, to reduce the ringing effects due to the edge of the images<fnref n="1" />.

In the case of the interferometric observations, we padded the region not covered by the mosaic to reduce the effects from the map coverage on the power spectrum. The maps are filled with Gaussian noise outside the mosaic with a sigma level of the value listed in Table A.1. In this analysis we use spatial frequency, k = 1/λ, where λ is the spatial scale probed.

## 3.2. Stitching different scales and power-law fits

The individual PS₂ from each molecular line has a different absolute amplitude, due to the difference in relative abundances and volumes traced by each molecular transition. However, the power-law exponent is inter-comparable. Therefore, we applied a normalisation coefficient to the PS₂ derived for each map. We normalized the C¹⁸O power spectrum, PS₂, C¹⁸O(k), to 1 at a spatial frequency of 10⁻⁰.⁵ pc⁻<fnref n="1" />, while the normalization parameter, fₘol, is added to the spatial power spectrum, log PS₂, mol + fₘol,

of each molecular transition as a free parameter, where mol is ¹³CO, H¹³CO⁺, and HNC, respectively. In practice, we generate two different data sets for analysis: D₁ = {¹³CO, C¹⁸O}, and H¹³CO⁺} and D₂ = {¹³CO, C¹⁸O}, and HNC. We fit a single power law to each of the combined data sets as

$$ \log PS_{2,i}(k) = A_i - \beta_i \log k, \quad (1) $$

where k is the spatial frequency in units of pc⁻<fnref n="1" />, βᵢ is the power spectrum power-law index, and Aᵢ sets the amplitude of the relation fitted to the dataset Dᵢ. The power laws and all the normalization parameters are fitted simultaneously using EMCEE (Foreman-Mackey et al. 2013) (see Appendix C for details). The results are shown in Fig. 2, while the fitted parameters are listed in Table 1.

# 4. Discussion

The physical interpretation of the intensity power spectrum relies on the assumption that the emission is optically thin, and therefore it is the power spectrum of the column density. Then the power spectrum is divided into three different regimes (e.g., Federrath et al. 2010). At large scales the transition between a plateau and a power law is associated with the injection scale of turbulence (usually expected to be similar to the cloud’s physical scale, McKee & Ostriker 2007). At intermediate scales, the power law is associated with the inertial range, which is the range of scales that characterize self-similarity in the flow. However, features at smaller scales could appear due to the influence of outflows or other feedback mechanisms (Nakamura & Li 2005; Padoan et al. 2009; Boyden et al. 2016; Xu et al. 2020). At small scales, the change in the slope of the power law is associated with the dissipation scale in the turbulence cascade. Our data appears to be beyond the range of the injection scales, but in what follows we discuss the other two features of the power spectra.

## 4.1. Exponent of the power spectra, β

We obtain a power spectrum with an average exponent of β = 2.9 ± 0.1, with the values of 2.92 ± 0.11 and 2.94 ± 0.11 for H¹³CO⁺ and HNC, respectively. This value is in agreement with the exponent found by Padoan et al. (2009) on the same region using ¹³CO (1–0) (see also Padoan et al. 2006), and it is also in rough agreement with the exponents found in H1, ¹²CO (1–0), ¹³CO (1–0), and dust extinction toward the Perseus molecular cloud (Pingel et al. 2018). However, Pingel et al. (2018) matches the pixel to the beam size and does not apply a beam correction in the power spectrum. This does not have an effect at large scales; however, it is important at scales closer to the beam and complicates a detailed comparison with the indexes obtained by Pingel et al. (2018). In other molecular clouds, studies with single-dish observations also recover an exponent comparable to the value determined here (Stutzki et al. 1998; Sun et al. 2006; Xu et al. 2020). Similar efforts focused on the power spectrum of the intensity taking advantage of interferometric mosaics combined with single-dish observations to recover all scales. In L1551, Swift & Welch (2008) found an exponent of 2.8±0.1 using C¹⁸O (1–0); while in different subregions of Orion A, Feddersen et al. (2019) found exponents between 3 and 4 using ¹²CO and ¹³CO (1–0), while slopes between 2.2 and 3.4 are determined in the case of C¹⁸O (1–0), see also Beuther et al. (2022) for analysis in a more distant high-mass region. In these cases, the most optically thin transition used is C¹⁸O (1–0),

<footnote n="1">¹ https://turbustat.readthedocs.io/en/latest/tutorials/applying_apodizing_functions.html</footnote>
<figure><img src="image_3.png" /><figcaption>**Fig. 1.** Integrated intensity map of the different tracers used to observe the NGC1333 region. Top left: Large scale traced with $^{13}$CO (1–0). The footprint of the C$^{18}$O (3–2) map is marked by the dashed box. Top right: Intermediate scale traced with C$^{18}$O (3–2). The footprint of the interferometric data is marked by the dashed box. Bottom left and right: Highest angular resolution maps of H$^{13}$CO$^+$ and HNC (1–0), respectively. The beam size and scale bar are shown in the top right and bottom right corners, respectively.</figcaption></figure>

which still should become optically thick close to dense cores or embedded YSOs, and gives exponent values comparable to those found in this work. We note that the wide range of exponents reported in the literature complicates our comparison. This variation could suggest that the power-law exponent depends on the local environment, or it may indicate that large uncertainties in previous studies, possibly due to optical depth issues, are affecting the results.

Our derived slopes match those derived by Miville-Deschênes et al. (2016) from their combined spectrum ranging

from scales between 0.01 and 50 pc (2.9 ± 0.1). However, the region considered in Miville-Deschênes et al. (2016) is a Galactic cirrus, with densities that are orders of magnitude lower than those characteristic of the gas in NGC 1333. Their 2.9 ± 0.1 exponent is comparable to that found in the infrared emission from cirrus observed with the Infrared Astronomical Satellite (IRAS) at 100 µm and Herschel at 250, 350, and 500 µm (Gautier et al. 1992; Miville-Deschênes et al. 2010). However, other Galactic cirrus clouds present different exponents in the power spectra derived from neutral atomic hydrogen (Hα), with
<table><caption>**Table 1.** Results of the power-law fits.</caption><thead><tr><th>Parameter</th><th>Value</th></tr></thead><tbody><tr><td>A<sub>H<sup>13</sup>CO+</sub></td><td>1.42<sup>+0.10</sup><sub>-0.10</sub></td></tr><tr><td>β<sub>H<sup>13</sup>CO+</sub></td><td>2.92<sup>+0.11</sup><sub>-0.11</sub></td></tr><tr><td>A<sub>HNC</sub></td><td>1.43<sup>+0.10</sup><sub>-0.10</sub></td></tr><tr><td>β<sub>HNC</sub></td><td>2.94<sup>+0.11</sup><sub>-0.11</sub></td></tr><tr><td>f<sub><sup>13</sup>CO</sub></td><td>2.51<sup>+0.08</sup><sub>-0.08</sub></td></tr><tr><td>f<sub>H<sup>13</sup>CO+</sub></td><td>-0.25<sup>+0.11</sup><sub>-0.11</sub></td></tr><tr><td>f<sub>HNC</sub></td><td>-1.50<sup>+0.10</sup><sub>-0.11</sub></td></tr></tbody></table>

exponents between 2.59 and 3.7 on the largest scales (e.g., Pingel et al. 2013; Martin et al. 2015). Similarly, the power spectra of the dust thermal emission observed with Herschel toward the Large Magellanic cloud show exponents between 1.0 and 2.43 (Colman et al. 2022). This variability is found in multiple studies across tracers (Szotkowski et al. 2019; Koch et al. 2020), implying that the power-law exponent at large scales is not universal but is related to local physical conditions.

Numerical simulations of hydrodynamic (HD) and magnetohydrodynamic (MHD) turbulence are commonly used to link the power spectrum exponent to physical quantities, such as the sonic Mach number ($M_s$), the Alfvén Mach number ($M_A$), and the mixture between solenoidal and compressive turbulent modes (e.g., Federrath & Klessen 2013; Burkhart et al. 2015a). Kim & Ryu (2005) used numerical simulations of turbulence in a compressible and isothermal medium to show that trans-sonic turbulence ($M_s \approx 1$) produces power spectra with exponents around 1.73, close to the expected values from Kolmogorov turbulence, and shallower for increasing $M_s$. A power spectrum slope of $\approx 3.0$ corresponds to their results for $M_s \approx 7$.

Dedicated studies of MHD turbulence in self-gravitating giant molecular clouds show power spectrum exponents between 1.5 and 2.5 for simulations with $0.5 < M_s < 20$ and $0.7 < M_A < 2.0$ (Burkhart et al. 2015a). These slopes are hard to reconcile with our observational results and strongly suggest the effect of additional physics. Boyden et al. (2016) include the effect of stellar winds and magnetic fields in the MHD simulations aimed at reproducing the physical conditions in the Perseus molecular cloud and produce synthetic $^{12}$CO (1–0) emission observations. However, the slopes of their resulting power spectra are lower than in our results, falling in the range between 2.64 and 2.77. Boyden et al. (2018) showed that larger variations in temperature and abundance tend to flatten the power spectrum slope, and synthetic observations of MHD turbulence produced using full astrochemical networks yield slopes of $\approx 2.8$.

Optical depth has a strong effect on the derived slope, such that the power spectrum slope of optically thick turbulent gas limits to 3.0 (Lazarian & Pogosyan 2004; Boyden et al. 2018). However, this is not the case in our observations, where the $C^{18}O$ and the higher-density tracers do not show evidence of being optically thick. Given the large degeneracy between the physical parameters that can result in our observed slope, we focus on other aspects of our observations, namely the similarity between the results for neutral and ionized species and the absence of a cut-off frequency in the power spectra.

Figures 2 and 3 show that there is a feature at $log k \approx 0.9$ pc$^{-1}$ ($\approx 150$ mpc) there is a feature in the power spectra. However, this feature is not in the $C^{18}O$ power spectrum, which covers

the same scales and with more independent samples. Therefore, this feature could be related to the shape of the interferometric mosaic since the 90" scale (corresponding to $log k \approx 0.9$ pc$^{-1}$) is comparable to the width of the mosaic covering the southeast filamentary structure.

## 4.2. Lack of difference between ions and neutrals

The neutral and ionized species are expected to decouple at certain scales due to ambipolar diffusion (Mestel & Spitzer 1956; Kulsrud & Pearce 1969; Mouschovias et al. 2011). This decoupling should introduce differences in the emission power spectra between neutral and ionized tracers (e.g., Houde et al. 2000).

The comparison of the complete power spectra between ions and neutrals (see Fig. 2) shows that they present the same power law and that neither has a turnover. This suggests that ions and neutrals do not present a substantial difference down to 20 mpc. Although both spectra agree well within $1\sigma$, there is a systematic difference between the power spectra. Figure 3 shows this possible offset appears between the power spectra of HNC and $H^{13}CO^+$ for the highest angular resolution data. This is opposite to the prediction from ambipolar diffusion theories (e.g., Li & Houde 2008). Similarly, we cannot rule out that this discrepancy in the power spectra is related to systematic differences in the abundance and/or excitation conditions at higher angular scales (Gaches et al. 2015; Pineda et al. 2022; Tritsis et al. 2023).

## 4.3. MHD wave cutoff

In a molecular cloud, where the ionization level is low, Alfvén waves cannot propagate when the collision frequency between ions and neutrals is comparable to or smaller than the MHD wave frequency (Kulsrud & Pearce 1969; Mouschovias et al. 2011; Hennebelle & Lebreuilly 2023). Therefore, there is a critical length for wave propagation beyond which waves do not transmit. The critical length for wave propagation (Stahler & Palla 2005; Mouschovias et al. 2011) is obtained from

$$ \lambda_{\min} = \sqrt{\frac{\pi}{4\mu m_H n(H_2)}} \frac{B_o}{n(H_2) X(e) \langle \sigma v \rangle}, \quad (2) $$

where $B_o$ is the unperturbed magnetic field and $\langle \sigma v \rangle$ is the rate coefficient for elastic collisions. The rate coefficient term is approximated by the Langevin term,

$$ \langle \sigma v \rangle = 1.69 \times 10^{-9} \text{ cm}^3 \text{ s}^{-1}, \quad (3) $$

for $HCO_+ - H_2$ collisions (McDaniel & Mason 1973). We estimate the magnetic field strength, $B_o$, using the relation of Crutcher et al. (2010) (see also Myers & Basu 2021; Pattle et al. 2023), a typical volume density of $10^{3.6} \text{ cm}^{-3}$, and the median value of $X(e) = 10^{-6.5}$ reported by Pineda et al. (2024). With these values, we obtain an MHD wave cutoff scale of 0.034 pc (7000 au). Therefore, we predict that if the turbulence is MHD in nature, then a dissipation scale should be observed at $\approx 0.034$ pc (or $\approx 37''$ at the distance of Perseus).

Previously, this scale had been proposed as a possible origin for the transition to coherence between the supersonic cloud and subsonic cores ($\sim 0.1$ pc; Goodman et al. 1998; Caselli et al. 2002; Pineda et al. 2010). Recent observations of dense gas with different tracers (e.g., $N_2H^+$ and $NH_3$) and angular resolutions show that the filaments present at the southern end of NGC 1333 present subsonic levels of turbulence at scales of $\approx 40''$ (Friesen et al. 2017; Hacar et al. 2017; Dhabal et al. 2019;
<figure><img src="image_4.png" /><figcaption>**Fig. 2.** Combined power spectrum of the three different datasets. The left and right panels respectively show the power spectrum with ion and neutrals at the smallest scales. The estimated MHD cutoff scale, 34 mpc, is marked by a vertical line. The fitted power law is shown by the dotted line, with exponents of $\beta = 2.92 \pm 0.11$ and $2.94 \pm 0.11$ for the left and right panels, respectively.</figcaption></figure>

<figure><img src="image_5.png" /><figcaption>**Fig. 3.** Comparison of power spectra of ions and neutrals with our highest angular resolution data. The spectra are normalized as $PS_2(k) = 1$ at the spatial frequency of $10 \text{ pc}^{-1}$. The ions have a subtle but systematic excess of power at the smallest scales when compared to the neutrals. This is opposite to the expectations from ambipolar diffusion.</figcaption></figure>

Sokolov et al. 2020), which corresponds to 0.034 pc or 7000 au at the distance of Perseus.

Finally, an improved estimate of the density or magnetic fields strength are needed to improve the obtained constraints. We used the volume density map of the region (Pineda et al. 2024) to determine a mean density $10^{3.6 \pm 0.2} \text{ cm}^{-3}$, and with this uncertainty we derived an uncertainty on the MHD wave cutoff of $\pm 13 \text{ mpc}$. This would push the $\lambda_{\min}$ right below the scales sampled by these observations. Therefore, higher angular resolution and tighter constraints on the density and magnetic field strength in the region are required to make more progress in this topic.

### 4.4. Interpretation of single power-law spectrum

The combination of previous evidence for ions presenting a higher level of turbulence than neutrals at small scales suggests that a more exotic physical process is at play. Pineda et al. (2021) suggested that MHD waves could penetrate within dense cores and perturb the magnetic field lines. As a result, this process would inject kinetic energy into the ions at smaller scales, and therefore it could remove the scale for dissipation and possible differences between ions and neutrals.

Another possibility proposed by Hennebelle & Lebreuilly (2023) involves the effect of the dust grain inertia on the transfer of energy at smaller scales. This mechanism involves the interaction between dust particles and gas at small scales, which could inject energy in the denser regions remove the scale for dissipation.

On the other hand, different two-fluid simulations have studied the nature of the turbulence of both ions and neutrals (Oishi & Mac Low 2006; Burkhart et al. 2015b; Hu et al. 2024); however, these results show that the energy can be transported across the ambipolar diffusion (AD) scale. These results suggest that our understanding of the AD process in more realistic conditions is incomplete.

Unfortunately, all these possible explanations have not provided synthetic observations to better compare with the different observations, leaving the door open for improved comparisons. This includes possible effects due to radiative transfer and/or chemistry that must be taken into account to make the comparison.

## 5. Conclusions

We studied the molecular cloud structure across different densities and scales by combining three different data sets. We derived the combined power spectrum covering more than two orders of magnitude (from $\approx 3$ pc to 20 mpc) in linear scale, with the smallest scales probed by two tracers: $H^{13}CO^+$ and HNC. The combined power spectrum fitting across all scales is well fit with a $PS_2(k) \propto k^{-\beta}$ function, where the exponents are $2.92 \pm 0.11$.
and 2.94 ± 0.11 for H$^{13}$CO$^{+}$ and HNC, respectively. The power spectrum shows no evidence of a turnover (dissipation) down to 20 mpc.

The power spectra of the ions and neutrals are not substantially different. However, a systematic offset (extra) power is present in the ions at small scales when compared to the neutrals, which is opposite to the expectations of ambipolar diffusion models. These suggestive results, combined with the previous observations showing broader line widths for ions compared to neutrals (Pineda et al. 2021), imply that, contrary to the classical picture of turbulence dissipation, our understanding of AD processes is incomplete. More theoretical and observational work is needed (e.g., Hennebelle & Lebreuilly 2023; Hu et al. 2024) to explain these results and provide a general framework for the turbulence of ion and neutral, and to further confirm these differences.

## Data availability

The code to reproduce the results of the paper are hosted at github.com/jpinedaf/NGC1333_NOEMA_turbulence

**Acknowledgements.** Part of this work was supported by the Max-Planck Society. We kindly thank the anonymous referee for the comments that helped improve the manuscript. EWK acknowledges support from the Smithsonian Institution as a Submillimeter Array (SMA) Fellow and the Natural Sciences and Engineering Research Council of Canada. DMS is supported by an NSF Astronomy and Astrophysics Postdoctoral Fellowship under award AST-2102405. JEP thanks M. A. Miville-Deschênes and F. Boulanger for stimulating discussion when starting the project, and P. Hennebelle, A. Traficante, and S. Molinari for valuable discussions. The authors acknowledge Interstellar Institute's program With Two Eyes and the Core2disk-III residential program of Institut Pascal at Université Paris-Saclay, with the support of the program "Investissements d'Avenir" ANR-11-INDEX-0003-01. This work made use of Astropy (http://www.astropy.org): a community-developed core Python package and an ecosystem of tools and resources for astronomy (Astropy Collaboration 2013, 2018, 2022). This research has made use of NASA's Astrophysics Data System Bibliographic Services.

## References

* Astropy Collaboration (Robitaille, T. P., et al.) 2013, A&A, 558, A33

* Astropy Collaboration (Price-Whelan, A. M., et al.) 2018, AJ, 156, 123

* Astropy Collaboration (Price-Whelan, A. M., et al.) 2022, ApJ, 935, 167

* Beuther, H., Wyrowski, F., Menten, K. M., et al. 2022, A&A, 665, A63

* Boyden, R. D., Koch, E. W., Rosolowsky, E. W., & Offner, S. S. R. 2016, ApJ, 833, 233

* Boyden, R. D., Offner, S. S. R., Koch, E. W., & Rosolowsky, E. W. 2018, ApJ, 860, 157

* Burkhart, B., Collins, D. C., & Lazarian, A. 2015a, ApJ, 808, 48

* Burkhart, B., Lazarian, A., Balsara, D., Meyer, C., & Cho, J. 2015b, ApJ, 805, 118

* Caselli, P., Benson, P. J., Myers, P. C., & Tafalla, M. 2002, ApJ, 572, 238

* Colman, T., Robitaille, J.-F., Hennebelle, P., et al. 2022, MNRAS, 514, 3670

* Crutcher, R. M., Wandelt, B., Heiles, C., Falgarone, E., & Troland, T. H. 2010, ApJ, 725, 466

* Curtis, E. I., Richer, J. S., & Buckle, J. V. 2010, MNRAS, 401, 455

* Dhabal, A., Mundy, L. G., Chen, C.-Y., Teuben, P., & Storm, S. 2019, ApJ, 876, 108

* Elmegreen, B. G., & Scalo, J. 2004, ARA&A, 42, 211

* Feddersen, J. R., Arce, H. G., Kong, S., Ossenkopf-Okada, V., & Carpenter, J. M. 2019, ApJ, 875, 162

* Federrath, C., & Klessen, R. S. 2013, ApJ, 763, 51

* Federrath, C., Roman-Duval, J., Klessen, R. S., Schmidt, W., & Mac Low, M. M. 2010, A&A, 512, A81

* Foreman-Mackey, D., Hogg, D. W., Lang, D., & Goodman, J. 2013, PASP, 125, 306

* Friesen, R. K., Pineda, J. E., Rosolowsky, E., et al. 2017, ApJ, 843, 63

* Gaches, B. A. I., Offner, S. S. R., Rosolowsky, E. W., & Bisbas, T. G. 2015, ApJ, 799, 235

* Gautier, T. N. I., Boulanger, F., Perault, M., & Puget, J. L. 1992, AJ, 103, 1313

* Goodman, A. A., Barranco, J. A., Wilner, D. J., & Heyer, M. H. 1998, ApJ, 504, 223

* Hacar, A., Tafalla, M., & Alves, J. 2017, A&A, 606, A123

* Hennebelle, P., & Lebreuilly, U. 2023, A&A, 674, A149

* Heyer, M. H., & Brunt, C. M. 2004, ApJ, 615, L45

* Heyer, M., & Dame, T. M. 2015, ARA&A, 53, 583

* Heyer, M., Krawczyk, C., Duval, J., & Jackson, J. M. 2009, ApJ, 699, 1092

* Houde, M., Bastien, P., Peng, R., Phillips, T. G., & Yoshida, H. 2000, ApJ, 536, 857

* Hu, Y., Xu, S., Arzamasskiy, L., Stone, J. M., & Lazarian, A. 2024, MNRAS, 527, 3945

* Kim, J., & Ryu, D. 2005, ApJ, 630, L45

* Koch, E. W., Rosolowsky, E. W., Boyden, R. D., et al. 2019, AJ, 158, 1

* Koch, E. W., Chiang, I.-D., Utomo, D., et al. 2020, MNRAS, 492, 2663

* Kulsrud, R., & Pearce, W. P. 1969, ApJ, 156, 445

* Larson, R. B. 1981, MNRAS, 194, 809

* Lazarian, A., & Pogosyan, D. 2004, ApJ, 616, 943

* Li, H.-B., & Houde, M. 2008, ApJ, 677, 1151

* Martin, P. G., Blagrave, K. P. M., Lockman, F. J., et al. 2015, ApJ, 809, 153

* McDaniel, E. W., & Mason, E. A. 1973, The Mobility and Diffusion of Ions in Gases, Wiley Series in Plasma Physics (John Wiley & Sons)

* McKee, C. F., & Ostriker, E. C. 2007, ARA&A, 45, 565

* Mestel, L., & Spitzer, Jr., L. 1956, MNRAS, 116, 503

* Miville-Deschênes, M. A., Martin, P. G., Abergel, A., et al. 2010, A&A, 518, L104

* Miville-Deschênes, M. A., Duc, P. A., Marleau, F., et al. 2016, A&A, 593, A4

* Mouschovias, T. C., Ciolek, G. E., & Morton, S. A. 2011, MNRAS, 415, 1751

* Myers, P. C., & Basu, S. 2021, ApJ, 917, 35

* Nakamura, F., & Li, Z.-Y. 2005, ApJ, 631, 411

* Oishi, J. S., & Mac Low, M.-M. 2006, ApJ, 638, 281

* Padoan, P., Juvela, M., Kritsuk, A., & Norman, M. L. 2006, ApJ, 653, L125

* Padoan, P., Juvela, M., Kritsuk, A., & Norman, M. L. 2009, ApJ, 707, L153

* Pattle, K., Fissel, L., Tahani, M., Liu, T., & Ntormousi, E. 2023, ASP Conf. Ser., 534, 193

* Pineda, J. E., Goodman, A. A., Arce, H. G., et al. 2010, ApJ, 712, L116

* Pineda, J. E., Schmiedeke, A., Caselli, P., et al. 2021, ApJ, 912, 7

* Pineda, J. E., Harju, J., Caselli, P., et al. 2022, AJ, 163, 294

* Pineda, J. E., Arzoumanian, D., Andre, P., et al. 2023, ASP Conf. Ser., 534, 233

* Pineda, J. E., Sipilä, O., Segura-Cox, D. M., et al. 2024, A&A, 686, A162

* Pingel, N. M., Stanimirović, S., Peek, J. E. G., et al. 2013, ApJ, 779, 36

* Pingel, N. M., Lee, M.-Y., Burkhart, B., & Stanimirović, S. 2018, ApJ, 856, 136

* Ridge, N. A., Di Francesco, J., Kirk, H., et al. 2006, AJ, 131, 2921

* Roman-Duval, J., Federrath, C., Brunt, C., et al. 2011, ApJ, 740, 120

* Sokolov, V., Pineda, J. E., Buchner, J., & Caselli, P. 2020, ApJ, 892, L32

* Stahler, S. W., & Palla, F. 2005, The Formation of Stars (Wiley-VCH)

* Steer, D. G., Dewdney, P. E., & Ito, M. R. 1984, A&A, 137, 159

* Stutzki, J., Bensch, F., Heithausen, A., Ossenkopf, V., & Zielinsky, M. 1998, A&A, 336, 697

* Sun, K., Kramer, C., Ossenkopf, V., et al. 2006, A&A, 451, 539

* Swift, J. J., & Welch, W. J. 2008, ApJS, 174, 202

* Szotkowski, S., Yoder, D., Stanimirović, S., et al. 2019, ApJ, 887, 111

* Tritsis, A., Basu, S., & Federrath, C. 2023, MNRAS, 521, 5087

* Virtanen, P., Gommers, R., Oliphant, T. E., et al. 2020, Nat. Methods, 17, 261

* Xu, D., Offner, S. S. R., Gutermuth, R., & Oort, C. V. 2020, ApJ, 890, 64

<footnote n="1">1. Max-Planck-Institut für extraterrestrische Physik, Giessenbacherstrasse 1, D-85748 Garching, Germany</footnote>

<footnote n="2">2. Istituto di Astrofisica e Planetologia Spaziali (IAPS), INAF, Via Fosso del Cavaliere 100, 00133 Roma, Italy</footnote>

<footnote n="3">3. Department of Astronomy, The University of Texas at Austin, 2500 Speedway, Austin, TX 78712, USA</footnote>

<footnote n="4">4. Center for Astrophysics | Harvard & Smithsonian, 60 Garden St, Cambridge, MA 02138, USA</footnote>

<footnote n="5">5. Institut de Radioastronomie Millimétrique (IRAM), 300 rue de la Piscine, F-38406 Saint-Martin d'Hères, France</footnote>

<footnote n="6">6. Department of Astronomy, University of Virginia, Charlottesville, VA 22904, USA</footnote>

<footnote n="7">7. SKA Observatory, Jodrell Bank, Lower Withington, Macclesfield SK11 9FT, UK</footnote>

<footnote n="8">8. IPAG, Université Grenoble Alpes, CNRS, F-38000 Grenoble, France</footnote>

<footnote n="9">9. Green Bank Observatory, PO Box 2, Green Bank, WV 24944, USA</footnote>

<footnote n="10">10. Department for Physics, Engineering Physics and Astrophysics, Queen's University, Kingston, ON K7L 3N6, Canada</footnote>
# Appendix A: Processing of IRAM data

The HNC (1–0) and H$^{13}$CO$^{+}$ (1–0) lines were observed with IRAM 30-m and NOrthern Extended Millimeter Array (NOEMA) and the details are described in Sec. A.1 and A.2, respectively. The combination procedure is described in Sec. A.3.

## A.1. IRAM 30m telescope

The observations were carried out with the IRAM 30m telescope at Pico Veleta (Spain) on 2021 November 9, 10, and 11; and 2022 February 19, 20, under project 091-21. The EMIR E090 receiver and the FTS50 backend were employed. We used two spectral setups to cover the H$^{13}$CO$^{+}$ (1–0) and HNC (1–0) lines at 86.7 and 90.6 GHz (see Table A.1). We mapped a region of ≈150′′×150′′ with the On-the-Fly mapping technique, and using position switching. Data reduction was performed using the CLASS program of the GILDAS package<fnref n="2" />. The beam efficiency, $B_{eff}$, is obtained using the Ruze formula (available in CLASS), and it is used to convert the observations into main beam temperatures, $T_{mb}$.

## A.2. NOEMA interferometer

The observations were carried out with the IRAM NOEMA interferometer within the S21AD program using the Band 1 receiver. The observations were carried out on 2021 July 18, 19, and 21; August 10, 14, 15, 19, 22, and 29; and September 1 in the D configuration. We observed a total of 96 pointings, which were separated on four different scheduling blocks. The mosaic's center is located at $\alpha_{J2000}=03^h29^m10.2^s$, $\delta_{J2000}=31^\circ13'49.4''$. We use the PolyFix correlator with a LO frequency of 82.505 GHz and an instantaneous bandwidth of 31 GHz spread over two sidebands (upper and lower) and two polarisations. The centers of the two 7.744 GHz wide sidebands are separated by 15.488 GHz. Each sideband is composed of two adjacent basebands of ~3.9 GHz width (inner and outer basebands). In total, there are thus eight basebands, which are fed into the correlator. The spectral resolution is 2 MHz throughout the 15.488 GHz effective bandwidth per polarization. Additionally, a total of 112 high-resolution chunks are placed, each with a width of 64 MHz and a fixed spectral resolution of 62.5 kHz. Both polarizations (H and V) are covered with the same spectral setup, and therefore the high-resolution chunks provide 66 dual polarisation spectral windows. The high spectral resolution windows used here are listed in Table A.1.

<table><caption>**Table A.1.** Spectral information for each of the spectral lines analyzed.</caption><thead><tr><th>Transition</th><th>Rest Freq.<br>(MHz)</th><th>Beam (PA)<br>("x")</th><th>rms</th><th>rms unit</th><th>V<sub>min</sub><br>(km s<sup>-1</sup>)</th><th>V<sub>max</sub><br>(km s<sup>-1</sup>)</th><th>Padding noise<br>(mJy beam<sup>-1</sup> km s<sup>-1</sup>)</th></tr></thead><tbody><tr><td><sup>13</sup>CO (1–0)</td><td>110201.35430</td><td>46</td><td>0.12</td><td>K</td><td>-1.0</td><td>11.5</td><td>—</td></tr><tr><td>C<sup>18</sup>O (3–2)</td><td>329330.55250</td><td>17.7</td><td>0.18</td><td>K</td><td>5.5</td><td>10</td><td>—</td></tr><tr><td>HNC (1–0)</td><td>90663.568</td><td>4.9×4.7(−38°)</td><td>23</td><td>mJy beam<sup>-1</sup></td><td>3.5</td><td>11</td><td>30</td></tr><tr><td>H<sup>13</sup>CO<sup>+</sup> (1–0)</td><td>86754.2884</td><td>5.0×4.9(−43°)</td><td>15</td><td>mJy beam<sup>-1</sup></td><td>5.3</td><td>10</td><td>20</td></tr></tbody></table>

## A.3. Image combination

The original IRAM 30m data is resampled to match the spectral setup of the NOEMA observations. We use the task uvshort to generate the pseudo-visibilities from the 30-m data for each NOEMA pointing. The imaging is done with natural weighting, a support mask, and using the SDI deconvolution algorithm (Steer et al. 1984).

The noise level of the combined images is reported in Table A.1. The integrated intensity maps are calculated using the velocity ranges listed in Table A.1, this velocity range covers all the emission seen in both molecules. The bottom left and right panels of Fig. 1 show the integrated intensity maps for H$^{13}$CO$^{+}$ and HNC.

# Appendix B: Median spectra

We show the median spectra for across the four different data cubes used in Fig. B.1. They are determined as the median value of all the values in the cubes.

<figure><img src="image_6.png" /><figcaption>**Fig. B.1.** Typical spectra for the different tracers used calculated over the regions used in this work.</figcaption></figure>

Given the average spectra, we estimate the optical depth of the brightest line in the sample: HNC (1–0). We use the radiative transfer solution to estimate the optical depth as

$$ \tau(\text{HNC}) = -\ln\left(1.0 - \frac{T_p(\text{HNC})}{J_v(T_{ex}) - J_v(T_{cmb})}\right), \quad (\text{B.1}) $$

where $T_p(\text{HNC})$ is the observed line peak brightness, and $J_v(T)$ is the brightness temperature of a black body with temperature $T$ at frequency $v$. For HNC (1–0) we assume the same excitation

<footnote n="2">² http://www.iram.fr/IRAMFR/GILDAS</footnote>
temperature of $T_{ex}$ = 12 K as used by Pineda et al. (2024), and obtain $\tau$(HNC) = 0.23, which shows the emission is optically thin. We perform the same calculations with $T_{ex}$ = 10 K and obtain an optical depth of 0.31, confirming that the optically thin approximation is reasonable in this region.

## Appendix C: Power-law fit with EMCEE

The emcee fit is initialized starting from the result of the linear fit using the minimize function within the scipy python package (Virtanen et al. 2020). In the emcee run, we use uniform priors $U(a, b)$, which is constant between *a* and *b*. We use $U(1, 2.5)$ for $A_{H^{13}CO^+}$ and $A_{HNC}$, $U(2.0, 3.4)$ for $\beta_{H^{13}CO^+}$ and $\beta_{HNC}$, $U(2, 3)$ for $f_{H^{13}CO}$, $U(-2, -1)$ for $f_{H^{13}CO^+}$, and $U(-1, 0.5)$ for $f_{HNC}$, where the ranges are set to cover the best fit from the initial linear fit and without reaching the limits with the chains.

We used 56 random walkers, with 50 000 steps. We estimate an autocorrelation time of 80 steps, and therefore before analyzing, we discarded the first 600 steps and then thinned the samples by 40 to obtain better estimates. The corner plots for both fits are shown in Fig. C.1 and summarized in Table 1.

<figure><img src="image_7.png" /><figcaption>Fig. C.1. Corner plot for all parameters fitted. See Table 1 for a summary of the fits.</figcaption></figure>