© The Authors 2024

1 Introduction

For star formation to occur, sufficient accretion of material by the protostar must take place, and for this, part of the initial angular momentum must be redistributed away, so that it is not part of the protostar system. If it was not the case, the rotational velocity would be too large for the protostar to accrete sufficient material. The specific angular momentum (SAM) can be reduced by one order of magnitude by means of magnetic tension that carries it outward through a process known as “magnetic braking” (Allen et al. 2003; Galli et al. 2006). From a theoretical standpoint, it has been shown that in ideal magnetohydrodynamic (MHD), spherical collapse, and nonturbulent conditions, magnetic braking catastrophically redistributes most of the angular momentum during the collapse, thereby modifying the kinematics of the gas in the inner envelopes and potentially preventing the formation of centrifugally supported disks if no diffusive effects or misalignment of the magnetic field with the rotational axis are considered (e.g., Li et al. 2014). However, more refined magnetized models have been developed that allow for angular momentum redistribution but include nonideal MHD to avoid the so-called catastrophe (e.g., Masson et al. 2016), or misalignment of the rotation axis of the system and the global magnetic field to reduce the magnetic braking efficiency (Joos et al. 2012; Gray et al. 2018). The exact role of magnetic fields and magnetic braking in setting the properties of protostellar disks is still debated in the literature of protostellar models. For example, Seifried et al. (2012) have shown that the collapse of turbulent cloud cores makes it possible to avoid the catastrophe even in the ideal MHD limit, while Santos-Lima et al. (2012); Kuffmeier et al. (2017) propose that disk formation is a natural consequence of the inherited gas turbulence from large scales, even in magnetized conditions. However, Wurster et al. (2019) suggest that it is the density and velocity structure and dynamical interactions between protostars that regulate the properties of disks, and not the magnetic braking. Finally, Lebreuilly et al. (2021, 2024) conclude that the inclusion of magnetic field is critical to reproducing the observed protostellar disk sizes formed consistently, even from models including large-scale turbulence.

Observations of Class 0 protostars have shown that magnetic fields are ubiquitous in the inner envelopes at scales of ~100–5000 au (Galametz et al. 2018). Moreover, indications of the coupling between the infalling material and the magnetic field have been observed, such as pinched magnetic field lines forming hourglass patterns (e.g., Hull et al. 2014; Maury et al. 2018; Redaelli et al. 2019; Liu 2021), and streamers of gas aligned with magnetic field lines (Sadavoy et al. 2018; Le Gouellec et al. 2019). It is expected that in the region where the magnetic field is being pinched there will be a decrease in angular momentum, and thus a variation in rotational velocity, and therefore a more severe magnetic braking around the equatorial plane. Phe-nomenologically, the effect of the magnetic field has also been discussed in relation to the sizes of disks formed as a result of the angular momentum of the gas (Keplerian protostellar disks). While few observations are able to measure the gas kinematics at disk scales, millimeter dust continuum emission was used to estimate the size of the dusty disks and show that while most Class 0 protostars harbor a disk, the vast majority are only found at radii smaller than 60 au, which is difficult to reconcile with purely hydrodynamical models (Maury et al. 2019; Sheehan et al. 2022). Moreover, tentative correlations were found between the misalignment of the rotation axis of the system and the magnetic field direction in the inner envelopes of protostars and the kinetic energy of the gas measured from observed velocity gradients (Galametz et al. 2020; Huang et al. 2024), suggesting a possible causality between these two quantities due to a change in magnetic braking efficiency.

Establishing a connection between observation and modeling is necessary to understand what mechanisms govern the formation of stars and disks. In addition to the observation of thermal dust and polarized emission at (sub)millimeter wavelengths, the study of molecular line emission can give us information about the gas kinematics on inner envelope and disk scales (see Belloche 2013, for a review of observational techniques used to measure gas angular momentum). It is important to investigate how variation in the initial conditions, for example the magnetic field, affects the rotational properties of the inner core -specifically, how magnetic field variations in terms of intensity affect the real angular momentum of the system, but also how this fact is transferred to the projected angular momentum in the plane of the sky and, finally, how these possible features are manifested in observations. The influence of magnetic braking in models of protostellar evolution has been investigated in the past by paying attention to the misalignment between the rotation axis and the magnetic field, to the ratio between thermal and gravitational energy, and to different conditions affecting the MHD resistivities (for a review, see Hirano et al. 2020). However, to our knowledge, there are no studies in the literature that focus on the ability of state-of-the-art observations to capture the kinematic signatures of angular momentum spatial distribution and evolution in different magnetic braking conditions. In this work, we have analyzed the kinematics of the gas in nonideal MHD models of low-mass Class 0 protostar formation, carried out with RAMSES (Teyssier 2002). First, we analyzed the gas kinematics directly from the models, then we investigated how the gas kinematics can be probed with observations from molecular line emission, post-processing the models with the radiative transfer (RT) tool POLARIS (Reissl et al. 2016; Brauer et al. 2017) and analyzing the resulting synthetic observations. The models and its analysis are presented in Sect. 2. We calculated the three-dimensional (3D) SAM at different envelope scales and evolutionary stages. Section 3 shows the C¹⁸O (2−1) synthetic observations. We apply state-of-the-art techniques used to derive angular momentum profiles from observational data, and show the resulting radial profiles. We discuss the impact of magnetic fields in the angular momentum transport in nonideal MHD models and the ability of observational methods to detect the braking of gas motions due to magnetic fields in Sect. 5. In Sect. 6, we present our conclusions.

2 Magnetohydrodynamic models of Class 0 protostars

2.1 Description of the models

We selected four numerical realizations modeling the collapse and evolution of 1 M_⊙ protostellar cores, which will ultimately form low-mass stars. Two of these models were published in Hennebelle et al. (2020), and two are new realizations performed specifically for the study presented here. All these numerical models use the adaptive mesh refinement code RAMSES (Teyssier 2002) to solve the MHD equations through the finite volume method (Fromang et al. 2006) (for more details on the numerical setup and physics included in the models, we refer to Hennebelle et al. 2020): Table 1 summarizes the initial conditions (Cols. 2, 3) and properties (Cols. 4–6) of the outputs analyzed in our study.

Two models, μ3ℳ0 and μ10ℳ0, have the same initial conditions, but different ratios of gravitational to magnetic energies (a.k.a., the mass-to-flux ratio, μ, expressed in units of the critical value, as it is defined in Eq. (1)),

$$\mu =\frac{(M/\text{Φ})}{{(M/\text{Φ})}_{\text{crit }}},$$(1)

where (M/Φ)_crit is the critical ratio when the magnetic energy equals the gravitational energy (for more details, see Nakano & Nakamura 1978). The models have, respectively, μ of 3.33 and 10 (see Col. 2 of Table 1). These models develop disks, the orientations and sizes of which are estimated according to the criterion of Lebreuilly et al. (2021, i.e., the azimuthal component of velocity dominating the rest of the components and considering 1×10⁻¹³ g cm⁻³ as the density threshold). This is well suited to capture the disk at this evolutionary stage, minimizing contamination from the envelope and streamers. The difference in the initial magnetic energy results in different disk radii; for example, at the most advanced stage examined here, the disks are ~22 au (μ = 3.33, stronger braking) and ~42 au (μ = 10, weaker braking), respectively. The similar initial conditions and different magnetization allow us to discuss the effect of the magnetic field on the gas velocity during the disk formation stage. The two other models, μ3ℳ1 and μ10ℳ1, are similar to μ3ℳ0 and μ10ℳ0 respectively, but contain more turbulent energy with an initial Mach of 1.

We extracted the 3D models (gas density, magnetic field lines, and the three gas velocity components) and examined them in a projection along the z-x plane that corresponds to a quasi-edge-on orientation (i.e.,with the outflow axis in the plane of sky). We made this choice because it is well suited to trace the rotational motions of the gas in synthetic observations, since their projection along the line of sight is maximum. Out of the full model box of 15 676 au, we analyzed a square field of view of side 8000 au centered on the protostar position.

From each of these four numerical simulations, we analyzed three model outputs at different evolutionary stages (except for μ10ℳ1, the execution of which did not reach the third evolutionary stage), selected by the accreted mass from the envelope into the central protostellar embryo, modeled by a sink particle. Class 0 protostars are defined as objects that have accreted less than half of their final stellar mass, and for which the envelope mass is still larger than the star+disk mass. Thus, for models considering an initial 1 M_⊙ core, we would expect that objects with stellar masses of M_star <0.3 M_⊙ are still in this early phase of vigorous accretion: we consider here three outputs spanning this stage. These models still show a significant envelope to be present, but we stress that the most evolved M_star = 0.3 M_⊙ output may be more representative of the transition Class 0/I. Finally, we note that models with larger initial turbulent energy have a slower time evolution; that is, the protostar takes longer to reach the same accreted mass than nonturbulent models. Also, models with larger magnetic energy build stellar mass faster than models with a weaker initial magnetic field. Hence, comparing models at similar stellar masses does not mean comparing similar times. Figure A.1 shows the total gas column density maps integrated along the line of sight, for the 11 outputs studied here. The protostar is located approximately at the geometric center of the maps shown.

2.2 Rotational velocity of the gas in protostellar models

In order to study the distribution and evolution of angular momentum in collapse models, we first looked at the rotational component of the gas velocity. Figures 1 and 2 present the radial rotational velocity profiles (ν_ϕ) computed in the azimuthally averaged models, along the equatorial plane. From these profiles, it is striking that the global radial evolution of ν_ϕ is mostly determined by the level of magnetization set by initial parameters: the models with stronger magnetic fields (dashed lines, μ ~ 3) show rotational velocities increasing much less steeply as the gas is collapsing toward the star-disk system, and as a result gas rotational velocities in the inner envelopes become up to an order of magnitude lower than in the less magnetized models (plain lines, μ ~ 10). We also stress that the increase in ν_ϕ from large scales to small scales is overall smooth and monotonous for the less magnetized models, while the profiles show more oscillations at radii of <1000 au in the more magnetized models (although the log-scale amplifies the variations at lower rotational velocity levels).

Finally, some variations in the low values of ν_ϕ as is shown for example by the most magnetized models, cannot be captured by the typical spectral resolutions available in observations of molecular line emission (see the two horizontal lines in Figs. 1 and 2 corresponding to the velocity resolution of 0.12 and 0.06 km s⁻¹). Going to a higher spectral resolution is prohibitive in terms of observing time, and is difficult because velocity differences of less than a few 0.01 km s⁻¹ would be embedded in the thermal broadening of the line. Thus, to perform the synthetic observations (described in Sect. 3), we adopted a spectral resolution of 0.06 km s⁻¹.

Fig. 1 Radial profiles of the rotational velocity of the gas in models without initial turbulence (i.e., μ3ℳ0 and μ10ℳ0). The radial profiles were computed from azimuthally averaged models, along the equatorial plane (disk plane).

Fig. 2 Radial profiles of the rotational velocity of the gas in models with initial turbulence (i.e., μ3ℳ1 and μ10ℳ1). The radial profiles were computed from azimuthally averaged models, along the equatorial plane (disk plane).

Fig. 3 Maps of the azimuthally averaged SAM in the models, $j_{V\phi }^{\text{model }}$. The white star and dashed blue line depict the location of the protostar and the equatorial plane, respectively. The leftmost panels, from (a) to (f), show the nonturbulent models, while the rightmost panels, from (g) to (k), present the turbulent models’ outputs.

2.3 Specific angular momentum of the gas in protostellar models

The next step of the analysis consisted of setting the ground truth by characterizing the SAM of the gas from the modeled proto-stellar envelopes. We computed the local SAM in the 3D model from Eq. (2), where j is the local SAM, $r_{\text{cyl}}^{\delta }$ is the radial distance to the rotational axis in cylindrical coordinates, and $v_{\phi }^{\delta }$ is the local gas rotational velocity component, within δ planes that contain the rotational axis. In the following, we always refer to this quantity as SAM, but we stress that the mathematically correct definition of SAM should take the mass of the particles into account, as well as all the cross products between r and v. However, observational estimates of the SAM cannot use this mass-weighted definition of the SAM, and are focused on measuring angular momentum in the equatorial plane using the velocity projected on the line of sight, where only ν_ϕ contributes. Hence, we used this quantity, $j_{{\text{V}}_{\phi }}^{\text{model }}$, as the estimator of the SAM in order to discuss whether synthetic observations of protostellar envelopes can distinguish signatures of magnetically braked gas. We provide in Appendix B further information regarding the expression used, as well as a brief discussion of the angular momentum and SAM radial profiles from the models.

$$j_{{\text{V}}_{\phi }}^{\mathrm{mod}el}=\frac{1}{n}{\displaystyle \sum _{\delta =0}^{2\pi }j},\text{\hspace{0.17em}\hspace{0.17em}with\hspace{0.17em}\hspace{0.17em}}j=r_{\text{cyl}}^{\delta }{\upsilon }_{\phi }^{\delta }$$(2)

As a second step, we built maps of the azimuthally averaged SAM: we averaged the 3D SAM along the azimuth, δ, sampling all 2π radians with steps of 0.017 rad. The azimuthally averaged maps, $j_{V_{\phi }}^{\text{model }}$, are shown in Fig. 3 for ℳ= 0 (panels from a to f) and ℳ=1 models (panels from g to k), where the grid has been rotated by 30° to show the rotation axis in the vertical direction. The values of the $j_{V_{\phi }}^{\text{model }}$ associated with rotational motions of the protostellar gas are typically a few 10⁻³ km s⁻¹ pc, with variations at envelope scales of usually less than an order of magnitude. While all models start with strictly identical rotational motions, ℳ0 model maps (μ3ℳ0 and μ10ℳ0) show that already very early on (the first snapshot at M_star < 0.1 M_⊙), the spatial distribution of the $j_{V_{\phi }}^{\text{model }}$ is quite different between the two models in the inner 4000 au that we study here. The μ3ℳ0 model (panels a–c) evolves in producing its largest values of $j_{V\phi }^{\text{model }}$ at large radii, with rotational motions largely associated with the development of the protostellar outflow, due to enforced corotation of gas out to some distance along the field line (Alfvén radius), as is seen in many MHD models (e.g., Machida & Hosokawa 2013; Marchand et al. 2020) and as was recently observed (Tabone et al. 2017; Lee et al. 2021). On the other hand, the μ10ℳ0 model (panels d–f of Fig. 3) concentrates most of its initial envelope SAM on the equatorial plane around the protostar, in the form of a flattened structure, the axis of which is relatively well aligned with the rotation axis of the envelope. Both models show the development of non-axisymmetric features despite the absence of an initial turbulent field.

Turbulent models, μ3ℳ1 and μ10ℳ1, also show very different spatial distribution of the SAM early on. The more magnetized μ3ℳ1 model (panels g–i) starts from a highly asymmetric distribution of its SAM, which is subsequently quickly redistributed outward by the development of outflow lobes at M_star< 0.2 M_⊙ (panel h), and finally redistributed outside the region examined here (radii >4000 au) when M_star reaches 0.3 M_⊙ (panel i). The less magnetized μ10ℳ1 model (panels j and k) evolves more similarly to μ10ℳ0, redistributing the initial SAM mostly within the equatorial plane, although with a less symmetric distribution than in the nonturbulent case.

As a third step, we built the radial profiles of $j_{V_{\phi }}^{\text{model }}$ along the equatorial plane, indicated by the blue line in Fig. 3. These profiles are shown in the first column of Fig. 4: panels a and e for ℳ0 and ℳ1 models, respectively. Regardless of the amount of turbulent energy and the magnetization level, the SAMs have similar values in the outer regions of the core, beyond 3000 au, in the eleven outputs of the models, as was expected from the common initial conditions. Altogether, the μ3 models show a larger decrease in the SAM toward smaller radii than the μ10 models, but while these more magnetized models really differentiate themselves from the μ10 at 1000-to-100 au in the nonturbulent case, most of the difference in the models with turbulence comes from larger scales where almost one order of magnitude decrease is observed between 3000-to-500 au scales. At smaller envelope radii, ≲500 au, the SAM of the more magnetized models has small values of ~10⁻⁴ to ~10⁻⁵ km s⁻¹ pc and remains overall constant in this range, while exhibiting large relative variations, at all evolutionary stages. In the case of the μ10ℳ0 and μ10ℳ1 models (panels a and e, solid lines), we see a much flatter profile with a weak slope that translates into a decrease in the SAM by a factor of ~5 in the range from ~4000 au to ~40 au.

As is shown in Fig. 4, the magnetization level influences the radial dependence as well as the magnitude of the SAM of the gas in the modeled envelopes. These radial profiles suggest that the more intense the magnetic field is, the less SAM is contained in the rotational motions (the dashed curves below the solid curves) and the more efficient the transport of the SAM toward the outer part of the core (the dashed curves exhibit steeper slopes than the solid ones), provided that the SAM contained in the initial core-scale rotation is not transmitted to the forming stellar embryo. This holds true even when we significantly increase the level of turbulence in the system (see Fig. 4, panel e).

To conclude, the model SAM profiles show distinguishable differences both in magnitude and in their evolution with envelope radii: weakly magnetized models produce large SAM values and a monotonous decrease toward smaller radii, while strongly magnetized models are associated with lower SAM values in the inner envelope and profiles with both large gradients of the SAM at some intermediate envelope radii and larger relative variations around the lowest values.

3 Synthetic observations of the gas kinematics with molecular tracers

In this section, we produced synthetic observations of the molecular line emission from the gas in the models, to assess the impact of source geometry and RT effects on recovering the SAM radial profiles from observations.

3.1 Synthetic observations: post-processing method

RAMSES outputs (gas density, gas pressure, or temperature, three components of the velocity field, and three components of the magnetic field) can be post-processed with POLARIS (POLArized RadIation Simulator), a 3D RT code (Reissl et al. 2016; Brauer et al. 2017) that solves the RT problem self-consistently on the basis of the Monte Carlo method. In the local thermodynamic equilibrium (LTE) approximation, the radiation of the gas is determined by the local kinetic temperature and its internal properties. The atomic or molecular level populations are dominated by particle collisions that obey a Maxwell-Boltzmann distribution, and used together with Kirchhoff’s law for thermal radiation to integrate the RT equation with isotropic ray tracing crossing the grid, and also to calculate the gas temperature. The source of photons is placed at the center of the model, where the sink particle traces the position of the protostellar object. The resulting specific intensity is then used to compute the observed molecular line intensity from the physical source, at the chosen spectral resolution. The resulting output FITS file is a cube; in other words, a series of maps, where each one corresponds to a certain channel or velocity. We placed the source at a distance of 250 pc, a typical distance to the closest star-forming regions (Lada & Lada 2003; Zucker et al. 2019), and proceeded with the RT from a central source with an accretion luminosity of 1 L_⊙, a typical value for bolometric luminosity observed toward solartype Class 0 and I protostars (Maury et al. 2011; Dunham et al. 2014). To facilitate the propagation of photons from the central source through the high-density optically thick material immediately surrounding it, and reduce the computational times, we emptied a small area of radius 4 au around the sink particle (as in Valdivia et al. 2019).

We chose to trace the gas kinematics with synthetic observations of the spectral emission from the C^l8O molecule. This molecule is a good tracer of high-column-density material because of its low abundance (the emission remains optically thin). While its low critical density (~8.4 × 10³ cm⁻³) would make it quite an abundant tracer of gas at high density, the gas-phase CO freezes out onto dust grains at low temperatures and is released only in the inner envelopes (i.e., at temperatures ≳25 K). Finally, the upper-level energy of the (2−1) transition (E_up = 15.81 K) makes the C^l8O(2−1) a good tracer of the lukewarm inner protostellar envelope, down to disk scales, as has been confirmed by observations (e.g., Gaudel et al. 2020; Tychoniec et al. 2021). We assumed a relative to H₂ abundance of 3 × 10⁻⁸ for the molecule, and LTE to calculate the level populations. The resulting spectral cubes contain 119 velocity channels of ~0.06 km s⁻¹, over a velocity range of ±~3.5 km s⁻¹, which covers the range of gas velocities in the models. We looked at the model putting the z−x plane in the plane of the sky, viewing the equatorial plane edge-on, as was done when exploring the models in Sect. 2.2. We focused our analysis on spatial scales of ≲8000 au around the central protostar; for example, the inner half of the core. This allowed us to avoid numerical effects at large radii caused by the periodic boundary conditions of the model.

3.2 C¹⁸O (2−1) synthetic observations: Velocity maps

We computed the intensity-weighted velocity maps ($M_1=\frac{{{\displaystyle \sum }}^{\text{}}{I}_i{v}_i}{{{\displaystyle \sum }}^{\text{}}{I}_i}$; a.k.a., moment 1 maps) using the immoments task in the CASA software¹.

The moment 1 maps for the nonturbulent models are shown in Fig. 5 (μ3ℳ0 in the upper row and μ10ℳ0 in the bottom row). The μ10ℳ0 model shows a prominent velocity gradient in the northeast-southwest direction, perpendicular to the rotational axis of the core, suggestive of highly organized motions in the equatorial plane. On the other hand, the μ3ℳ0 model presents more complex spatial distribution of rotational motions than μ10ℳ0. In the two earliest stages, the highest gas velocities are associated with the cavity walls of the southeast-northwest outflow, while the latest stage shows the development of highvelocity motions at small scales around the protostar, with supersonic streamer structures connecting the inner envelope to the disk scales.

As is shown in Fig. 6, the addition of an initial turbulent field produces more disorganized first moment maps. This is especially true for the μ3ℳ1 model, which shows much smaller gas velocities recovered in the first moment map, and much weaker signatures of a bipolar outflow than in the nonturbulent case. The least magnetized conditions (the μ10ℳ1 model) keep an organized velocity gradient from northeast to southwest at large scales, but the highly symmetrical equatorial features are replaced by a more compact and warped structure at high velocities around the protostar.

Fig. 4 SAM radial profiles computed through the equatorial plane for nonturbulent (upper row) and turbulent (lower row) models. The first column shows the profiles from the models, while the second to last columns show the profiles from synthetic observations of the models, using different techniques to measure the rotational velocities (see Sect. 3). More (μ = 3.33) and less (μ = 10) magnetized models are depicted with dotted and solid lines, respectively. Different evolutionary states are represented on a color scale: black, blue, and red for protostellar masses equal to ~0.06/0.07, 0.18, and 0.30 M⊙, respectively. The shaded gray area represents the region that cannot be sampled because it lies below the velocity resolution of synthetic observations (0.06 km s−1).

Fig. 5 First moment map (color scale) derived from the C18O (2−1) line for μ = 3.33, ℳ = 0 and μ = 10, ℳ = 0 models (upper and bottom rows, respectively). The dashed black line depicts the equatorial plane.

Fig. 6 First moment map (color scale) derived from the C18O (2−1) line for μ = 3.33, ℳ = 1 and μ = 10, ℳ = 1 models (upper and bottom rows, respectively). The dashed black line depicts the equatorial plane.

3.3 Specific angular momentum profiles from synthetic observations

The main goal of this study is to test whether analysis methods commonly applied to molecular line observations manage to capture the magnitude of SAM and its different trends depending on the efficiency of magnetic braking. In the following, we present our implementation of the three main approaches found in the literature, applied to our synthetic observational data-cubes.

First, we built position-velocity (PV) diagrams, computed by scanning along the equatorial plane in the data-cubes, from which we extracted the profile of maximum velocity with distance to the protostar. We applied the method showed in Seifried et al. (2016) to extract the maximum rotation velocity radial profile, considering an intensity threshold of 5σ, where the adopted σ RMS value is 0.01 mJy pixel⁻¹. This method makes the hypothesis that the gas moving at the fastest velocities is dominated by rotational motions in the equatorial plane, and is commonly used as an approximation of the Keplerian velocity in observations at disk scales. The $j_{\text{Vmax }}^{\text{obs}}$ SAM profiles in the equatorial plane obtained from this method are shown in panels b and f of Fig. 4. Comparing panels b and f to panels a and e, we note that the SAM profiles built from this method tend to overestimate the SAM of all models at large envelope radii. Moreover, because of its poor sensitivity to capturing motions at very small velocities, such as those produced in the inner envelope in the most magnetized case, this method systematically overestimates the SAM and fails to identify the decrease in rotational motions due to magnetic braking in those cases.

In the second method, the rotational velocity used to build the SAM profile was extracted directly from the intensity-weighted velocity maps presented in Sect. 3.2. This method makes the hypothesis that rotational motions are dominating the line-of-sight velocity in the equatorial plane, and thus small variations in the intensity-weighted peak velocity can be used as a direct measure of the rotational velocity to derive the SAM $j_{\text{Vpeak }}^{\text{obs }}$. Figure 4 shows these $j_{\text{Vpeak }}^{\text{obs}}$ radial profiles (panels c and g presents ℳ = 0 and ℳ = 1 models, respectively). As was expected, this approach produces observed SAM values smaller than the method using $j_{\text{Vmax }}^{\text{obs}}$. Moreover, because it is tailored to measure small velocity offsets from the material at rest, this method produces very small values for the rotational velocities, and a somewhat smaller velocity resolution than that adopted for the synthetic observations. As a result, part of the $j_{\text{Vpeak }}^{\text{obs }}$ profile falls below the detectability limit of the synthetic observations (gray zones in Fig. 4), especially at large scales where rotational velocities are small. In the quantitative analyses of profiles provided in Sect. 5, we will consider that $j_{\text{Vpeak }}^{\text{obs }}$ values lower than 0.06 [km s⁻¹] × r [pc] are undetected motions (consistent with null velocity). However, this method is able to capture the global magnitudes and trends of the models $j_{V_{\phi }}^{\text{model }}$ radial profiles better than the previous method, and models with different magnetizations are observed to occupy a different parameter space in the profiles from synthetic observations obtained from the first moment maps.

Finally, we also followed a third method that is close in spirit to the second one and that is also used in the literature (as in, e.g., Sai et al. 2023), where we minimized a Gaussian line profile to each spectrum along the equatorial plane, directly from the spectral cube, to estimate the peak velocity, from which we computed the SAM $j_{\text{vGauss}}^{\text{obs}}$. We present $j_{\text{vGauss}}^{\text{obs}}$ radial profiles for nonturbulent and turbulent models in panels d and h of Fig. 4, respectively. The profiles are very similar to the profiles obtained from first moment maps. The trends of models with different magnetization are mostly captured and rendered as occupying a different parameter space when using this method of computing the SAM.

In the following, we discuss the kinematic signatures of magnetic braking from SAM in profiles computed in the equatorial plane of the models, then we quantitatively test the ability of different observational techniques to capture the trends of SAM profiles and identify different magnetization regimes.

4 Accuracy of observational methods in capturing the specific angular momentum of gas in protostellar envelopes

4.1 Accuracy in recovering the absolute values of the specific angular momentum in protostellar envelopes

We first estimated the accuracy of the observational methods in estimating the amount of SAM in protostars by computing the distance between model SAM profiles and observed SAM profiles. Table D.1 shows the distance (d) between the models and synthetic observations for the three methods described in Sect. 3.3, and between models with different magnetization levels (see Appendix D).

On the other hand, Fig. 7 and Cols. 2–4 in Table D.1 show that the methods based on peak velocity (both moment 1 and Gaussian fitting) are the most efficient at reproducing the intrinsic SAM of the analyzed models. These methods are thus shown to be more prone to capturing radial velocity variations occurring on the envelope scale, with the distance to the model ranging from 2.1 to 36. Moreover, both methods produce radial profiles with similar magnitudes and trends. As Col. 8 shows, their difference in magnitude from the equatorial profile (measured as their distance, d) is under ~12 (where 0 is obtained for identical curves). On average, $j_{\text{Vpeak }}^{\text{obs}}$ and $j_{\text{VGauss }}^{\text{obs }}$ show distances to the models of 13.9 and 11.8, respectively, which makes it difficult to favor one method. The method based on maximum velocity systematically performs less well in capturing the model SAM profiles (as is illustrated by the high values in Col. 2 of Table D.1).

Column 5 of Table D.1 provides the distance between radial profiles stemming from models with different magnetization levels, which measure how different the radial profiles of SAM are in the models. The differences between these SAM profiles are significantly lower than the typical error of observational methods in capturing the magnitude of SAM profiles (e.g., the values in Col. 5 are on average lower than the values in Cols. 2–4) for the nonturbulent models. This shows that the observational method is most of the time not accurate enough to discern the difference in SAM magnitudes between two differently magnetized scenarios, in the nonturbulent case. For turbulent models, Table D.1 shows larger differences between models of different magnetization, with larger distance values reported in Col. 5 than the error of the observational methods reported in Cols. 2–4. Especially if using the $j_{\text{Vpeak }}^{\text{obs }}$ or $j_{\text{VGauss }}^{\text{obs}}$ methods to build the radial profiles, this suggests that in some cases observational methods could identify the difference in the magnitude of the SAM, depending on the level of magnetization. In practice, and to conclude, observational methods seem mostly unable to distinguish envelopes with different magnetization by measuring differences in the magnitude of SAM in the equatorial plane.

We note also that Col. 6 of Table D.1 shows that the more magnetized model exhibits larger differences between the SAM radial profiles estimated from the peak velocity, and the one estimated from the maximum velocity, than the less magnetized one – 60 versus 40%, respectively – regardless of the turbulent energy. Hence, in highly magnetized protostellar envelopes where magnetic braking is efficient, one expects to detect a large difference between the SAM radial profiles computed from the maximum velocity and from the peak velocity, probably because the SAM probed with the maximum gas velocity is only efficient at capturing the rotational component of velocity in environments where rotation dominates over infall, as is the case in the μ10ℳ0 and μ10ℳ1 models.

4.2 Accuracy in recovering the radial variations of the SAM

Second, we tested the accuracy of the methods in measuring the variation of the SAM at different envelope radii by adjusting power law functions to the synthetic observations profiles, including a change of power law indices (a broken power law) when a clear break could be identified by eye. Tables C.1 and C.2 present the power law parameters. Figure 7 shows the SAM radial profiles from synthetic observation and models, together with the best-fit power law: it shows that none of the observational methods of inferring the SAM manages to perfectly capture the magnitude or trend of the models’ SAM (the red, green, and blue profiles do not approach the black one). We stress that our analysis of synthetic observations suggests that while the break in power law indices could be detected, at least in the most evolved stages, no such large index could be recovered (see Table C.2), as the “true” 3D SAM index of 1.32 in μ3ℳ0 is recovered as 0.83 and 0.93 from the $j_{\text{VGauss }}^{\text{obs}}$ and $j_{\text{Vpeak }}^{\text{obs }}$, respectively, and an index of 3.96 in μ3ℳ1 is recovered as ~0.86 from $j_{\text{Vpeak }}^{\text{obs }}$ and $j_{\text{VGauss }}^{\text{obs}}$. However, our synthetic observations show that the shape of observed SAM profiles can still provide clues as to whether protostars experience strong magnetic braking. Specifically, the radial profiles of SAM built along the equatorial plane show different properties, exhibiting a rather monotonous power law decline in the weak braking scenario, and complex profiles when the magnetic braking is strong.

Fig. 7 SAM of the gas: radial profiles computed along the equatorial plane. They are computed from the gas velocity in the model’s outputs, and from the velocity recovered with synthetic observations of the C18O (2−1) line emission. Model: black circles show 3D SAM computed considering the azimuthally averaged rotational (νϕ) velocity component. Synthetic observation: Red crosses show the SAM computed from the maximum velocity. Green stars show the SAM computed from the intensity-weighted velocity (moment 1). Blue circles show the SAM computed from the Gaussian-fitted peak velocity. The shaded gray area represents the region beyond the velocity resolution. Red, green, and blue shadows show error bars on the synthetic data.

5 Discussion

5.1 Magnetic braking: Whether it can be identified and measured from observations of protostars

Different methods of extracting the rotational component of the gas velocity from observational data are used in the literature. The use of PV diagrams, at disk scales where the high-velocity motions mostly originate from Keplerian rotational motions, has been discussed (e.g., Seifried et al. 2016; Aso & Machida 2020, and references inside). However, at envelope scales, only a few modeling works have been published with the aim of testing the methods. This is key, as the envelope geometry is more complex, and it is unclear whether the largest velocities extracted from PV diagrams are solely due to rotational motions in such conditions, for example.

None of the observed profiles perfectly captures the exact profiles of SAM in the models. We note however that SAM profiles from a weakly magnetized environment will show a smoother decline of the SAM, although these differences may be too subtle to favor one scenario or another on their own in real observations. We summarize our conclusions in the scheme shown in Fig. 8, where we show the typical broken power law radial profiles of SAM obtained by the best observational method, the one in which the rotational velocity is estimated from the peak velocity in the C¹⁸O emission lines ($j_{\text{Vpeak }}^{\text{obs}}$). Altogether, our analysis suggests that larger slopes will be found for j(r) in the outer envelopes, ≳ 1000 au, in the more magnetized case, although these differences can be difficult to detect observationally.

Furthermore, our analysis suggests that the SAM profiles will show very low values, and velocities close to the most widespread velocity resolution allowed by observations of molecular line emission, in environments where magnetic braking is important. In less magnetized environments, the greater amount of SAM accumulated in the equatorial plane makes it easily detectable at scales of a few thousand au, also showing a more monotonous evolution of the SAM with equatorial radii.

Finally, we stress that the spatial distribution of the velocity field in the moment 1 map can be discriminant in identifying objects subject to different magnetic braking efficiencies. In the most magnetized case, the velocity field will present a greater asymmetry, and most of the high-velocity gas will be found around outflow cavities and in the outermost region of the core. In objects where magnetic braking is inefficient, the velocity field will present a higher symmetry with respect to the central object and the high-velocity gas will be found distributed in the equatorial plane.

Fig. 8 Sketch illustrating the influence of the magnetic field on the SAM in protostellar envelopes, in the case of strong magnetic braking (top row) and weak magnetic braking (bottom row). Left: diagram illustrating the distribution of the rotational, νϕ, and radial, υϱ, velocity components in the case of strong or weak magnetic braking in nonideal MHD models. Thin black arrows depict the radial velocity component. The thick black arrow shows the initial orientation of the magnetic field. Red and blue shapes show the rotational component of the velocity field. The dashed magenta line shows the equatorial plane, perpendicular to the rotational axis, along which the SAM profiles were computed. Right: radial profiles of the observed SAM. The blue, red, and black stars show the SAM radial profiles computed from synthetic observations of the C18O (2−1) moment 1 maps, and dashed lines show the best-fit broken power law models for the three evolutionary stages of the nonturbulent (left) and turbulent (right) models. Values of the χ2 per degree of freedom are reported in Table 2.

5.2 Comparison with observational works

Here we discuss whether signatures of magnetic braking have been identified in observational works, and compare our findings to constraints from the observations of embedded protostars.

Observations have revealed a relationship between the kinetic energy of protostellar gas, measured from the velocity gradients at scales a few thousands au in Class 0 envelopes, and the misalignment of the core-scale magnetic field (inferred by dust polarization observations) and the large-scale rotation axis (inferred from the outflow direction) (Galametz et al. 2020). Moreover, comparing the dust polarization from the JCMT BISTRO-1 survey with the gas kinematics traced by the C¹⁸O emission in SMA maps of embedded protostars at ~1000 au scales, Gupta et al. (2022) find that the ratio of rotation motions to infall motions (measured as gradients perpendicular to the outflow axis w.r.t. gradients parallel to the outflow axis) increases in the envelopes with increasing misalignment angles between the core rotation axis and the magnetic field orientation. Since models suggest that a global magnetic field aligned with the initial rotation axis of the core will produce a stronger magnetic braking of envelope motions (Hennebelle & Ciardi 2009; Joos et al. 2012), these results could be tied to the efficiency of magnetic braking in the surveyed protostars. In both studies, the gas velocity is estimated from first moment maps of the molecular emission and the only hypothesis is that velocity gradients perpendicular to the outflow direction preferentially trace rotation, while velocity gradients parallel to the outflow trace infall motions. We stress that estimating the SAM radial profiles from these observations using the methods explored in our work, for example the peak and the maximum velocity, would allow one to test further whether a more aligned magnetic field at the core scale is responsible for braking the gas in inner envelopes. As is explained in Sect. 4.1, using one method or the other should produce large differences (e.g., typically Euclidian distances between the radial profiles obtained >60; see Table D.1) between the SAM profiles when probing a highly magnetized environment.

Several observational works have measured the magnitude of SAM associated with the gas in envelopes of embedded protostars at scales ≲5000 au. Gaudel et al. (2020) analyzed the CALYPSO molecular gas emission maps to estimate the rotational velocity and SAM over three orders of magnitude in scale, in the envelopes of 12 Class 0 protostars. The SAM values computed from the Gaussian fit of molecular line emission profiles is in the range of 10⁻³–10⁻⁴ km s⁻¹ pcat scales between 1500 and 100 au. Yen et al. (2015a) also finds values of specific angular momenta in embedded sources of a few 10⁻⁴ km s⁻¹ pc at 1000 au scales. Pineda et al. (2019) have measured SAM profiles in three Class 0 envelopes, at large scales from 10 000 au down to 1000 au, and find values in the range 10⁻² to a few 10⁻⁴ km s⁻¹ pcfrom large to small scales. The range of values measured by observations are similar to those found in our models (~10⁻⁴−10⁻³ km s⁻¹ pc, or ~10¹⁹−10²⁰ cm² s⁻¹ ). We stress that the magnitudes of SAM found in these observations are in better agreement with the predictions from the more magnetized model as no large values of SAM > 2 × 10⁻³ km s⁻¹ pc, such as the ones measured with synthetic observations of the μ10ℳ0 or μ10ℳ1 models, have been measured at envelope radii of 100–2000 au in any of these Class 0 protostars.

Pineda et al. (2019); Heimsoth et al. (2022); Pandhi et al. (2023) have observed the molecular line emission from samples of cores and find SAM scales with a radius following a power law with the index j ∝ r^~1.8: such large radial variation is more pronounced than that observed in our models. However, the scales probed in these observations are also somewhat larger, reaching >10 000 au envelope radii, and the method used is also slightly different, as the values are averaged over slices parallel to the axis of rotation, while we focus on profiles along the equatorial plane. The limited spatial resolution of these two last observational studies, however, does not allow one to detect a change in SAM profiles in the inner envelopes in which magnetic braking could be at work. Thanks to the high angular resolution of their data, Gaudel et al. (2020) could identify that the SAM in the envelopes shows a broken power law profile, more prominent in some sources from their sample. In sources showing the break, a steep dependency with the envelope radius is observed at radii ≳1600 au (index of 1.6), and a flattening of the radial profile in the inner regions (index of ≳0.3 down to radii ~100 au). While the large slope found at large envelope radii is not present in any of the synthetic observations performed out of the models explored, the transition from a large slope at large envelope radii to smaller variations at small envelope radii is suggestive of the models with efficient magnetic braking. Models with inefficient magnetic braking, μ10ℳ0 and μ10Μ1, do not show such clear broken power law profiles, nor a steep evolution of SAM with radius (see Table C.1 and panels a and e of Fig. 4). Finally, the complex spatial distribution of the velocity field from moment 1 maps, in the CALYPSO data, also suggests that most of the sources may be subject to magnetic braking in their inner envelopes.

Sai et al. (2023) observed the SAM profiles toward sources of various evolutionary stages. Among them, the youngest is the Class 0 protostar IRAS15398-3359: the gas velocities probed in the inner 1000 au of the envelope are extremely small, less than 0.1 km s⁻¹, and the profile has a slope, j ∝ r^−0.7, with very disrupted profiles exhibiting large local variations. Such profiles are more consistent with the profiles found from synthetic observations of models with efficient magnetic braking. Ohashi et al. (2014) computed rotation velocity profiles from the C¹⁸O emission map in the L1527 envelope, using the Gaussian fit method. The slope they find is very shallow (j(r) ~ r^0.3) at envelope radii of 100–800 au, suggestive of the slopes found in models with efficient magnetic braking. However, they find no break in the velocity profile at envelope scales, large values of rotational velocities (between 0.5 and 1 km s⁻¹), and relatively monotonous variations in the rotational velocities with envelope radii, which are features more commonly found in our models with less efficient magnetic braking.

In the B335 Class 0 protostar, polarization of the infrared emission of background stars by the dust contained in the envelope has been used to infer a magnetically supercritical state of the core at scales of ≳ 10 000 au, with a mass-to-flux ratio of around three (Kandori et al. 2020). Moreover, to explain the small upper limit on the disk size (e.g., 10 au, see Yen et al. 2015b), the observed magnetic field geometry has been confronted with predictions from MHD models, which suggest that a rather strong magnetic braking is at work in this object (Maury et al. 2018). Such efficient magnetic braking could be possible also because of the large gas ionization fractions in the inner envelope of B335, (Cabedo et al. 2023) leading to a very good coupling between the B field and the gas. Interestingly, the analysis of molecular line emission has revealed that the SAM radial profile follows r^1.6 at large scales (≳ 5000 au), and flattens at radii of ≳5000 au (Yen et al. 2011, and references inside). Yen et al. (2015b) show that the SAM is 1.5 × 10⁻³ km s⁻¹ pc at a radius of 9000 au; they also estimate an upper limit of ≲7 × 10⁻⁵ km s⁻¹ pc at 1000 au scales, and find an SAM at 100 au scales ~(3–5) × 10⁻⁵ km s⁻¹ pc. These characteristics are similar to those present in the μ3ℳ0 and μ3ℳ1 models (Fig. 7) that we analyzed: the models feature a SAM redistribution up to one order of magnitude in the range of ~ 1000–4000 au, together with a flattening of the radial profile at ≲1000 au, at least in the turbulent case (μ3ℳ1). Therefore, several observed characteristics of the B335 protostar are in agreement with a protostellar evolution and disk formation scenario, including rather efficient magnetic braking, and we note that the signatures seen in our models could be related to the features observed in its gas kinematics.

To summarize, a good fraction of protostars, observed in relatively large samples, exhibit gas kinematics features that could be consistent with the expected signatures of magnetic braking: a flat inner (0.23 and 0.47) profile and a steep increase (0.93 and 0.83) in the SAM at envelope radii of > 1000–2000 au. These signatures, however, remain poorly explored so far in observations of the gas kinematics. Beyond B335, detailed studies of other embedded protostars combining the characterization of the magnetic field topology and of the envelope kinematics may help us to tell in the future whether a magnetized scenario is required to explain the spatial distribution of SAM and disk sizes in Class 0 protostars.

5.3 Caveats

We stress that only models including turbulence and strong magnetic fields can reproduce the very steep profiles observed at large scales (j(r) ∝ r^1.8). Gaudel et al. (2020) suggest that the angular momentum observed at large envelope radii could be inherited from the turbulence cascading down from the interstellar medium (ISM) surrounding the envelopes, as it follows quite closely the Kolmogorov relationship, j ∝ r^5/3, found at larger scales in molecular clouds, both in observations (Hily-Blant et al. 2008; Heyer & Brunt 2012; Heyer & Dame 2015) and in models (Dib et al. 2010; Arroyo-Chávez & Vázquez-Semadeni 2022). Pandhi et al. (2023) also favor a similar interpretation, suggesting that the velocity gradients observed at envelope radii >l000 au have significant contributions from both solid body rotation and turbulent motions. Interestingly, while no large-scale turbulent driving is included in our μ3ℳ0 model, a steep power law behavior consistent with an average index of 5/3 is observed developing in the more magnetized μ3ℳ1 model. However, these large values of the slope are not correctly captured by synthetic observations, which always measure power laws with indices smaller than one.

Some numerical works have included the large-scale physics, such as turbulence (Kuffmeier et al. 2017), external accretion (Kuznetsova et al. 2020), or the influence of filament fragmentation (Misugi et al. 2024), and provide some ideas as to how the properties of the surrounding environment may change the angular momentum of protostellar cores. For example, Misugi et al. (2024) study the evolution of the angular momentum of cores formed through filament fragmentation. They find that the outer angular momentum of cores is inherited from the initial Kolmogorov turbulent velocity field of the filament only for massive cores, while for low-mass cores SAM radial profiles depend on the magnetization of the parent filament (with their strongly magnetized case having an initial μ ~ 6). However, they do not discuss SAM profiles obtained from synthetic observations of the gas kinematics. Kuznetsova et al. (2020) find that the accretion of gas with large angular momentum produces larger average values than in the present work, by one order of magnitude approximately, for both MHD in the ideal limit (2.04×10² km s⁻¹ pc) and HD models (1.72×10² km s⁻¹ pc). A deeper discussion of the contribution of large-scale conditions, such as turbulence, to synthetic observations of SAM profiles at core scales, is beyond the scope of this study, however, as it would require one to analyze synthetic observations of the gas kinematics in nonideal MHD models forming self-coherently cores out of a turbulent magnetized medium. We stress that an extensive model exploration would be required to check whether only strongly magnetized models develop such broken power law features of j(r), as our study is limited to four realizations only, as it intends to showcase the detectability of kinematic features due to magnetic braking in the most favorable conditions.

In our study, we use only protostellar evolution models that include magnetic fields, and no pure hydrodynamical models, which would feature no magnetic braking. However, we point out that our nonturbulent, least magnetized model produces SAM values of about 10⁻³ km s⁻¹ pc, and a flat distribution between ~ 1000 and ~4000 au, similar to the values and trends reported by Takahashi et al. (2016) in their hydrodynamic simulations of collapsing cores. This suggests that the μ10ℳ0 models are a good approximation of situations with near-negligible magnetic braking. Our study is based on the realization of four models of similar envelope mass, and we analyzed three evolutionary stages for each of these three models. These realizations cannot grasp the full extent of possible parameter combinations for the initial conditions, so our work is not easily applicable to all cores observed in the literature. Moreover, the SAM profiles that we extracted and analyzed may not be unique to the conditions of the models presented, and similar profiles may be achieved with different initial conditions. However, an exhaustive analysis exploring all possible combinations of initial model parameters is both prohibitive in computing time and beyond the scope of this work.

We also point out the difficulty in detecting velocity gradients in the outer envelopes (radii ≳1000 au) due to the limited velocity resolution of observations. The magnitude of this limitation can be appreciated in Fig. 7, where the shaded gray areas enclose the region where the observational analysis fails to capture the velocity gradient because they stand below the spectral resolution; that is, under 0.06 km s^−l.

In addition, we emphasize that our analysis focuses on the equatorial plane since we have followed the observers’ strategy: the analysis of SAM radial profiles in the plane where rotation is expected to dominate the gas motions on the line of sight. Observations of non-axisymmetrical gas features in embedded protostars (e.g., Flores et al. 2023; Pineda et al. 2023; Gupta et al. 2024; Cacciapuoti et al. 2024) suggest that some streamers may be playing an additional role in funneling angular momentum from other directions and contributing to the angular momentum responsible for building the disk. As is proposed in Cabedo et al. (2021), such gas streamers could develop as a natural consequence of infall preferentially flowing along magnetic field lines, and indeed such streamers are observed developing in the MHD models that we explore (as in other models presented in the literature). However, from an observational perspective, estimating the angular momentum carried by the streamers to the disk scales is complex and strongly depends on the assumption of an underlying model (e.g., Ulrich 1976; Mendoza et al. 2009). As such, there are very limited observational constraints on the angular momentum carried by these structures; hence, we did not investigate the angular momentum that they may carry inward in the models explored here, and believe that such complementary studies will be better addressed in dedicated future works using multiscale simulations from clouds to disks.

Finally, we stress that, to facilitate the comparison between the SAM of the models and synthetic observations, we have used an approximation of the SAM where the mass is neglected, as is typically computed in real observations. In the Appendix B, we present the actual SAM radial profiles through the equatorial plane. In these profiles, the division between the models prevails according to their magnetic field, showing a similar trend to that of the previously presented profiles (see Fig. B.2).

6 Conclusions

We have used nonideal MHD models of protostellar formation and evolution with different magnetization conditions to identify kinematic features due to magnetic fields braking the rotational motions of the gas. We then performed RT of these models to obtain the expected spectral line emission, and analyzed these synthetic observations following state-of-the-art methods used to interpret observations of gas kinematics in protostellar envelopes. Our analysis shows the following main results:

• The magnetization level of the models produces strong (μ = 3) or weak (μ = 10) magnetic braking, as has already been discussed extensively in the literature. Our analysis shows that some signatures of magnetic braking in the models are present in the profiles of angular momentum computed in the equatorial plane, with two families of profiles being generated depending on the strength of the initial magnetic field. The strongly magnetized model produces a profile with a steep decrease in the angular momentum, of around one order of magnitude, from the outer envelope (>1000 au) to the inner one (<1000 au), while the weakly magnetized model exhibits a softer trend, with a weak decrease in the angular momentum, of around a factor of five, from the outer envelope to the inner one. The effect of the magnetic field is also reflected in the two-dimensional distribution of the SAM. While in more magnetized environments the SAM tends to be redistributed preferentially toward the outer envelope and possibly associated with the formation of outflows, less magnetized environments have most of their SAM associated with the equatorial plane;

• Analyzing the synthetic observations of the C¹⁸O (2−1) molecular line emission, we show that analyzing the gas kinematics can provide clues to distinguish objects experiencing strong magnetic braking from ones where it is absent. First, envelopes experiencing weak magnetic braking exhibit velocity maps from molecular line emission that are organized along symmetrical spatial distributions around the equatorial plane, even in the presence of turbulence: this leads to larger velocity gradients being measured in the equatorial plane. Second, the radial profiles of SAM built along the equatorial plane show different properties, exhibiting a rather monotonous power law decline in the least magnetized cases, and complex profiles when the magnetic braking is strong. The low velocities in the outer envelopes make the detection of rotation difficult with commonly used observational techniques, especially in the most magnetized environments. Therefore, in magnetized environments our work predicts SAM profiles close to the observational limit in terms of velocity resolution and thermal broadening of the line, while less magnetized environments are associated with higher velocities even in the outermost parts of the envelope.

We have presented evidence that the initial magnetization of protostellar cores produces different radial profiles of the SAM of the gas, when measured at core to disk scales in the equatorial plane. We show that observational methods are not able to measure accurately the fast decrease in the angular momentum at large envelope radii, in the case of strong magnetic braking, and that they also have issues in capturing the purely rotational component of the velocity, and hence the magnitudes of the SAM. However, our synthetic observations show that the shapes of observed SAM profiles can still provide clues as to whether protostars experience strong magnetic braking. Further dedicated studies will have to investigate the role of large-scale gas kinematics from the dense ISM in nonideal MHD models, and also examine the role of accretion proceeding from directions other than throughout the equatorial plane.

Acknowledgements

We thank Nagayoshi Ohashi for the discussion and to the Core2disk-III residential program of Institut Pascal at Université Paris-Saclay, with the support of the program “Investissements d’avenir” ANR-11-IDEX-0003-01 for making it possible. We would like to thank the anonymous referee for constructive comments that helped to improve the paper. This project has received funding from the ERC under the European Union’s Horizon 2020 research and innovation program (MagneticYSOS, grant agreement no. 679937). P. Hennebelle, U. Lebreuilly and N. Aflez-Lopez acknowledge financial support from the European Research Council (ERC) via the ERC Synergy Grant ECOGAL (grant 855130).

Appendix A Total number density

Figure A.1 shows the total gas number density for μ3ℳ0 and μ10ℳ0 (upper panel), as is presented in Hennebelle et al. (2020), and μ3ℳ1, and μ10ℳ1 models (bottom panel) developed for the present work. Specifically, the panels show the density in the z−x plane integrated along the y plane.

Fig. A.1 Total gas number density for μ3ℳ0 and μ10ℳ0 models (upper panel: Top and bottom rows, respectively) and for μ3ℳ1 and μ10ℳ1 (bottom panel: Top and bottom rows, respectively) as they evolve (from left to right).

Appendix B Expressions used for the angular momentum and the specific angular momentum

The angular momentum within a core, described as a continuum distribution, is given by the following equation:

$$L={\displaystyle {\int }_{core}r}\times PdV={\displaystyle {\int }_{core}r}\times \rho (r)\upsilon dV\text{, }$$(B.1)

where r stands for the position vector, ρ stands for density, m for the mass, and υ for the velocity of each differential volume. Then, we can approximate this equation for a discrete distribution of particles without a loss of generality:

$$L={\displaystyle \sum _{core}^i{r}_i}\times m_i{\upsilon }_i={\displaystyle \sum _{core}^i{r}_i}\times {\rho }_i{V}_i{\upsilon }_i\text{, }$$(B.2)

where the i index run for every particle in the core. It is equivalent to describe the system as a function of unit density, which allows for easier computation from the adaptative mesh refinement (AMR) model, instead:

$$l\equiv L/V={\displaystyle \sum _{core}^i{r}_i}\times {\rho }_i{\upsilon }_i,$$(B.3)

where ρ_i is the mass per unit of volume for each particle in the system.

In order to compute this quantity from models, and to build radial profiles similar to those analyzed in observational works (which focus on measuring the SAM in the equatorial plane where rotational motions are expected to dominate over infall motions for the line-of-sight component of the velocity vector), we averaged azimutally the 3D angular momentum as follows:

$$l=\left|l\right|=\frac{1}{n}{\displaystyle \sum _{\delta =0}^{2\pi }{r}_{cyl}^{\delta }}{\upsilon }_{\phi }^{\delta }{\rho }^{\delta },$$(B.4)

where stands for each of the n planes that contains the rotational axis. To compute this scalar quantity, we assumed that the radial distance toward the rotational axis, r_cyl, is the magnitude of the position vector, and ν_ϕ is the rotational velocity component perpendicular to the radial vector. Figure B.1 presents the radial profiles of angular momentum L^model (along the equatorial plane) for all the models explored in this work.

Finally, we computed the SAM described as the angular momentum per unit density as follows:

$$l^{\text{model }}=l/\rho =\frac{1}{n}\frac{{\displaystyle \sum _{\delta =0}^{2\pi }{r}_{cyl}^{\delta }}{\upsilon }_{\phi }^{\delta }{\rho }^{\delta }}{{\displaystyle \sum _{\delta =0}^{2\pi }{\rho }^{\delta }}},$$(B.5)

where we divided the angular momentum, l, by the azimutally averaged density. Figure B.2 presents the radial profiles of SAM l^model (along the equatorial plane) for all the models explored in this work.

We stress that, to make comparisons with observable quantities, in the present work we have been using an approximation of the SAM, neglecting the mass of the gas particles, which we refer as SAM j throughout the text. This quantity, j= r × v, is thequantity that can be probed by observations of the protostellar gaseous envelopes, and is the quantity that we use throughout the paper to discuss the ability of observations to distinguish signatures of magnetic braking. This approximation is reasonable, as profiles from the l in models (Figure B.2) show very similar trends to the profiles from this approximated j (Figure 7, panels a and e) used to make comparisons with synthetic observations. This quantity, j, as well as the formulations in Eq. B.4 and B.5, neglect the angular momentum due to the velocity components in the radial and poloidal directions, which contribute a negligible amount to the global angular momentum, as well as the product v_ϕr_z,which is negligible in the equatorial plane that we are examining (where r_z is very small).

Fig. B.1 Radial profiles of the gas angular momentum in the models (defined as L in the description above; see Eq. B.4). The colors and lines are presented as in Fig.4. The left and right panels present nonturbulent and turbulent models, respectively.

Fig. B.2 Radial profiles of the SAM (defined as lmodel in the description above; see Eq.B.5). The colors and lines are presented as in Fig.4. The left and right panels present nonturbulent and turbulent models, respectively.

Appendix C Broken power law fitting

Appendix D Distances

The normalized Euclidean distance was used in the present work to quantify the goodness of approximations by observational methods with respect to the 3D SAM. Using the following equation, we measured the distances between the different profiles (the results are shown in Table D.1):

$$d=100\frac{\sqrt{{\displaystyle {\sum }_{i=1}^n{\left(p_i-q_i\right)}^2}}}{{\displaystyle {\sum }_{i=1}^n\left(p_i+q_i\right)}},$$(D.1)

where p and q are the SAM values for each profile. The distance (d) ranges from 0 to 100 for perfect agreement with the maximum respective difference.

References

All Tables

All Figures

	Fig. 1 Radial profiles of the rotational velocity of the gas in models without initial turbulence (i.e., μ3ℳ0 and μ10ℳ0). The radial profiles were computed from azimuthally averaged models, along the equatorial plane (disk plane).
In the text

	Fig. 2 Radial profiles of the rotational velocity of the gas in models with initial turbulence (i.e., μ3ℳ1 and μ10ℳ1). The radial profiles were computed from azimuthally averaged models, along the equatorial plane (disk plane).
In the text

	Fig. 3 Maps of the azimuthally averaged SAM in the models, $j_{V\phi }^{\text{model }}$. The white star and dashed blue line depict the location of the protostar and the equatorial plane, respectively. The leftmost panels, from (a) to (f), show the nonturbulent models, while the rightmost panels, from (g) to (k), present the turbulent models’ outputs.
In the text

Fig. 4 SAM radial profiles computed through the equatorial plane for nonturbulent (upper row) and turbulent (lower row) models. The first column shows the profiles from the models, while the second to last columns show the profiles from synthetic observations of the models, using different techniques to measure the rotational velocities (see Sect. 3). More (μ = 3.33) and less (μ = 10) magnetized models are depicted with dotted and solid lines, respectively. Different evolutionary states are represented on a color scale: black, blue, and red for protostellar masses equal to ~0.06/0.07, 0.18, and 0.30 M⊙, respectively. The shaded gray area represents the region that cannot be sampled because it lies below the velocity resolution of synthetic observations (0.06 km s−1).

In the text

	Fig. 5 First moment map (color scale) derived from the C18O (2−1) line for μ = 3.33, ℳ = 0 and μ = 10, ℳ = 0 models (upper and bottom rows, respectively). The dashed black line depicts the equatorial plane.
In the text

	Fig. 6 First moment map (color scale) derived from the C18O (2−1) line for μ = 3.33, ℳ = 1 and μ = 10, ℳ = 1 models (upper and bottom rows, respectively). The dashed black line depicts the equatorial plane.
In the text

Fig. 7 SAM of the gas: radial profiles computed along the equatorial plane. They are computed from the gas velocity in the model’s outputs, and from the velocity recovered with synthetic observations of the C18O (2−1) line emission. Model: black circles show 3D SAM computed considering the azimuthally averaged rotational (νϕ) velocity component. Synthetic observation: Red crosses show the SAM computed from the maximum velocity. Green stars show the SAM computed from the intensity-weighted velocity (moment 1). Blue circles show the SAM computed from the Gaussian-fitted peak velocity. The shaded gray area represents the region beyond the velocity resolution. Red, green, and blue shadows show error bars on the synthetic data.

In the text

Fig. 8 Sketch illustrating the influence of the magnetic field on the SAM in protostellar envelopes, in the case of strong magnetic braking (top row) and weak magnetic braking (bottom row). Left: diagram illustrating the distribution of the rotational, νϕ, and radial, υϱ, velocity components in the case of strong or weak magnetic braking in nonideal MHD models. Thin black arrows depict the radial velocity component. The thick black arrow shows the initial orientation of the magnetic field. Red and blue shapes show the rotational component of the velocity field. The dashed magenta line shows the equatorial plane, perpendicular to the rotational axis, along which the SAM profiles were computed. Right: radial profiles of the observed SAM. The blue, red, and black stars show the SAM radial profiles computed from synthetic observations of the C18O (2−1) moment 1 maps, and dashed lines show the best-fit broken power law models for the three evolutionary stages of the nonturbulent (left) and turbulent (right) models. Values of the χ2 per degree of freedom are reported in Table 2.

In the text

	Fig. A.1 Total gas number density for μ3ℳ0 and μ10ℳ0 models (upper panel: Top and bottom rows, respectively) and for μ3ℳ1 and μ10ℳ1 (bottom panel: Top and bottom rows, respectively) as they evolve (from left to right).
In the text

	Fig. B.1 Radial profiles of the gas angular momentum in the models (defined as L in the description above; see Eq. B.4). The colors and lines are presented as in Fig.4. The left and right panels present nonturbulent and turbulent models, respectively.
In the text

	Fig. B.2 Radial profiles of the SAM (defined as lmodel in the description above; see Eq.B.5). The colors and lines are presented as in Fig.4. The left and right panels present nonturbulent and turbulent models, respectively.
In the text