Magnetic field properties in star formation: A review of their analysis methods and interpretation

1 Introduction

Star formation within molecular clouds, the densest part of the interstellar medium (ISM), is regulated by the complex interplay among gravity, turbulence, magnetic fields, and other factors (e.g, protosteller feedback and feedback from previous generations of stars) at different scales. Magnetic fields interact with the other two major forces (gravity and turbulence) by providing supports against gravitational collapse (Shu et al., 1987) and generating anisotropic turbulence (Goldreich and Sridhar, 1995). Observational studies of magnetic fields are crucial to distinguish between strong-field star formation theories in which magnetically sub-critical clouds slowly form super-critical substructures that subsequently collapse (Mouschovias et al., 2006), and weak-field star formation theories where large-scale supersonic turbulent flows form overdense intersecting regions that dynamically collapse (Mac Low and Klessen, 2004).

Polarized thermal dust emission observations have been the most common way to trace the plane-of-sky (POS) component of magnetic field orientations with the assumption that the shortest axis of irregular dust grains is aligned with magnetic field lines (Davis and Greenstein, 1949; Lazarian, 2007; Lazarian and Hoang, 2007). The Goldreich-Kylafis (GK) effect provides an alternative way to trace the POS field orientation (with a 90° ambiguity) with molecular line polarization observations (Goldreich and Kylafis, 1981). The recently developed Velocity Gradient Technique (VGT) proposed another way to trace the POS field orientation with spectra line observations based on the notion that the gradient of velocity centroids (VCG, González-Casanova and Lazarian, 2017) or thin velocity channels (VChG, Lazarian and Yuen, 2018) is perpendicular to the magnetic field due to the intrinsic properties of magneto-hydrodynamic (MHD) turbulence (Goldreich and Sridhar, 1995).

Several analysis techniques have been developed to infer the properties of magnetic fields based on their orientations: the Davis-Chandrasekhar-Fermi (DCF) method was proposed by Davis (1951) and Chandrasekhar and Fermi (1953) approximately 70 years ago and has been the most widely used method to indirectly derive the magnetic field strength with statistics of field orientations. A new tool, the polarization-intensity gradient method (here after the KTH method, Koch et al., 2012a), was proposed about one decade ago and can also be used to assess the significance of magnetic fields based on ideal MHD equations. The Histogram of Relative Orientations (HRO) analysis (Soler et al., 2013), which was proposed right after the KTH method, measures the relative orientation between the magnetic field and density structures and can be used to link the observed magnetic morphology with the physics of simulations. These methods provide information on the magnetic properties in star-forming molecular clouds and allow us to investigate both qualitatively and quantitatively the dynamical role of magnetic fields in the collapse and fragmentation of dense molecular structures.

In this chapter we review the concept and recent developments of these techniques and discuss their limitations. We also summarize the application of these methods to observations of star formation regions and discuss the role of magnetic fields at different spatial scales. In particular, we focus on the relative importance of the magnetic field as compared to gravity and turbulence at different scales of star-forming clouds. In Section 2, we review the DCF method. In Section 3, we review the HRO analysis. In Section 4, we review the KTH method. In Section 5, we summarize this chapter.

2 The DCF method

In the middle 20th century, Davis (1951) and Chandrasekhar and Fermi (1953) proposed the DCF method to estimate the mean¹ magnetic field strength (B^m) of the interstellar medium (ISM) in the spiral arm based on the first interstellar magnetic field orientation observation made by Hiltner (1949). Since then, the method has been improved and adopted by the community to estimate the field strength in star-forming regions. In this section, we present a review of the original DCF method and its modifications.

2.1 Basic assumptions

2.1.1 Energy equipartition

The original DCF method assumes an equipartition between the transverse (i.e., perpendicular to the underlying field B^u) turbulent magnetic and kinetic energies (i.e., the Alfvén relation, hereafter the DCF53 equipartition assumption):

in SI units², where $B_{\perp }^{\text{t}}$ is the transverse turbulent magnetic field, δv_⊥ is the transverse turbulent velocity dispersion, μ₀ is the permeability of vacuum, ρ is the gas density. The DCF53 assumption is also adopted by the recently proposed Differential Measure Approach (DMA, Lazarian et al., 2022). In the POS, the DCF53 assumption yields

where “pos” stands for the POS.

Alternatively, Federrath (2016) assumed an equipartition between the coupling-term magnetic field (${\mathit{B}}_{\Vert }^{\text{t}}\cdot {\mathit{B}}^{\text{u}}$, where “‖” denote the direction parallel to B^u) and the turbulence (hereafter the Fed16 equipartition assumption) when the underlying field is strong. Skalidis et al. (2021) further proposed that only the transverse velocity field is responsible for ${\mathit{B}}_{\Vert }^{\text{t}}\cdot {\mathit{B}}^{\text{u}}$ and the transverse velocity field for the POS coupling-term field can be approximated with the line-of-sight (LOS) velocity dispersion δv_los (see their Section 4.2). Thus they obtained

where the POS transverse velocity dispersion is neglected.

2.1.1.1 Sub-alfvénic case

Conventionally, a sub-Alfvénic state means that the underlying magnetic energy is greater than the turbulent kinetic energy when comparing the magnetic field with the turbulence. It is widely accepted that the DCF53 equipartition assumption is satisfied for pure incompressible sub-Alfvénic turbulence due to the magnetic freezing effect where the perturbed magnetic lines of force oscillate with the same velocity as the turbulent gas in the transverse direction (Alfvén, 1942). However, the star-forming molecular clouds are highly compressible (Hennebelle and Falgarone, 2012). For compressible sub-Alfvénic turbulence, there are still debates on whether the DCF53 or Fed16 equipartition assumptions is more accurate (e.g., Skalidis et al., 2021; Skalidis and Tassis, 2021; Lazarian et al., 2022; Li et al., 2022; Liu et al., 2022).

Observational studies usually adopt the local underlying field within the region of interest instead of the global underlying field at larger scales. In this case, the volume-averaged coupling-term magnetic energy is 0 by definition (Zweibel and McKee, 1995; Liu et al., 2022), which should not be used in analyses. Several numerical studies (Federrath, 2016; Beattie et al., 2020; Skalidis and Tassis, 2021; Beattie et al., 2022) have suggested that the volume-averaged RMS coupling-term magnetic energy fluctuation should be studied instead of the volume-averaged coupling-term magnetic energy. With non-self-gravitating sub-Alfvénic simulations, they found that the coupling-term field energy fluctuation is in equipartition with the turbulent kinetic energy within the whole simulation box. However, it is unclear whether the Fed16 equipartition assumption still holds in sub-regions of their simulations. Investigating the local energetics is very important because the local and global properties of MHD turbulence can be very different (Lazarian and Vishniac, 1999; Brandenburg and Lazarian, 2013). In small-scale sub-regions below the turbulent injection scale and without significant self-gravity, the local underlying magnetic field is actually part of the turbulent magnetic field at larger scales and the local turbulence is the cascaded turbulence (Federrath, 2016). Within self-gravitating molecular clouds, the gravity is comparable to or dominates the magnetic field and turbulence at higher densities (Crutcher, 2012; Kauffmann et al., 2013; Liu et al., 2022), which has a strong effect on both magnetic fields and turbulence. e.g., the gravity can compress magnetic field lines and amplify the field strength; the gravitational inflows can accelerate the gas and enhance turbulent motions. As observations can only probe the magnetic field in part of the diffuse ISM or molecular clouds, it is necessary to test the validity of the Fed16 assumption in sub-regions of simulations with or without self-gravity. Moreover, Li et al. (2022) and Lazarian et al. (2022) pointed out that the physical meaning of the RMS coupling-term energy fluctuation adopted by the Fed16 equipartition assumption is still unclear, which needs to be addressed in the future.

The traditional DCF53 equipartition assumption has been tested by many numerical works (e.g., Haugen and Brandenburg, 2004; Falceta-Gonçalves et al., 2008; Liu et al., 2021; Skalidis et al., 2021; Chen et al., 2022). For non-self-gravitating simulations, Haugen and Brandenburg (2004) found the DCF53 equipartition is violated throughout the inertial range (i.e., between the turbulence injection scale L_inj and dissipation scale) in initially very sub-Alfvénic $({\mathcal{M}}_{A0}=0.05-0.5)$ simulations, while Falceta-Gonçalves et al. (2008) found an exact equipartition between turbulent magnetic and kinetic energies for initially slightly sub-Alfvénic (initial Alfvénic Mach number ${\mathcal{M}}_{A0}=0.7$) models. Another numerical work by Skalidis et al. (2021) with non-self-gravitating simulations found that the DCF53 assumption is approximately fulfilled for initially trans-Alfvénic $({\mathcal{M}}_{A0}=0.7-2)$ models, but the turbulent magnetic energy is much smaller than the turbulent kinetic energy in initially sub-Alfvénic $({\mathcal{M}}_{A0}=0.1-0.5)$ models. Similarly, as these studies adopted the whole simulation-averaged field as the mean field, it is unclear whether these relations still holds in sub-regions where the local properties are dominant. Lazarian et al. (2022) suggested that the DCF53 equipartition in sub-Alfvénic and non-self-gravitating media can be established in the regime of strong turbulence at small scales $(<L_{\mathit{inj}}{\mathcal{M}}_A^2)$. For self-gravitating simulations, a numerical study of clustered star-forming clumps found that the DCF53 energy equipartition assumption is approximately fulfilled in trans-Alfvénic $({\mathcal{M}}_A=0.7-2.5)$ clumps/cores at 1–0.1 pc scales (Liu et al., 2021). Another numerical study by Chen et al. (2022) with self-gravitating and initially trans-Alfvénic $({\mathcal{M}}_{A0}=1.04-1.45)$ simulations found that the DCF53 equipartition approximately holds in the whole simulations box. Chen et al. (2022) further suggests that the DCF53 equipartition breaks in sub-blocks with insufficient cell numbers (e.g., $<41^3$ cells blocks). It is unclear whether the DCF53 energy equipartition assumption still holds in very sub-Alfvénic $({\mathcal{M}}_A\ll 0.7)$ clumps/cores, although the real question may be whether there are very sub-Alfvénic clumps/cores in gravity-accelerated gas (see discussions in Sections 2.4.2 and Section 3.3).

In summary, the DCF53 equipartition assumption is valid within trans- or slightly sub-Alfvénic self-gravitating molecular clouds, but its validity in very sub-Alfvénic self-gravitating regions still needs more investigations. The Fed16 equipartition assumption can be used as an empirical relation when studying the global underlying and turbulent magnetic field in the diffuse ISM beyond the turbulent injection scale if the ISM is sub-Alfvénic, but its physical interpretation, as well as its applicability in part of the ISM below the turbulent injection scale and within self-gravitating molecular clouds, are still unclear. The equipartition problem has only been investigated theoretically and numerically. We are unaware of any observational attempts to look into this problem within molecular clouds yet. There is also a lack of observational methods to study the energy equipartition.

2.1.1.2 Super-alfvénic case

A super-Alfvénic state of turbulence conventionally means that the underlying magnetic energy is smaller than the turbulent kinetic energy. In super-Alfvénic case, the magnetic field is dominated by the turbulent component. Numerical studies have shown that both the turbulent magnetic energy and the RMS coupling-term magnetic energy fluctuation are smaller than the turbulent kinetic energy in super-Alfvénic simulations (Falceta-Gonçalves et al., 2008; Federrath, 2016; Liu et al., 2021). This energy non-equipartition in super-Alfvénic cases could lead to an overestimation of the magnetic field strength. Lazarian et al. (2022) suggested that the super-Alfvénic turbulence transfers to sub-Alfvénic turbulence at $<L_{\mathit{inj}}{\mathcal{M}}_A^{-3}$ scales due to the turbulence cascade in the absence of gravity.

2.1.2 Isotropic turbulence

The original DCF method assumes isotropic turbulent velocity dispersion, thus the unobserved δv_pos⊥ in Equation 2 can be replaced by the observable δv_los, which implies

On the other hand, the modified DCF method in Skalidis and Tassis (2021) (hereafter ST21) requires an assumption of isotropic turbulent magnetic field, so that the $B_{\text{pos}\Vert }^{\text{t}}$ in Equation 3 can be replaced by $B_{\text{pos}\perp }^{\text{t}}$ (Skalidis et al., 2021). Then there is

However, both incompressible and compressible MHD turbulence are anisotropic in the presence of a strong magnetic field (Shebalin et al., 1983; Higdon, 1984; Goldreich and Sridhar, 1995). In particular, the fluctuations of the Alfvénic and slow modes are anisotropic, while only the fast mode has isotropic fluctuations (Lazarian et al., 2022). The anisotropic velocity field in low-density non-self-gravitating regions has been confirmed by observations (Heyer et al., 2008). Lazarian et al. (2022) suggested that the anisotropy of MHD turbulence is a function of ${\mathcal{M}}_A$ and relative fraction of MHD modes. The anisotropy of turbulence is also scale-dependent, which increases with decreasing scales (Lazarian and Vishniac, 1999). The recently developed DMA does not require the assumption of isotropic turbulence, but analytically considers the anisotropy of MHD turbulence in non-self-gravitating turbulent media. In particular, the DMA has written the POS transverse magnetic field structure function $({\tilde{D}}_{yy})$ and the POS velocity centroid structure function $({\tilde{D}}^v)$ as functions of the turbulent field fluctuation and the velocity dispersion, respectively. Both ${\tilde{D}}_{yy}$ and ${\tilde{D}}^v$ are complicated functions of the mean field inclination angle³ γ, ${\mathcal{M}}_A$, composition of turbulence modes, and distance displacements. We refer the readers to the original DMA paper (Lazarian et al., 2022) for the derivation of ${\tilde{D}}_{yy}$ and ${\tilde{D}}^v$.

With non-self-gravitating simulations, Heitsch et al. (2001) and Skalidis et al. (2021) found that the full cube turbulent magnetic field is approximately isotropic within a factor of 2 with their ${\mathcal{M}}_A$ = 0.1–14 simulations, but they did not investigate the property of turbulence in sub-regions at smaller scales where the anisotropy is expected to be more prominent. Both works found the turbulent field is more anisotropic for models with smaller ${\mathcal{M}}_A$ values.

Studies on the anisotropy of MHD turbulence in the high-density self-gravitating regime are rarer. Spectroscopic observations toward high-density regions did not find obvious evidences for velocity anisotropy (Heyer and Brunt, 2012). A recent numerical work by Liu et al. (2021) found both the turbulent magnetic field and the turbulent velocity dispersion are approximately isotropic within a factor of 2 at 1–0.01 pc scales in their trans-to-super-Alfvénic simulations of clustered star-forming regions. They did not find obvious relation between the anisotropy and levels of initial magnetization in their simulations. This may be due to the local super-Alfvénic turbulence at high densities and/or the complex averaging along the LOS for different local anisotropic turbulent field at various directions. Moreover, Otto et al. (2017) found that the velocity anisotropy due to magnetic fields disappears in high-density and sub-Alfvénic (local ${\mathcal{M}}_A=0.54$) regions of their self-gravitating simulations. Thus, the uncertainty brought by the anisotropic MHD turbulence should be a minor contributor for the DCF method in star-forming clumps/cores where self-gravity is important.

2.1.3 Tracing ratios between field components with field orientations

In the POS, the local total field at position i is a vector sum of the underlying field and the local turbulent field: ${\mathit{B}}_{\text{pos}}^{\text{tot,\hspace{0.17em}i}}={\mathit{B}}_{\text{pos}}^{\text{u}}+{\mathit{B}}_{\text{pos}}^{\text{t,\hspace{0.17em}i}}$. Equation 4 can be rewritten as

or

to derive the POS underlying magnetic field strength $B_{\text{pos}}^{\text{u}}$ or the POS total magnetic field strength $B_{\text{pos}}^{\text{tot}}$. Similarly for the ST21 approach, Equation 5 can be rewritten as

to derive the POS underlying magnetic field strength. Both the original DCF method and the ST21 approach have adopted the assumption that the turbulent-to-underlying field strength ratio $B_{\text{pos}\perp }^{\text{t}}/B_{\text{pos}}^{\text{u}}$ or turbulent-to-total field strength ratio $B_{\text{pos}\perp }^{\text{t}}/B_{\text{pos}}^{\text{tot}}$ in the above equations can be estimated from statistics of the POS magnetic field orientations, which is usually done by calculating the dispersion of magnetic field position angles or applying the angular dispersion function method (hereafter the ADF method, Hildebrand et al., 2009; Houde et al., 2009; Houde et al., 2016). On the other hand, the DMA assumes that the POS polarization angle structure function D^ϕ can be expressed as a function of the ratio between the POS transverse and total magnetic field structure functions $({\tilde{D}}_{yy}/{\tilde{D}}_{\mathit{Bpos}})$, where the term ${\tilde{D}}_{yy}/{\tilde{D}}_{\mathit{Bpos}}$ can be expressed as a function of the turbulent-to-underlying field strength ratio.

2.1.3.1 Angular dispersions

The underlying magnetic field reflects the intrinsic property of an unperturbed magnetic field. There are different approaches to relate the dispersion of magnetic field position angles with $B_{\text{pos}\perp }^{\text{t}}/B_{\text{pos}}^{\text{u}}$. Table 1 summarises these angular relations and the corresponding DCF formulas for $B_{\text{pos}}^{\text{u}}$.

TABLE 1. Angular relations and corresponding formulas to derive $B_{\text{pos}}^{\text{u}}$.

All the relations listed in Table 1 have assumed that $B_{\text{pos}\Vert }^{\text{t}}\ll B_{\text{pos}}^{\text{u}}$, so the $B_{\text{pos}\Vert }^{\text{t}}$ can be neglected and the dispersion on the field angle ⁴ of the magnetic field can be approximated with the dispersion of the position angle. The angular relation $B_{\text{pos}\perp }^{\text{t}}/B_{\text{pos}}^{\text{u}}\sim \delta \varphi$ has adopted an additional small angle approximation (i.e., δϕ ∼ tan δϕ ∼ δ(tan ϕ) ≪ 1 or $B_{\text{pos}\perp }^{\text{t}}\ll B_{\text{pos}}^{\text{u}}$). Ostriker et al. (2001) and Liu et al. (2021) found that the angle limit for the small angle approximation is δϕ ≲ 25° in their simulations. Liu et al. (2021) also suggested that δϕ can significantly underestimate $B_{\text{pos}\perp }^{\text{t}}/B_{\text{pos}}^{\text{u}}$ for large δϕ values. For the relation $B_{\text{pos}\perp }^{\text{t}}/B_{\text{pos}}^{\text{u}}\sim \delta (\tan \varphi )$, the simulations by Heitsch et al. (2001) and Liu et al. (2021) suggested that δ(tan ϕ) can show significant scatters due to large values of δ(tan ϕ) when ϕ ∼ 90°. Thus, δ(tan ϕ) is valid to trace the $B_{\text{pos}\perp }^{\text{t}}/B_{\text{pos}}^{\text{u}}$ when ϕ is small, and δ(tan ϕ) reduces to δϕ in such cases. In addition, Chen et al. (2022) found that tan δϕ does not trace $B_{\text{pos}\perp }^{\text{t}}/B_{\text{pos}}^{\text{u}}$ well in their simulations.

The total magnetic field is the sum of the underlying magnetic field and the turbulent magnetic field. There are also different approaches trying to relate the dispersion of magnetic field position angles with $B_{\text{pos}\perp }^{\text{t}}/B_{\text{pos}}^{\text{tot}}$. Table 2 summarises these angular relations and the corresponding DCF formulas for $B_{\text{pos}}^{\text{tot}}$.

TABLE 2. Angular relations and corresponding formulas to derive $B_{\text{pos}}^{\text{tot}}$.

Similarly, all the relations listed in Table 2 have approximated the dispersion on the field angle of the magnetic field with the dispersion of the position angle, which requires $B_{\text{pos}\Vert }^{\text{t}}<B_{\text{pos}}^{\text{u}}$. The numerical study by Liu et al. (2021) found that δϕ correlates with $B_{\text{pos}\perp }^{\text{t}}/B_{\text{pos}}^{\text{u}}$ or $B_{\text{pos}\perp }^{\text{t}}/B_{\text{pos}}^{\text{tot}}$ in small angle approximation, but δϕ estimates $B_{\text{pos}\perp }^{\text{t}}/B_{\text{pos}}^{\text{tot}}$ rather than $B_{\text{pos}\perp }^{\text{t}}/B_{\text{pos}}^{\text{u}}$ for large δϕ values. Liu et al. (2021) also suggested that δ(sin ϕ) provides a better estimation of $B_{\text{pos}\perp }^{\text{t}}/B_{\text{pos}}^{\text{tot}}$ than δϕ. Due to the scatters of δ(tan ϕ), the formula $\sqrt{{\mu }_0\rho }\delta v_{\text{los}}\times \frac{{(1+3\delta {(\tan \varphi )}^2)}^{1/2}}{\delta (\tan \varphi )}$ does not correctly estimate the total magnetic field strength (Heitsch et al., 2001; Liu et al., 2021). The angular relation $B_{\text{pos}\perp }^{\text{t}}/B_{\text{pos}}^{\text{tot}}\sim \tan \delta \varphi$ has not been tested by simulations yet. In addition, Chen et al. (2022) found that sin δϕ and $B_{\text{pos}\perp }^{\text{t}}/B_{\text{pos}}^{\text{tot}}$ are well correlated in regions where the polarization percentage is larger than 20% of its maximum value in the synthetic polarization maps.

2.1.3.2 The ADF method

Structure functions and correlation functions have been widely used in astrophysical studies. Falceta-Gonçalves et al. (2008) introduced the structure function in the study of polarization position angles. Later, the ADF method was developed by Hildebrand et al. (2009) to estimate the POS turbulent-to-ordered field strength ratio $B_{\text{pos}\perp }^{\text{t}}/B_{\text{pos}}^{\text{o}}$ based on the structure function of magnetic field position angles, where B^o is the ordered field. The ADF approach developed initially in Hildebrand et al. (2009) (Hil09) only corrects for the large-scale ordered field structure (see Section 2.2.1 for more discussions about the ordered field). Later, the ADF technique was extended by Houde et al. (2009) (Hou09) and Houde et al. (2016) (Hou16) to be applied to single-dish and interferometric observational data by additionally taking into account the LOS signal integration over multiple turbulent cells, the beam-smoothing effect, and the interferometer filtering effect. The Hou09 and Hou16 variants of the ADF are in the form of the auto-correlation function, which transform to the structure function in the small angle limit. Basically, the ADF of magnetic field orientations are fitted to derive $B_{\text{pos}\perp }^{\text{t}}/B_{\text{pos}}^{\text{o}}$ or $B_{\text{pos}\perp }^{\text{t}}/B_{\text{pos}}^{\text{tot}}$. The simplest ADF only accounting for the ordered field structure has the form (Hildebrand et al., 2009; Houde et al., 2009; Houde et al., 2016)

where ΔΦ(l) is the POS angular difference of two magnetic field segments separated by a distance l, $a_2^{\prime }{l}^2$ is the second-order term of the Taylor expansion of the ordered component of the correlation function, ${(\frac{B_{\text{pos}\perp }^{\text{t}}}{B_{\text{pos}}^{\text{tot}}})}_{\text{adf}}^2=1/(1+1/{(\frac{B_{\text{pos}\perp }^{\text{t}}}{B_{\text{pos}}^{\text{o}}})}_{\text{adf}}^2)$, and the subscript “adf” indicates ADF-derived parameters. The variation of the ordered field is characterised by a scale l_d and it is assumed that the higher-order terms of the Taylor expansion of the ordered field do not have significant contributions to the ADF at l < l_d. The ADF method also assumes that the local turbulent correlation scale (see Section 2.2.2) is smaller than l_d. Because ΔΦ(l) is constrained to be within [ − 90, 90] degrees (i.e., the field angular dispersion approximated with the position angular dispersion) and $a_2^{\prime }$ is defined as positive values (i.e., positive ordered field contribution), the maximum values of the derivable integrated turbulent-to-ordered and -total strength ratio from the simplest ADF are ∼0.76 and ∼0.6 (Liu et al., 2021), respectively. Due to the space limit of this review, we refer readers to the original ADF papers for more complicated forms of ADFs and their detailed derivations. The validity of the ADF method is further discussed in Section 2.2.

2.1.3.3 The DMA

Combining Eq. 1 and the ratio between the polarization angle structure function D^ϕ and the velocity centroid structure function ${\tilde{D}}^v$, the DMA estimates the field strength as

where the factor f is a function of γ, ${\mathcal{M}}_A$, and composition of turbulence modes. Note that the DMA assumes the velocity and magnetic field have the same scaling, therefore, f does not depend on the distance displacement. Lazarian et al. (2022) has listed the formula of f in different physical conditions (see their Table 3). Note that their structure function $D^{\varphi }(l)=\frac{1}{2}⟨1-\cos [2\mathrm{\Delta }\mathrm{\Phi }(l)]⟩$ is different from the one adopted by the ADF method (the left term of Equation 9). Lazarian et al. (2022) claims that D^ϕ(l) is applicable to cases of large angle fluctuations while the 1 − ⟨ cos[ΔΦ(l)] adopted by the ADF method is not.

TABLE 3. Correction factors and corresponding simulation parameters.

2.2 Uncertainties in the statistics of field orientations

As stated in Section 2.1.3, the turbulent-to-underlying or -total magnetic field strength ratio is assumed to be traced by statistics of magnetic field position angles. Other than the uncertainty on the assumption itself, there are various effects that could introduce uncertainties in the statistics of position angles. Here we describe these effects and summarize on how they are treated in different approaches. Note that the estimation of gas density from dust emission is associated with uncertainties on the dust-to-gas ratio, temperature, dust opacity, and source geometry (e.g., Hildebrand, 1983; Ossenkopf and Henning, 1994), whereas the statistics of turbulent velocity field using spectroscopic data is affected by the chemical processes and excitation conditions of particular line tracers (e.g., van Dishoeck et al., 1993), density fluctuations (e.g., Yuen et al., 2021; Yuen et al., 2022), and ordered velocity fields due to star-forming activities (e.g., infall, rotation, and outflows). We do not discuss the uncertainty on the gas density and turbulent velocity field in this paper as it is beyond the scope of this review.

2.2.1 Ordered field structure

The original DCF method was derived assuming the large-scale ordered field lines are straight. For a non-linear large-scale field, the contribution from the ordered field structure can overestimate the angular dispersion that should be only attributed to turbulence. For highly ordered magnetic fields, the underlying field structure can be fitted with simple hourglass models (hereafter the HG technique. e.g., Lai et al., 2002; Girart et al., 2006; Rao et al., 2009; Qiu et al., 2014) or even more complex models (e.g., Myers et al., 2018).

Pillai et al. (2015) (the spatial filtering technique. Hereafter the Pil15 technique) and Pattle et al. (2017) (the unsharp masking technique. Hereafter the Pat17 technique) tried to universally derive the ordered field orientation at each position by smoothing the magnetic field position angle among neighboring positions. Pattle et al. (2017) tested the Pat17 technique with a set of Monte Carlo simulations and found that this technique does correctly recover the true angular dispersion if the measurement error is small. By varying the smoothing size, the Pil15 and Pat17 techniques can account for the ordered structure at different scales.

The ADF method analytically takes into account the ordered field structure (see Section 2.1.3.2 for detail) and has been the most widely used method to remove the contribution from the ordered field to the angular dispersion. With a set of Monte Carlo simulations, Liu et al. (2021) found that the ADF method works well on accounting for the ordered field. Figure 1 shows the overestimation factor of the angular dispersion due to the POS ordered field structure, which is quantified by the ratio R_o between the directly estimated angular dispersion δϕ_obs and the ADF-derived integrated (i.e., without corrections for the LOS signal integration. See Section 2.2.3) turbulent-to-ordered magnetic field strength ratio⁵ ${(\frac{B_{\text{pos}\perp }^{\text{t}}}{B_{\text{pos}}^{\text{o}}})}_{\text{adf,int}}$ from a compilation of previous DCF estimations (Liu et al., 2022). Figure 1 does not show strong relations between R_o and n_{H_{2}}, which contradicts with the expectation that the ordered field structure is more prominent in higher-density regions where gravity is more dominant (e.g., Hull and Zhang, 2019). However, most of the estimations shown in Figure 1 were from the simplest Hil09 variant of ADF, which considers less effects and gives more uncertain results. A revisit of the data in the literature with the more refined Hou09/Hou16 approaches should give a more reliable relation between the ordered field contribution and the gas density. There is a group of data points with R_o < 1 values derived from the Hou09 approach, which implies that the contribution from the ordered field is negative (i.e., $a_2^{\prime }<0$) and is not physical. This is because the original studies for those data points applied the Hou09 variant of ADF to interferometric data and/or fitted the ADF within an upper limit of l that is too large. With Monte Carlo simulations, Liu et al. (2019) found that the sparse sampling of magnetic field detections can generate an artificial trend of decreasing ADF with increasing l at large l, which can explain the unphysical $a_2^{\prime }<0$ values. The average R_o for R_o > 1 values is 2.2 ± 1.1.

FIGURE 1. Overestimation factor of angular dispersions (R_o) due to contributions from the ordered field structure as a function of gas number density based on the DCF catalogue from Liu et al. (2022) which collected the estimations in the literature. Here δϕ_obs is the directly estimated angular dispersion for the simplest DCF version as mentioned in Section 2.1.3.1. Different colors indicate different ADF variants used to derive ${(\frac{B_{\text{pos}\perp }^{\text{t}}}{B_{\text{pos}}^{\text{o}}})}_{\text{adf,int}}$.

A deficiency of the original ADF methods (Hil09/Hou09/Hou16) is a lack of a universal way to define the fitting upper limit of l for the $a_2^{\prime }{l}^2$ term. Therefore, users usually perform the fitting of the ADF within an arbitrary range of l, which can give very different results depending the adopted l range. Lazarian et al. (2022) suggested that the ordered field contribution to the ADF can be removed with multi-point ADFs (Cho, 2019), which has the advantage of not requiring fitting the ordered field contribution but it is at the expense of an increased noise.

It should be noted that the concept of the “ordered” field is vague and is not well defined. The referred entity of an local ordered field depends on the range of scales (i.e., resolution to maximum size) of the region of interest. An example is that the simple hourglass-like magnetic field in G31.41 at a lower resolution (Girart et al., 2009) show complex poloidal-like structures at a higher resolution (Beltrán et al., 2019). It should also be noted that the non-linearly ordered field structure is not only due to non-turbulent processes such as gravity, shocks, or collisions, but can also result from larger-scale turbulence, where the curvature of the ordered field generated by pure turbulence depends on ${\mathcal{M}}_A$ (e.g., Yuen and Lazarian, 2020). The above mentioned techniques (except for the HG technique, where the hourglass shapes are often associated with structures under gravitational contraction) remove the contribution from non-turbulent ordered field structures as well as the contribution from the low spatial frequency turbulent field.

2.2.2 Correlation with turbulence

It is assumed that the turbulent magnetic field is characterized by a correlation length l_δ (Myers and Goodman, 1991; Hildebrand et al., 2009). The turbulent magnetic fields within a turbulent cell of length l_δ are mostly correlated with each other, while the turbulent fields separated by scales larger than l_δ are mostly independent.

Hou09 and Hou16 assumed a Gaussian form for the turbulent autocorrelation function and included it in the ADF analysis. Figure 2 shows the relation between l_δ derived from the ADF method and the resolution l_beam of the observations from the DCF catalogue compiled by Liu et al. (2022). There is an overwhelming trend that the l_δ and l_beam are correlated with each other within a factor of 2. At l_beam ∼ 5 mpc, there is a group of data points with l_δ > 2l_beam, which mostly correspond to the estimations from Palau et al. (2021).

FIGURE 2. The correlation length (l_δ) of the turbulent magnetic field derived from the ADF method versus the beam size (l_beam) of observations based on the Liu et al. (2022) catalogue which collected the estimations in the literature. Solid line corresponds to 1:1 relation. Dashed lines correspond to 1:2 and 2:1 relations.

The smallest observed l_δ is l_δ,min ∼0.6 mpc estimated by Hull et al. (2017a) in Serpens-SMM1. This scale is consistent with the lower limit of the observed ambipolar diffusion scale⁶ (∼1.8–17 mpc, Li and Houde, 2008; Hezareh et al., 2010; Houde et al., 2011) of ion-neutral decoupling, although recently (Tang et al., 2018) suggested that the previous observational estimates of ambipolar diffusion scales were biased towards small values due to the imperfect method used to analyse the turbulent spectra and the true ambipolar diffusion scale may be much larger (e.g., ∼0.4 pc in NGC 6334). The estimated l_δ,min is an order of magnitude smaller than the scale of the observed lower end of the Larson’s law (∼10 mpc, Padoan et al., 2009; Kauffmann et al., 2013; Yuen et al., 2022).

What can we learn from the correlation between l_δ and l_beam? One may intuitively think that there is an intrinsic turbulent correlation scale, which is overestimated at insufficient beam resolution. However, in this scenario, the observed angular dispersions at larger scales should be smaller than those at smaller scales due to the signal integration across more turbulent cells along the LOS (see Section 2.2.3), which contradicts observational results (Liu et al., 2022). Thus, it is reasonable to think that the magnetic field and turbulence are correlated at different scales instead of only at the smallest scale (Heitsch et al., 2001).

Alternatively, we propose that the local turbulence confined by the range of scales of the observations (from the size of the considered region, the effective depth along the LOS, or the filtering scale of interferometric observations to the beam resolution) is responsible for generating the local turbulent magnetic field at the particular scales. Note that the local ordered magnetic field partly consists of the turbulent magnetic field at larger scales and the two are not distinguishable in observations of limited scale ranges. i.e., the contribution of low spatial frequency turbulent fields can be removed by the ordered field term in the ADF analysis. In addition, coarser resolutions cannot resolve the high spatial frequency turbulence. Thus, the turbulent correlation scale derived from the ADF method actually corresponds to the local turbulent power spectra with cutoffs at the resolution and the maximum recoverable scale of particular observations. This local turbulence correlation scale in observations of molecular clouds should be much smaller than the driving scale of interstellar supersonic turbulence (∼100 pc, e.g., Heyer and Brunt, 2004; Haverkorn et al., 2008).

Densities may be another factor that affects the local turbulence. The property of the local turbulence relies on the gas component probed by the observations. At smaller scales, the higher densities decrease the mean free path of gas particles and thus the turbulent correlation length, which can also explain the scaling relation seen in Figure 2.

2.2.3 LOS signal integration

The contribution of the turbulent field responsible for the observed polarization is the mean of the N independent samples of the turbulent field along the LOS (Zweibel, 1990; Myers and Goodman, 1991). The integrated turbulent field seen by observations should be $1/\sqrt{N}$ times of the intrinsic turbulent field. Note that the measured angular dispersion in polarization maps is an approximation of $B_{\text{pos}\perp }^{\text{t}}/B_{\text{pos}}^{\text{u}}$ or $B_{\text{pos}\perp }^{\text{t}}/B_{\text{pos}}^{\text{tot}}$. If both the intrinsic and integrated turbulent field dominate the underlying field, the observed and intrinsic angular dispersions should be close to the value expected for a random field and there is little difference between them. Thus, the factor $\sqrt{N}$ is only an upper limit of the underestimation factor for the angular dispersion.

Lazarian et al. (2022) suggested that the underestimation of angular dispersion is dependent on the mixture of turbulence modes and the inclination of the mean field. With pure Alfvénic simulations and simulations of equal Alfvénic and slow modes, they found that the integrated angular dispersion decreases slower than N^−1/2 when the LOS integration depth L is smaller than the turbulence injection scale L_inj, but decreases faster than N^−1/2 when L ≳ L_inj. They also suggested that pure Alfvénic fluctuations decrease as N⁻¹ instead of N^−1/2 when the mean field is on the POS, which can quickly vanish if the integration depth is greater than the scale of the studied region. Liu et al. (2021) tested the LOS signal integration effect of angular dispersions with self-gravitating simulations. They found that the angular dispersion is only underestimated by a factor of $<$2 at scales $>$0.1 pc, while the angular dispersion can be significantly underestimated at scales $<$0.1 pc. However, both works only investigated the underestimation of the angular dispersion, which does not necessarily reflect the underestimation of the turbulent magnetic field. Future numerical studies should investigate the effect of LOS signal integration on the turbulent magnetic field as well as on the angular dispersion.

The ADF method derives the POS turbulent correlation length by fitting the ADFs and uses this information to derive the number of turbulent cells N_adf along the LOS under the assumption of identical LOS and POS turbulent correlation lengths. Note that this assumption is not satisfied in anisotropic MHD turbulence. Liu et al. (2021) applied the ADF method to simulations and found a significant amount of observed angular dispersions corrected by $\sqrt{N_{\text{adf}}}$ exceed the value expected for a random field, which they interpreted as N_adf being overestimated by the ADF method. However, as mentioned above, $\sqrt{N_{\text{adf}}}$ is only an upper limit of the underestimation factor for angular dispersions. Thus their results do not necessarily mean that the ADF method inaccurately accounts for the LOS signal integration. More numerical tests of the ADF method is required to address its validity on the signal integration effect.

Cho and Yoo (2016) proposed an alternative approach (here after the CY16 technique) to estimate the number of turbulent cells along the LOS

where δV_c is the standard deviation of centroid velocities. The CY16 technique was developed to account for the LOS signal integration at scales larger than the injection scale in non-self-gravitation media, which does not naturally extend to small-scale and high-density media where self-gravity is important. Liu et al. (2021) tested the CY16 technique with self-gravitating simulations. They found that the observed angular dispersions in synthetic polarization maps corrected by $\sqrt{N_{\text{cy}16}}$ agrees with the angular dispersions in simulation grids at scales $>$0.1 pc, but the correction fails at $<$0.1 pc. The failure of the CY16 technique at $<$0.1 pc scales may be related to the complex velocity fields associated with star formation activities. Similar to the removal of ordered magnetic fields, Lazarian et al. (2022) suggests that the contribution from non-turbulent velocity fields to δV_c can be removed by multi-point structure functions.

The non-uniform and complex density and magnetic field structures along the LOS tend to reduce the measured angular dispersion (e.g., Zweibel, 1990). The distribution of the LOS grain alignment efficiency can also affect the derived angular dispersion (e.g., Padoan et al., 2001). The reduction factor of angular dispersions due to these effects is highly dependent on the physical conditions of individual sources and cannot be solved universally.

2.2.4 Observational effects

The observed sky signals are limited by the angular resolution and the sampling of the particular observations. Interferometric observations are further affected by the filtering of large-scale signals.

As discussed in Section 2.2.2, the magnetic field is likely perturbed by the turbulence of different wavenumbers at different scales. In this case, the loss of turbulent power due to beam-smoothing can underestimate the angular dispersion of magnetic field position angles and thus overestimate the magnetic field strength (Heitsch et al., 2001). Falceta-Gonçalves et al. (2008) investigated the resolution effect with numerical simulations and suggested that the ratio between the derived field strength at different spatial resolutions (R_s) and the intrinsic field strength at an infinite resolution follows an empirical relation $(1+C/R_s^{0.5})$, where C is a constant obtained by fitting the relation. It is unclear whether this empirical relation is applicable to observations or other simulations. The ADF method is another approach trying to universally solve the beam-smoothing effect, which analytically introduces a Gaussian profile to describe the telescope beam. Liu et al. (2021) tested the ADF method with simulations and found that this method does correctly account for the beam-smoothing effect. Liu et al. (2021) also suggested that the information of the turbulent magnetic field is not recoverable if the polarization source is not well resolved (i.e., size of the polarization detection area smaller than 2–4 times of the beam size).

The minimal separation of antenna pairs in interferometric observations limits the largest spatial scale that the observations are sensitive to. This filtering effect also introduces uncertainties in the estimation of the turbulent magnetic field. The ADF method attempts to analytically solve this problem by modeling the dirty beam of interferometers with a twin Gaussian profile. With their numerical test, Liu et al. (2021) found the ADF method correctly accounts for this large-scale filtering effect as well.

For observations of polarized starlight extinction or low signal-to-noise-ratio dust polarization emission, the polarization detection is sparsely sampled. Soler et al. (2016) has suggested that this sparse sampling effect can introduce jitter-like features in the ADFs and affect the accuracy of the ADF analysis. Although there are no universal solutions to correct for the effect of sparse sampling, the uncertainties of the ADFs due to the pure statistics can be modeled with Monte Carlo analyses (e.g., Liu et al., 2019). With a simple Monte Carlo test, Liu et al. (2019) found that the ADF of sparsely sampled magnetic field orientations is underestimated at large distance intervals and has larger scatters compared to the ADF of uniformly sampled field orientations. The sparse sampling not only affects the ADF, but can also affect the velocity structure function (VSF). Ha et al. (2021) found similar jitter-like features on the VSF of sparsely sampled stars and estimated the uncertainty of the VSFs with a Monte Carlo random sampling analysis.

2.2.5 Projection effects

The angular relations in Section 2.1.3 are originally derived in 3D spaces where the orientation of the 3D magnetic field is known. Dust polarization observations only traces the POS field orientations, thus the DCF method can only measure the POS magnetic field strength. There are different attempts to reconstruct the 3D magnetic field from the DCF estimations.

The 3D magnetic field can be derived by including the inclination angle of the 3D field in the DCF equations (e.g., Houde, 2004; Lazarian et al., 2022). Note that the inclination angle of the 3D field is only meaningful when the underlying field is prominent (i.e., ${\mathcal{M}}_A\lesssim 1$). Houde et al. (2002) proposed a technique to derive the inclination angle with the combination of the polarimetric data and ion–to–neutral line width ratios. However, Houde (2011) later suggested that this technique cannot be widely applied due to the degeneracy between the inclination angle and the POS field strength. Chen et al. (2019) developed a technique to estimate the inclination angle from the statistical properties of the observed polarization fraction, but this technique is subject to large uncertainties (up to 30°). On the other hand, Hu et al. (2021) suggested that the field inclination angle can be derived with the anisotropy of the intensity fluctuations in spectroscopic data. However, this approach suggested by Hu et al. (2021) requires sophisticated multi-parameter fittings of idealized datasets, which limits its application to observations. The applicability of their technique in self-gravitating regime is also questionable since the gravity-induced motion can significantly affect the pure velocity statistics (Hu et al., 2022).

The 3D field strength can also be estimated by combining the POS and LOS components of the magnetic field. This can be done by combining DCF estimations and Zeeman estimations of the same material, where Zeeman observations are the only way to directly derive the LOS magnetic field strength. Recently, Tahani et al. (2018) proposed a new method to estimate the LOS magnetic field strength based on Faraday rotation measurements, which can also be used to reconstruct the 3D magnetic field in combination with the POS magnetic field estimations (e.g., Tahani et al., 2019).

In most cases, the information on the inclination angle and the LOS magnetic field is missing, or the measured POS and LOS magnetic fields do not correspond to the same material. One may still obtain an estimate of the 3D field from the POS field based on statistical relations. Crutcher et al. (2004) suggested that the statistical relation between the 3D and POS underlying magnetic field is $B_{3\text{D}}^{\text{u}}=\frac{4}{\pi }{B}_{\text{pos}}^{\text{u}}$ for a randomly inclined field between 0 and π/2. Note that this statistical relation only applies to a collection of DCF estimations where the 3D field orientation is expected to be random, but should not be applied to individual observations. For an isotropic turbulent magnetic field, the relation between the 3D and POS turbulent field is $B_{3\text{D}}^{\text{t}}=\sqrt{\frac{3}{2}}{B}_{\text{pos}}^{\text{t}}$. Therefore, the total field has the statistical relation $B_{3\text{D}}^{\text{tot}}=f{B}_{\text{pos}}^{\text{tot}}$, where the factor f should be between $\frac{4}{\pi }\sim 1.27$ and $\sqrt{\frac{3}{2}}\sim 1.22$ depending on whether the underlying or turbulent field is more dominant. For an anisotropic turbulent field, the statistical relation between 3D and POS turbulent or total field should depend on the extent of the anisotropy.

2.3 Simulations and correlation factors

Due to the uncertainties in the assumptions of the DCF method (see Section 2.1) and on the statistics of field orientations (see Section 2.2), a correction factor is required to correct for the magnetic field strength estimated from different modified DCF formulas. Several studies (Heitsch et al., 2001; Ostriker et al., 2001; Padoan et al., 2001; Falceta-Gonçalves et al., 2008; Liu et al., 2021; Skalidis and Tassis, 2021; Chen et al., 2022) have numerically investigated the correction factor Q_c at different physical conditions with 3D ideal compressible MHD simulations. Table 3 summarizes these numerically derived correction factors. Recently, Lazarian et al. (2022) made the attempt to analytically derive the correction factor⁷ f for their proposed DMA formula.

For the most widely used formula $B_{\text{pos}}^{\text{u}}\sim \sqrt{{\mu }_0\rho }\frac{\delta v_{\text{los}}}{\delta {\varphi }_{\text{obs}}}$, Ostriker et al. (2001) made the first attempt to quantify its uncertainty and derived a correction factor of ∼0.5 for an initially slightly sub-Alfvénic (${\mathcal{M}}_{A0}$ = 0.7) model with n₀ ∼ 10² cm⁻³ and a box length of 8 pc. Later, Padoan et al. (2001) found Q_c ∼ 0.4 for three selected ∼1 pc and n < 10⁵ cm⁻³ clumps in a initially slightly super-Alfvénic model with a box length of 6.25 pc. Recently, Liu et al. (2021) has expanded the analysis to high-density ($\sim 10^4-10^6$ cm⁻³) trans-Alfvénic clumps/cores at 1–0.2 pc scales and obtained Q_c ∼ 0.28 for several strong-field models. Liu et al. (2021) also found Q_c ≪ 0.28 at $<$0.1 pc or n > 10⁶ cm⁻³ regions. Both Ostriker et al. (2001) and Liu et al. (2021) proposed that their correction factors are only valid when δϕ_obs < 25°. Applying the same δϕ_obs < 25° criteria to the results in Padoan et al. (2001), their correction factor changes to Q_c ∼ 0.31. From the three studies, it is very clear that there is a decreasing trend of Q_c with increasing density and decreasing scale in self-gravitating regions. This could be due to reduced turbulent correlation scales and more field tangling at higher densities, which leads to a more significant LOS signal averaging effect.

Heitsch et al. (2001) proposed a formula to estimate the geometric mean of $B_{3\text{D}}^{\text{u}}$ and $B_{3\text{D}}^{\text{tot}}$ ($\sqrt{B_{3\text{D}}^{\text{u}}{B}_{3\text{D}}^{\text{tot}}}$, see Table 3). However, their formula includes δ(tan ϕ_obs), which was found to have large scatters (Liu et al., 2021). Also, the term $\sqrt{B_{3\text{D}}^{\text{u}}{B}_{3\text{D}}^{\text{tot}}}$ does not have physical meaning. Thus, it is not suggested to use the geometric mean formula to estimate the magnetic field strength.

Falceta-Gonçalves et al. (2008) proposed a formula $B_{\text{pos}\Vert }^{\text{tot}}\sim \sqrt{{\mu }_0\rho }\frac{\delta v_{\text{los}}}{\tan \delta {\varphi }_{\text{obs}}}$ and derived the correction factor for the total magnetic field strength for the first time. They compared the estimated $B_{\text{pos}\Vert }^{\text{tot}}$ with the initial input $B_{3\text{D}}^{\text{u}}$ in their non-self-gravitating simulations and suggested that $B_{\text{pos}\Vert }^{\text{tot}}$ is slightly larger than $B_{3\text{D}}^{\text{u}}$ when the magnetic field is in the POS. They also found $B_{\text{pos}\Vert }^{\text{tot}}<B_{3\text{D}}^{\text{u}}$ for large inclination angles. It is unclear how accurate the estimated $B_{\text{pos}\Vert }^{\text{tot}}$ values are compared to the total field strengths in their simulation. Recently, Chen et al. (2022) adopted the same formula $\sqrt{{\mu }_0\rho }\frac{\delta v_{\text{los}}}{\tan \delta {\varphi }_{\text{obs}}}$ but suggested that it estimates the underlying field strength than the total field strength. They also suggested that $\sqrt{{\mu }_0\rho }\frac{\delta v_{\text{los}}}{\sin \delta {\varphi }_{\text{obs}}}$ gives a better estimation of the total field strength. They tested the formulas with self-gravitating and initially trans-Alfvénic $({\mathcal{M}}_{A0}=1.04-1.45)$ simulations. Instead of using the area-averaged parameters (ρ, δv_los, and δϕ_obs) for the calculation, Chen et al. (2022) suggested that the average of the local pseudo-field strength calculated using the local physical parameters gives a better estimation of the true field strength. They proposed that the local gas density can be estimated with the velocity fitting method (Storm et al., 2014) or the equilibrium method (Elmegreen and Elmegreen, 1978), which, in combination with the local velocity dispersion and angular difference measurements, gives accurate field strength estimations within a factor of 2 at 1–20 pc scales and n ∼ 10^1.5–10⁴ cm⁻³ when γ > 30°. They also found that the correction factor decreases with increasing density and smaller scales if the analytic gas density in the simulation is adopted in the field strength estimation, which agrees with the trend found in the numerical studies of the simplest DCF formula (see discussions above in the same section).

Liu et al. (2021) also investigated the correction factor for the total magnetic field strength. They found Q_c ∼ 0.62 for $B_{\text{pos}}^{\text{tot}}\sim \sqrt{{\mu }_0\rho }\frac{\delta v_{\text{los}}}{\delta {\varphi }_{\text{obs}}}$ and Q_c ∼ 0.53 for $B_{\text{pos}}^{\text{tot}}\sim \sqrt{{\mu }_0\rho }\frac{\delta v_{\text{los}}}{\delta (\sin {\varphi }_{\text{obs}})}$. Although $\sqrt{{\mu }_0\rho }\frac{\delta v_{\text{los}}}{\delta {\varphi }_{\text{obs}}}$ is often used to derive the underlying field strength, it is more likely correlated with the total field strength when the angular dispersion is large. Note that the correction factors in Liu et al. (2021) only apply to scales of clumps and cores when densities are greater than 10⁴ cm⁻³ and are not applicable at larger scales (e.g., ISM or clouds.). Also note that the simulations in Liu et al. (2021) have physical conditions similar to clustered star-forming regions. It is possible that the correction factor for isolated low-mass star formation regions could be larger than those reported by Liu et al. (2021) at the same scales due to lower densities. Liu et al. (2021) did not find significant difference among the correction factors for different inclination angles.

Skalidis and Tassis (2021) and Skalidis et al. (2021) proposed a new formula $B_{\text{pos}}^{\text{u}}=\sqrt{\frac{{\mu }_0\rho }{2\delta {\varphi }_{\text{obs}}}}\delta v_{\text{los}}$ and tested their formula with scale-free non-self-gravitating models. They claimed that their formula does not need a correction factor. However, several assumptions in the derivation of the ST21 formula are approximations (see Section 2.1) and there are also uncertainties in statistics of magnetic field position angles (see Section 2.2). Therefore correction factors should be still required for this method to compensate these uncertainties. More tests are needed to understand the validity of the ST21 formula under different physical conditions (e.g., at high-density self-gravitating medium).

Other than in situations when the ADF method is improperly applied to observations (e.g., obtaining negative ordered field contribution or only adopting $a_2^{\prime }{l}^2$ when l > l_d), the uncertainty of the ADF method may mainly come from the maximum derivable values of the integrated turbulent-to-ordered or -total field strength ratio (see Section 2.1.3.2), which underestimates the turbulent-to-ordered or -total field strength ratio and overestimates the field strength. Liu et al. (2021) has estimated the correction factor for the total field strength derived from the ADF method in trans-Alfvénic clumps/cores. They failed to derive the strength of the non-linearly ordered field in their simulations, thus they did not compare the ADF-derived ordered field strength with simulation values.

A recent and important modification of the DCF method is the DMA (Lazarian et al., 2022), which theoretically derives the formula $B_{\text{pos}}^{\text{u}}=\sqrt{{\mu }_0\rho }f\frac{{\tilde{D}}^v}{D^{\varphi }}$ in the non-self-gravitating regime. The analytical correction factor f is a function of γ, ${\mathcal{M}}_A$, and fraction of turbulence modes. Lazarian et al. (2022) has listed the asymptotic forms of f and the DMA formula in typical conditions of the ISM and molecular clouds (see their Table 3). They tested the DMA formula with a set of non-self-gravitating simulations and found the analytically and numerically derived correction factors agree well with each other in typical interstellar conditions. They suggested a pronounced dependence of f on the mean field inclination angle and fraction of turbulence modes in molecular clouds. However, both parameters are difficult to obtain observationally, which makes it difficult to apply the DMA to observational data. A further extension of the DMA by including self-gravity is essential to increase its accuracy in the determination of magnetic field strengths in self-gravitating molecular clouds.

2.4 Observational DCF estimations in star-forming regions

The original DCF method and its modified forms have been widely applied to observations of magnetic fields in star-forming regions to estimate the magnetic field strength. Statistical studies of DCF estimations are of significant value to extend our understanding of the role of magnetic fields in star formation (e.g., Pattle and Fissel, 2019; Myers and Basu, 2021). Recently, Liu et al. (2022) compiled all the previous DCF estimations published between 2000 and June 2021 from polarized dust emission observations within star-forming molecular clouds. Similarly, Pattle et al. (2022) made a compilation of all types of DCF measurements published between 2001 and May 2021. Here we briefly summarise the previous observational DCF estimations.

2.4.1 Comparing magnetic field with gravity

2.4.1.1 B − n relation

During the contraction of molecular clouds, the gravity can compress and amplify the magnetic field. The power-law index of the B − n relation (B ∝ n^j) characterise the dynamical importance of magnetic field during the cloud collapse. In the case of extremely weak magnetic fields where gas collapses isotropically due to magnetic freezing, there is relation B ∝ n^2/3 for a self-gravitating cloud core during its temporal evolution (Mestel, 1966). In such case, the radial component of the magnetic field also has the 2/3 scaling dependence on the gas density at any time snapshots, whereas the tangential component does not follow this scaling relation. For stronger fields, the density increases faster than the magnetic field due to ambipolar diffusion at higher densities, which results in shallower power-law slopes (e.g., j ≲ 0.5, Mouschovias and Ciolek, 1999).

However, the temporal B − n relation of a single cloud is not obtainable observationally due to the long evolutionary time scale. Studies of the spatial B − n relation for a single cloud (e.g., Li et al., 2015) are also rare. Instead, observational B − n studies usually focus on the spatial B − n relation for an ensemble of star-forming regions at different evolution stages and different scales. Crutcher et al. (2010) made a pioneering study of the spatial B − n relation based on the Bayesian analysis of a collection of Zeeman observations. They found that the magnetic field does not scale with density at n_{H_{2}} < 300 cm⁻³, but scales with density as B ∝ n^0.65 at n_{H_{2}} > 300 cm⁻³. Myers and Basu (2021) compiled the DCF estimation in 17 dense cores and reported B ∝ n^0.66. With compilations of DCF estimations, Liu et al. (2022) and Pattle et al. (2022) found a similar trend to the Zeeman results in that the magnetic field does not scale with density at lower densities, but scales with density at higher densities. Due to the large scatters and the uncertainty in correction factors, they did not report the critical density and magnetic field strength for the transition. Liu et al. (2022) reported B ∝ n^0.57 with a simple power-law fit for the high-density regime.

Despite the progress in the observational B − n studies, concerns have been raised on whether the B − n relation from a collection of different sources can be compared with model predictions for individual sources (Li, 2021). For the Zeeman observations, Tritsis et al. (2015) and Jiang et al. (2020) found that adopting different observational uncertainties of n other than the factor of 2 uncertainty adopted by Crutcher et al. (2010) can affect the fitted slope of the B − n relation, which questions the validity of the Bayesian analysis in Crutcher et al. (2010). Tritsis et al. (2015) further found that the samples collected in Crutcher et al. (2010) are preferentially non-spherical, which is inconsistent with the B ∝ n^2/3 scaling. The DCF-derived B − n relation is also very uncertain due to the scatters on the DCF estimations and the intrinsic B ∝ n^0.5 dependence of the DCF method. We do not aim to present a detailed review of the B − n relation, thus we refer readers to Tritsis et al. (2015), Crutcher and Kemball (2019), Hennebelle and Inutsuka (2019), Li (2021), and Pattle et al. (2022) for additional detailed discussions.

2.4.1.2 Mass-to-flux-ratio to critical value

The relative importance between the magnetic field and the gravity of individual sources is usually parameterized by the magnetic critical parameter λ (i.e., mass-to-flux-ratio in units of the critical value). The critical value of the mass-to-flux-ratio is model-dependent (e.g., Strittmatter, 1966; Mouschovias and Spitzer, 1976; Nakano and Nakamura, 1978). The magnetic critical parameter for an isothermal disk-like structure is given by (Nakano and Nakamura, 1978; Crutcher et al., 2004)

where G is the gravitational constant, μ_{H_{2}} is the mean molecular weight per hydrogen molecule, and m_H is the atomic mass of hydrogen. The magnetic critical parameter for a spherical structure with a density profile of n ∝ r⁻ⁱ is given by (Liu et al., 2022)

where k_i = (5 − 2i)/(3 − i). Equation 13 reduces to $\lambda \sim 8.7\times 10^{-21}\frac{N_{\text{H}2}/({\mathrm{c}\mathrm{m}}^{-2})}{B/(\mu G)}$ when i ∼1.83 (Liu et al., 2022). λ > 1 indicates the gravity dominates the magnetic field (i.e., magnetically super-critical), and vice versa. Alternatively, the importance between the magnetic field and the gravity can also be compared with the magnetic virial parameter α_B = 1/λ or with their energy ratios. The magnetic critical parameter or the magnetic virial parameter can also be expressed as a function of the number density n_{H_{2}} and radius r through the relation N_{H_{2}} = 4n_{H_{2}}r/3.

Statistical DCF studies have suggested that while molecular clouds are magnetically sub-critical (Planck Collaboration P. A. R. et al., 2016), the molecular dense cores within clouds are super-critical (Myers and Basu, 2021). The recent more complete DCF compilation by Liu et al. (2022) found a clear trend of increasing λ with increasing N_{H_{2}} (see Figure 3), where the average state transits from sub-critical to super-critical at $\sim 3.4\times 10^{21}$ cm⁻². This trend appears to agree with the prediction of magnetic field controlled star formation theories (Mouschovias et al., 2006), where magnetically sub-critical molecular clouds gradually form trans-to-super-critical substructures that collapse. The collapse is more dynamical in higher density regions. The dissipation of magnetic flux at higher densities may be due to ambipolar diffusion (Mouschovias and Ciolek, 1999) or magnetic reconnection (Lazarian and Vishniac, 1999). Mass accumulation along magnetic field lines can also increase the mass-to-flux ratio at higher densities. Despite the general trend seen in Figure 3, the samples collected by Liu et al. (2022) are mostly from different regions. Future multi-scale studies of the same region would be of great significance. High-mass star-forming regions tend to have higher column densities than low-mass star-forming regions at the same scales (Tan et al., 2014), thus high-mass star formation may be more magnetically super-critical than low-mass star formation within molecular clouds. There is only one DCF estimation of high-mass star-forming clouds at 10 pc scales so far (Orion, Planck Collaboration P. A. R. et al., 2016). More DCF estimations would be helpful for us to better understand the dynamical states of massive star formation at cloud scales. With better calibrated modified DCF methods (e.g., the extension of the DMA to the self-gravitating regime) and observational constraints for the physical parameters required for the DCF estimations, future magnetic criticality studies in molecular clouds could shed light on constraining the critical density in specific clouds and on comparing the mass-to-flux ratio of sources at different evolutionary stages.

FIGURE 3. The relation between the magnetic critical parameter (Equation 13) and column density from the previous DCF estimations (Liu et al., 2022). Different symbols represent different instruments (see Liu et al. (2022) and references therein). Different colors represent different methods. To avoid an assumption on whether the underlying or turbulent field is more dominant, the magnetic critical parameter is estimated with the 3D total magnetic field strength adopting the statistical relation $B_{3\text{D}}^{\text{tot}}=1.25B_{\text{pos}}^{\text{tot}}$ (i.e., $(\sqrt{\frac{3}{2}}+4/\pi )/2\sim 1.25$, see Section 2.2.5) instead of $B_{3\text{D}}^{\text{tot}}=\sqrt{\frac{3}{2}}{B}_{\text{pos}}^{\text{tot}}$ adopted in Liu et al. (2022).

2.4.2 Comparing magnetic field with turbulence

The relative importance between the underlying magnetic field and the turbulence of individual sources is usually parameterized by the Alfvénic Mach number

If there is an equipartition between the turbulent magnetic and kinetic energies, Equation 14 reduces to ${\mathcal{M}}_A=B_{3\text{D}}^{\text{t}}/B_{3D}^{\text{u}}$. The ratio $B_{3\text{D}}^{\text{t}}/B_{3D}^{\text{u}}$ can be derived with the statistics of the observed polarization angles if the mean field inclination angle and the anisotropy of the turbulent field are known (Lazarian et al., 2022). With the relations $B_{3\text{D}}^{\text{u}}=f_u{B}_{\text{pos}}^{\text{u}}$ and $B_{3\text{D}}^{\text{t}}=f_t{B}_{\text{pos}\perp }^{\text{t}}$, there is ${\mathcal{M}}_A\sim (f_t/f_u)B_{\text{pos}\perp }^{\text{t}}/B_{\text{pos}}^{\text{u}}$, where f_u and f_t are factors for the 3D to POS conversion. In small angle approximation, we obtain ${\mathcal{M}}_A\sim (f_t/f_u)\delta {\varphi }_{\text{obs}}$. The term $B_{\text{pos}\perp }^{\text{t}}/B_{\text{pos}}^{\text{u}}$ can also be derived from the ADF method (see Section 2.1.3). The relation between ${\mathcal{M}}_A$ and $B_{\text{pos}\perp }^{\text{t}}/B_{\text{pos}}^{\text{u}}$ or δϕ_obs can be regarded as an extension of the DCF formula, thus the correction factors for the DCF formula should also be applied to ${\mathcal{M}}_A$ under the same equipartition assumption. i.e., we have ${\mathcal{M}}_A\sim (f_t/f_u/Q_c){(B_{\text{pos}\perp }^{\text{t}}/B_{\text{pos}}^{\text{u}})}_{\text{adf}}$ or ${\mathcal{M}}_A\sim (f_t/f_u/Q_c)\delta {\varphi }_{\text{obs}}$. Adopting an additional correction factor f_o to account for the ordered field contribution to the angular dispersion, the relation between ${\mathcal{M}}_A$ and δϕ_obs becomes ${\mathcal{M}}_A\sim (f_t/f_u/Q_c/f_o)\delta {\varphi }_{\text{obs}}$.

Liu et al. (2022) has compiled the observed angular dispersion from previous DCF studies (see Figure 4). They suggested that the observed angular dispersion does not provide too much information on the Alfvénic states of molecular clouds due to the maximum angle limit of a random field (δϕ_obs < 52°) and the lack of an appropriate DCF correction factor at cloud scales. Without knowledge of the inclination angle and turbulence anisotropy for most observations in the compilation, they adopted the statistical correction factor f_u ∼ 4/π for a randomly distributed 3D mean field orientation (Crutcher et al., 2004), $f_t\sim \sqrt{3}$ for approximate isotropic turbulent fields in self-gravitating regions (Liu et al., 2021), and f_o ∼ 2.5 based on the numerical work by Liu et al. (2021). Additionally adopting the Q_c for different modified DCF approaches reported in Liu et al. (2021), Liu et al. (2022) found that the average ${\mathcal{M}}_A$ for the substructures in molecular clouds is ∼0.9, which indicates the average state is approximately trans-Alfvénic. They also suggested that both sub- and super-Alfvénic states exist for the cloud substructures and did not find strong relations between ${\mathcal{M}}_A$ and n. Pattle et al. (2022) did a similar study of ${\mathcal{M}}_A$ with their compilation of observed angular dispersions. Unlike Liu et al. (2022), Pattle et al. (2022) corrected the estimated ${\mathcal{M}}_A$ with the ratio between the DCF-derived POS field strength and Zeeman-derived LOS field strength at similar densities assuming that the Zeeman estimations are accurate references. They found that the turbulence is approximately trans-Alfvénic on average and ${\mathcal{M}}_A$ has no clear dependence on n, which agree with Liu et al. (2022). Note that both Liu et al. (2022) and Pattle et al. (2022) are statistical studies, where the adopted statistical relations may not apply to individual sources.

FIGURE 4. The relation between the directly estimated angular dispersion and the number density from the previous DCF estimations (Liu et al., 2022), reproduced with permission. Different symbols represent different instruments. Some data points with δϕ_obs > 52° are unphysical and are not plotted.

As the ADF method removes the contribution from the ordered field, the turbulent-to-ordered field strength ratio derived by the ADF method should be more suitable for the study of the Alfvénic state than the directly estimated angular dispersions. However, the applicability of the ADF method to determine the ${\mathcal{M}}_A$ value is limited by the maximum derivable turbulent-to-ordered field strength ratio (see Section 2.1.3.2), its uncertainty on the LOS signal integration (see Section 2.2.3), and the lack of appropriate numerically-derived correction factors for these uncertainties (see Section 2.3).

If an alternative assumption of the Fed16 equipartition (see Section 2.1.1) is adopted, the ${\mathcal{M}}_A$ should be correlated with $\sqrt{B_{\text{pos}\perp }^{\text{t}}/B_{\text{pos}}^{\text{u}}}$ or $\sqrt{\delta {\varphi }_{\text{obs}}}$ instead of $B_{\text{pos}\perp }^{\text{t}}/B_{\text{pos}}^{\text{u}}$ or δϕ_obs (Skalidis et al., 2021). However, the applicability of these alternative relations is limited by the lack of correction factors at different physical conditions.

In summary, the average state of star-forming substructures within molecular clouds may be approximately trans-Alfvénic, but the observed angular dispersions do not yield clues on the Alfvénic state of molecular clouds themselves. Note that the equipartition assumption (either the DCF53 or the Fed16 assumption), which should be independently confirmed, is a prerequisite for using the angular dispersion to determine the Alfvénic Mach number. If the equipartition assumption is not satisfied for some of the sources, the average state should be more super-Alfvénic.

2.4.3 Equilibrium state

The equilibrium state of a dense structure is usually parameterized by the virial parameter. Neglecting the surface kinetic pressure and the thermal pressure, the total virial parameter for a spherical structure considering the support from both the magnetic field and the turbulence is estimated as the ratio of the total virial mass and the gas mass

The total virial mass is given by (Liu et al., 2020)

where the magnetic virial mass is estimated with

and the turbulent virial mass is estimated with

Alternatively, the equilibrium state can also be studied by comparing 2E_turb + E_B with E_G, where E_turb is the turbulent kinetic energy, E_B is the magnetic energy, and E_G is the gravitational energy.

Figure 5 shows the total virial parameter as a function of column density for the dense substructures within molecular clouds based on the DCF compilation by Liu et al. (2022). Due to the lack of mass estimations, we do not show the α_turb+B for the molecular clouds observed by the Planck. Since the magnetic field can solely support the clouds (see Section 2.4.1.2), these clouds should have α_turb+B > 1 (i.e., super-virial). Low-mass star-forming regions with M ≤ M_crit = 870M_⊙(r/pc)^1.33 and high-mass star-forming regions with M > M_crit (Kauffmann et al., 2010) are indicated with different colors. In Figure 5, the high-mass regions with highest column densities (e.g., N_{H_{2}} > 10²⁴ cm⁻²) tend to be trans- or sub-virial, but both super- and sub-virial states exist at lower column densities. The median α_turb+B in low-mass and high-mass regions are ∼1.1 and ∼0.66, respectively, suggesting that the gravity may be more dominant in high-mass star formation. It may also indicate that high-mass star formation within molecular clouds tends to be more likely in non-equilibrium. It is possible that the magnetic field strength is overestimated for some sources due to the energy non-equipartition, which suggests a even more dynamical massive star formation. In summary, it appears that star-forming regions with higher column densities may have smaller total virial parameters due to more significant roles of gravity, but this trend is highly uncertain due to the large scatters.

FIGURE 5. The relation between the total virial parameter (Equation 15) and column density for the DCF sample catalogued by Liu et al. (2022). Different symbols represent different instruments. Red and blue colors indicate low-mass (LM) and high-mass (HM) star-forming regions, respectively.

3 The HRO analysis

3.1 Basics

Soler et al. (2013) developed the HRO analysis to characterize the relative orientation of the magnetic field with respect to the density structures, which can be used to establish a link between observational results and the physics in simulations. In 3D, the HRO is the relation between the 3D magnetic field orientation and the number density gradient ∇n. In the POS, the HRO is the relation between the POS magnetic field orientation and the column density gradient ∇N. The calculation of the density gradient is by applying Gaussian derivative kernels to the density structure (Soler et al., 2013). The density gradient at different scales can be studied by varying the size of the derivative kernels. The orientation of the iso-density structure is perpendicular to the direction of the density gradient. For observational data in the POS, the angle ϕ_B-N between the orientation of the iso-column density structure and the POS magnetic field orientation⁸ is given by

where E is the POS electric field pseudo-vector (Planck Collaboration P. A. R. et al., 2016). The angle ϕ_B-N is within ±90°. When ϕ_B-N = 0°, the magnetic field is parallel to the orientation of the column density structure and is perpendicular to the column density gradient. When |ϕ_B-N| = 90°, the magnetic field is perpendicular to the orientation of the column density structure and is parallel to the column density gradient. The preferential relative orientation between the magnetic field and the density structure is characterized by the histogram shape parameter (Planck Collaboration P. A. R. et al., 2016; Soler et al., 2017)

where A₀ is the percentage of pixels with |ϕ_B-N| < 22.5° and A₉₀ is the percentage of pixels with 67.5° < |ϕ_B-N| < 90° in a polarization map. Positive ξ values indicate the column density structure is more likely aligned with the magnetic field, and vice versa. The ξ in 3D can be derived similarly.

Using the parameter ξ to characterize the relative orientation has some drawbacks. For instance, the derivation of the parameter ξ completely ignores angles within 22.5° < |ϕ_B-N| < 67.5°. It also suffers from intrinsic deficiencies of binning of angles. To overcome these shortcomings, Jow et al. (2018) improved the HRO analysis with the projected Rayleigh statistic (PRS), which uses the PRS parameter Z_x instead of ξ to characterize the preferential relative orientation. For a set of n′ angles of relative orientation {ϕ_B−N,i}, Z_x is estimated with

Z_x > 0 indicates a preferential parallel alignment between the magnetic field and column density structure, and vice versa. Jow et al. (2018) suggested that the parameter Z_x is more statistically powerful than the parameter ξ, especially when the sample size is small or when the angles {ϕ_B−N,i} are more uniformly distributed. Equation 21 cannot be directly applied to 3D data. The formula for the 3D PRS parameter still needs more theoretical investigations (Barreto-Mota et al., 2021).

The VGT group (e.g., González-Casanova and Lazarian, 2017; Lazarian and Yuen, 2018) introduced an alignment measure (AM) parameter to study the relative orientation between magnetic field orientations and velocity gradients, where the AM can also be used to study the relative orientation between magnetic fields and density structures. The AM for ϕ_B-N can be expressed as

Similar to Z_x, the AM is also based on the Rayleigh statistic. AM is in the range of -1 to 1, which can be regarded as a normalized version of Z_x.

3.2 Observations

The relation between the cloud density structures and magnetic field orientations has been extensively studied observationally. For example, Li et al. (2013) found the orientation of Gould Belt clouds tends to be either parallel or perpendicular to mean orientation of the cloud magnetic field, which they interpreted as strong fields channeling sub-Alfvénic turbulence or guiding gravitational contraction. Toward the same sample, clouds elongated closer to the field orientation were found to have 1) higher star formation rates, which was suggested to be due to their smaller magnetic fluxes as well as weaker magnetic support against gravitational collapse (Li et al., 2017); 2) shallower mass cumulative function slopes (Law et al., 2020), i.e., shallower column density probability distribution functions (N-PDFs), or, in other words, more mass at high densities. In filamentary clouds, there is evidence that the magnetic field is parallel to the low-density striations and is perpendicular to the high-density main filament (e.g., Chapman et al., 2011; Cox et al., 2016), which implies that the main filament accrete gas through the striations along the field lines. Besides the success of those observational studies, the HRO analysis has enabled a way to perform pixel-by-pixel statistics for the local alignment between the column density structure and the magnetic field. Observational studies using the HRO analysis have been focused on studying this alignment at different column densities.

The first HRO analyses were made with observations from the Planck/HFI at large scales. With a smoothed resolution of 15′, Planck Collaboration R. et al. (2016) found that ξ is mostly positive and is anti-correlated with the column density⁹ N_{H_{2}} over the whole-sky at N_{H_{2}} ≲ 5 × 10²¹ cm⁻². The Planck observations toward 10 nearby (d < 450 pc) clouds (Planck Collaboration P. A. R. et al., 2016) with a smoothed resolution of 10′ (∼0.4–1.3 pc) have revealed a prevailing trend of decreasing ξ with increasing N_{H_{2}}, with ξ being positive at lower N_{H_{2}} and becoming negative at higher N_{H_{2}} in most clouds except for those with low column densities (e.g., CrA). The transition of ξ from positive to negative values was found to be at N_{H_{2},tr} ∼ 2.5 × 10²¹ cm⁻².

Subsequent studies have expanded the HRO analysis to compare the large-scale magnetic field observed with Planck/HFI or the BLASTPol with the smaller-scale column density structures revealed by the Herschel Space Observatory. Malinen et al. (2016) compared the Herschel dust emission structures at a 20′′ (∼0.01 pc) resolution and the large-scale magnetic field orientation revealed by Planck polarization maps at a 10′ (∼0.4 pc) resolution and found a trend of decreasing ξ with increasing N_{H_{2}} in the nearby cloud L1642, which transits at N_{H_{2},tr} ∼ 8 × 10²⁰ cm⁻². Soler et al. (2017) found the same decreasing ξ − N_{H_{2}} trend in the Vela C molecular complex by comparing the large-scale magnetic field orientation revealed by BLASTPol at a smoothed resolution of 3′ (∼0.61 pc) with the column density structures revealed by Herschel at a smoothed resolution of 35.2′′ (∼0.12 pc). The transition column density in Vela C is N_{H_{2},tr} ∼ 6 × 10²¹ cm⁻². They also found the slope of the ξ − N_{H_{2}} relation is sharper in sub-regions where the high-column density tails of N-PDFs are flatter. Soler (2019) compared the Herschel column density maps at a resolution of 36′′ (∼0.03–0.08 pc) with the large-scale magnetic field from Planck observations for the ten clouds studied by Planck Collaboration P. A. R. et al. (2016) and found the ξ (or Z_x) decreases with increasing N_{H_{2}} in most clouds, which is in agreement with the study by Planck Collaboration P. A. R. et al. (2016). They also found that regions with more negative cloud-averaged ξ (or Z_x) tend to have steeper N-PDF tails, but did not find a clear trend between the cloud-averaged ξ (or Z_x) and the star formation rate. In addition, Fissel et al. (2019) compared the magnetic field orientation revealed by BLASTPol with the integrated line-intensity structures of different molecular lines observed with Mopra in Vela C. They found that the line emission for low-density tracers are statistically more aligned with the magnetic field, while high-density tracers tend to be perpendicular to the magnetic field. The transition occurs at n_{H_{2},tr} ∼ 10³ cm⁻³.

At smaller scales, the polarization observations from the SOFIA/HAWC + or JCMT/POL2 can be used to study the relative orientation between the magnetic field and column density structures within elongated filamentary structures. Lee et al. (2021) compared magnetic field orientations from the HAWC + observations with the Herschel column density maps in Ophiuchus/L1688 at a smoothed resolution of 36.3′′ (∼0.02 pc). They found smaller ξ values at higher N_{H_{2}}, with ξ being mostly negative except for the first two column density bins. The transition column density in L1688 is N_{H_{2},tr} ∼ 5 × 10²¹ cm⁻². They also found that ξ increases from negative to ∼0 at higher column densities in the dense core Oph A within L1688, which suggests the magnetic alignment behavior is more complex at higher column densities. Li et al. (2022) compared the inferred magnetic field from HAWC+ with the Herschel column density map in Taurus/B211 and found that the magnetic field is more likely perpendicular to the column density structures in B211. Kwon et al. (2022) found negative ξ values in the Serpens Main cloud with JCMT observations at a resolution of 14′′ (∼0.03 pc), where they suggested that the first column density bin at N_{H_{2}} ∼ 9.3 × 10²¹ cm⁻² is approximately the transition column density. In Kwon et al. (2022), the ξ decreases from 0 to negative values from N_{H_{2}} ∼ 9.3 × 10²¹ cm⁻² to $\sim 4.6\times 10^{22}$ cm⁻², which mainly corresponds to the filamentary structures. Between N_{H_{2}} ∼ 4.6 × 10²² cm⁻² to $\sim 10^{23}$ cm⁻², ξ increases back to 0, suggesting that the magnetic fields are trending parallel to elongated structures. At N_{H_{2}} > × 10²³ cm⁻², ξ again decreases with increasing N_{H_{2}}. This behavior suggests a complex interplay between mass accumulation, hub-filament interaction, and gravitational collapse within filamentary structures.

Interferometric polarization observations are capable of revealing the magnetic field orientation within molecular dense cores. With ALMA observations at a resolution of 0.35′′ (∼140 AU), Hull et al. (2017b) found that the magnetic field is slightly more perpendicular than parallel to the column density structure around the class 0 protostellar source Ser-emb 8 in Serpens Main. Using ALMA observations at a resolution of 1.1′′ (∼200 AU), Hull et al. (2020) found that the magnetic field tends to be perpendicular to the column density structure in the class 0 protostellar source BHR 71 IRS1. They also found the magnetic field tends to be perpendicular to the column density structure but is parallel to outflow-cavity walls in another class 0 source BHR 71 IRS2, suggesting that the magnetic field is affected by the outflow activity. With ALMA observations at a resolution of 1.2′′ (∼0.02 pc) toward the high-mass star-forming region G327.3, Beuther et al. (2020) found that the relative orientation between the magnetic field and the dust emission structure changes from perpendicular to a random distribution when it is closer to the dust emission peak (i.e., ξ or Z_x increases with increasing N_{H_{2}}). This clearly contradicts the general decreasing ξ − N_{H_{2}} trend at larger scales found by studies of lower resolution data, and might be related to the rotational properties of potential disks.

In summary, there is a general trend that the relative orientation between the magnetic field and the column density structure changes from parallel to perpendicular with increasing gas column densities. The transition from parallel to perpendicular alignment occurs at N_{H_{2},tr} ∼ 10²¹–10²² cm⁻² and may be at n_{H_{2},tr} ∼ 10³ cm⁻³, which are comparable to the transition or critical densities of the B − n and B − N relations measured from DCF estimations (Liu et al., 2022) or Zeeman observations (Crutcher, 2012). Within filaments and molecular dense cores, the relative orientation could return to random alignment in some sub-regions due to the impact of the projection effect (see Section 3.3), accreting gas flows, outflows, and/or the disk rotation. More observations are needed to better understand the returning to random alignment at high densities.

3.3 Simulations

As indicated in Section 3.2, most clouds show a transition from a preferentially parallel alignment between the magnetic field orientation and density structures at low densities to a preferentially perpendicular alignment at higher densities. Numerical HRO studies have been focused on interpreting the contradictory alignment at different densities and understanding the implication of this transition.

Soler et al. (2013), Planck Collaboration P. A. R. et al. (2016), and Soler et al. (2017) have analysed 3 self-gravitating models with different initial magnetization levels (plasma β₀ = 100, 1, and 0.1). The 3 models are initially super-Alfvénic (${\mathcal{M}}_{A0}$ = 100, 10, and 3.16). At $\sim 2/3$ of the flow crossing time¹⁰, the 3 models become super-Alfvénic, trans-Alfvénic, and sub-Alfvénic on average (Planck Collaboration P. A. R. et al., 2016), respectively. The magnetic field in the 3 models are preferentially parallel to the density structures. Only in the trans- or sub-Alfvénic models, the relative orientation between the magnetic field and the density strucutre changes from parallel to perpendicular at high densities (Planck Collaboration P. A. R. et al., 2016). Soler et al. (2017) found that the transition occurs when the term A₂₃ changes its sign. The A₂₃ term is a coefficient in the time-varying equation describing the angle between the magnetic field and density gradient, which is a function of the velocity field, magnetic field, and density distribution. They also found that the local values of ${\mathcal{M}}_A>1$ and ∇v < 0 (i.e., converging gas flow) in high-density regions are associated with perpendicular alignment but are not sufficient for the relative orientation to become perpendicular.

Chen et al. (2016) analysed 3 self-gravitating models with different initial magnetic levels (${\mathcal{M}}_{A0}$ = 2.2, 4.4, and 8.8). The 3 initially super-Alfvénic models become trans-to-sub-Alfvénic (${\mathcal{M}}_A$ = 0.78, 0.81, and 0.84) on average due to shock compression when the most evolved core starts to collapse. They suggested that the change of relative orientation from parallel to perpendicular above a transition density happens when the sub-Alfvénic gas is gravitationally accelerated to become locally super-Alfvénic again in overdense regions. However, Otto et al. (2017) found that the change of relative orientation can still happen even when the high-density self-gravitating region is sub-Alfvénic (local ${\mathcal{M}}_A=0.54$) in their simulations. Hence ${\mathcal{M}}_A>1$ may be sufficient but is not necessary for the transition of relative orientation.

Seifried et al. (2020) performed a set of self-gravitating colliding flow simulations with ${\mathcal{M}}_{A0}$ = 0.6–4. At a time snapshot of 19 Myr, the mass-to-flux ratio μ of the simulations become 0.54–4.3, where the critical μ is 1. They found that only models with μ ≲1 shows a clear transition from parallel to perpendicular alignment for the relative orientation. Note that their models with μ ≲1 also have ${\mathcal{M}}_{A0}\lesssim 1$. Similar to Soler et al. (2017), they suggested that the observed transition can be well explained by the term A = A₁ + A₂₃, where A₁ is another coefficient describing the evolution of the angle between the magnetic field and density gradient. They found that ξ < 1 is associated with local ${\mathcal{M}}_A\gg 1$ in high density regions, which agrees with Soler et al. (2017). Seifried et al. (2020) also performed a set of SILCC-Zoom simulations with different physical conditions, but their output are less clear than those of the colliding flow simulations. Therefore no conclusions were drown from these SILCC-Zoom simulations.

Körtgen and Soler (2020) performed a set of initially supersonic (initial Mach number ${\mathcal{M}}_{s0}=7.5$) non-gravitating simulations with different initial magnetic levels (β₀ = 10, 1, and 0.01; ${\mathcal{M}}_{A0}$ = 16.7, 5.3, and 0.5), where the turbulence is driven either solenoidally or compressively. They found only models with β₀ = 0.01 (${\mathcal{M}}_{A0}$ = 0.5) and compressible turbulence clearly show the transition from parallel to perpendicular alignment at higher densities, which may suggest the main driver of the change of relative orientation is the compression of the gas. They also found the transition point for ξ ∼ 0 does not perfectly correspond to A = A₁ + A₂₃ ∼ 0, where A ∼ 0 happens at slightly lower densities than ξ ∼ 0. Their analysis suggests that A = A₁ + A₂₃ > 0 is not a sufficient condition and self-gravity is not a necessary condition for the transition of relative orientation. The change of relative orientation without self-gravity is also seen in the initially trans-to-sub-Alfvénic and very supersonic simulations by Hu et al. (2019a).

Girichidis (2021) also found a transition from parallel alignment to perpendicular alignment at higher densities in a 10 pc region of their large-scale simulation. They suggested that the transition happens when the local mass-to-flux ratio exceeds its critical value and the gravitational forces dominates the combination of thermal and magnetic pressures.

Barreto-Mota et al. (2021) studied the relative orientation between the magnetic field and density structure with simulations with a wide range of initial states (sub-Alfvénic and super-Alfvénic; with and without gravity). In non-gravitating simulations, only initially sub-Alfvénic models with the largest initial sonic number $({\mathcal{M}}_{s0}=7)$ show signs of perpendicular alignment between the magnetic field and column density structure in 2D at late stages of evolution. In simulations with gravity, the initial distribution of Z_x is similar to the simulations without gravity. At later stages, all gravitating sub-Alfvénic models show a change of Z_x from positive to negative values at high densities. They concluded that self-gravity may help to create structures perpendicular to the magnetic field.

Ibáñez-Mejía et al. (2022) analysed 3 clouds in a galactic-scale simulation with gravity. They found that the gas is trans-Alfvénic at n_{H_{2}} < 10² cm⁻³. At n_{H_{2}} > 10² cm⁻³, the gas becomes super-Alfvénic, which they suggested to be due to gravitational collapse. The relative orientation changes from predominantly parallel to much more random or perpendicular at even higher densities. i.e., they also found ${\mathcal{M}}_A>1$ at high density is associated with perpendicular alignment but is not a sufficient condition for the transition of relative orientation, which agrees with Soler et al. (2017) and Seifried et al. (2020).

Observations can only trace the projected quantities in the POS. The projection effect can make the 2D ξ − N_{H_{2}} relation different from the 3D ξ − n_{H_{2}} relation. Two parallel vectors in 3D are still parallel when projected in 2D, but 2 perpendicular vectors in 3D can have any relative orientations when projected in 2D depending on the viewing angle (e.g., Planck Collaboration R. et al., 2016). Soler et al. (2013) found that the relative orientation preserve well from 3D to 2D in sub-Alfvénic environments, except when the magnetic field orientation is close to the LOS. In contrast, the projected relative orientation does not necessarily reflect the relative orientation in 3D in super-Alfvénic environments. Other numerical studies found similar results for the projection effect (Seifried et al., 2020; Barreto-Mota et al., 2021; Girichidis, 2021).

The transition density at which the relative orientation changes its sign has also been studied by simulations. There is evidence that models with stronger magnetic field tend to have smaller transition densities in individual numerical studies (Soler et al., 2013; Chen et al., 2016), but this transition may not be solely dependent on the Alfvénic Mach number (Pattle et al., 2022). It should be noted that different numerical studies do not have a consistent way for measuring the local ${\mathcal{M}}_A$. Each study reports the local ${\mathcal{M}}_A$ averaged over arbitrary scales, which makes it challenging to compare their results.

Table 4 summarizes the simulations that show the transition of relative orientation from parallel to perpendicular alignment at higher densities. In most simulations (Soler et al., 2013; Planck Collaboration P. A. R. et al., 2016; Chen et al., 2016; Soler et al., 2017; Körtgen and Soler, 2020; Seifried et al., 2020; Barreto-Mota et al., 2021; Ibáñez-Mejía et al., 2022), the transition of relative orientation only occurs when the large-scale environment is trans- or sub-Alfvénic. For non-self-gravitating simulations, the transition only happens when the initial environment is sub-Alfvénic and supersonic (Körtgen and Soler, 2020; Barreto-Mota et al., 2021). Alternatively, Seifried et al. (2020) proposed that the transition only occurs in an magnetically sub-critical large-scale environment. It is still unclear what triggers the relative orientation to change from parallel to perpendicular. Plausible reasons include local ${\mathcal{M}}_A>1$¹¹, A₁ + A₂₃ > 0, and/or ∇v < 0 in high-density regions, which are associated with the perpendicular alignment but may not be sufficient conditions for the transition. Alternatively, local μ > 1 in high-density regions may also be responsible for the transition, where a dominant gravity can help to create perpendicular alignment between the magnetic field and density structure.

TABLE 4. Simulations that show transitions of relative orientation from parallel to perpendicular alignment at higher densities.

The analytical explanation for the parallel alignment between the magnetic field and density structure at low densities and the perpendicular alignment at high densities is still unclear. At low densities, the gravity is not dominant, so the alignment is mainly due to the interplay between the magnetic field and turbulence. Soler et al. (2013), Planck Collaboration R. et al. (2016), Planck Collaboration P. A. R. et al. (2016), and Soler et al. (2017) proposed that the initially super-Alfvénic compressive turbulent flows can stretch the magnetic field and density structures in the same direction due to magnetic freezing. However, the parallel alignment is also seen in initially sub-Alfvénic simulations (see Table 4). So an initially super-Alfvénic enviroment might not be a necessary condition for the parallel alignment. Alternatively, Chen et al. (2016) and Fissel et al. (2019) suggested that the anisotropic turbulent eddies in sub-Alfvénic gas naturally create elongated density structures parallel to the magnetic field. The anisotropic nature of sub-Alfvénic turbulence can explain the majority of the low-density gas with local ${\mathcal{M}}_A<1$ and ξ > 0, but cannot explain the faction of gas with local ${\mathcal{M}}_A>1$ and ξ > 0 in some simulations (Soler et al., 2017; Seifried et al., 2020; Ibáñez-Mejía et al., 2022). At high densities, the situation is more complicated as the gravity steps in. Soler et al. (2013) and Planck Collaboration P. A. R. et al. (2016) suggested that the gas compression along field lines in sub-Alfvénic turbulence can create density structures perpendicular to the magnetic field even without gravity. Alternatively, a magnetized gravitational collapse in the presence of a strong magnetic field would naturally cause an elongated density structure that is perpendicular to the magnetic field (Mouschovias, 1976a; Mouschovias, 1976b). However, both explanations are inconsistent with the results of some simulations where ${\mathcal{M}}_A>1$ is associated with the perpendicular alignment. Another alternative explanation is that the perpendicular alignment at high densities may be determined by the strong large-scale magnetic field at low densities rather than the small-scale magnetic field at high densities (Planck Collaboration R. et al., 2016).

4 The KTH method

4.1 Basics

Koch et al. (2012a) proposed the KTH method to determine the magnetic field strength with the relative orientations between the magnetic field, emission intensity gradient, and local gravity. This method is based on ideal MHD force equations, where the intensity gradient is assumed to trace the resulting direction of motions in the MHD equation (i.e., the inertial term, Koch et al., 2013). This method leads to maps of position-dependent magnetic field strengths and magnetic-to-gravitational force ratios.

Under the assumptions of negligible viscosity, infinite conductivity (ideal MHD), isotropic magnetic field pressure, small turbulent-to-ordered field strength ratio, smoothly and slowly varying field strength, stationarity, and that the intensity gradient indicates the resulting direction of motions, Koch et al. (2012a) considered the force balance among pressure, gravity, and the curvature term of magnetic field to derive the field strength

where ψ is the relative orientation between the gravity and the intensity gradient, α is the relative orientation between the magnetic field and the intensity structure¹², P is the hydrostatic dust pressure, ϕ_G is the local gravitational potential, and R_B is the local magnetic field line radius. Koch et al. (2012a) then introduced a parameter Σ_B to measure the local field significance. The field significance parameter Σ_B is given by the reformulation of Equation 23

which quantifies the relative importance of the magnetic force F_B = B²/(4πR_B) compared to the combination of the gravitational and pressure forces F_G + F_P = ∇P + ρ∇ϕ_G. Koch et al. (2012a) further demonstrated analytically that Eqs 23, 24 do not suffer too much from the geometry and projection effects. If local changes in temperatures and densities are small compared to gravity, the pressure terms (∇P and F_P) can be omitted from Eqs 23, 24. Thus, Σ_B can be used to quantitatively indicate whether the magnetic field in a region is strong enough to prevent gravitational collapse (Σ_B > 1) or not (Σ_B < 1). Later, Koch et al. (2012b) suggested that the global mass-to-flux ratio normalized to the critical value within a specific region can be estimated with

where they adopted the magnetic critical mass from Shu et al. (1999) in the derivation of λ_KTH.

Due the a series of assumptions, the KTH method is subjected to many uncertainties (Koch et al., 2012a). In contrast to the DCF method that requires the information on the turbulent field, the KTH method regards the turbulent field as insignificant and only requires the ordered field structure. The contribution from turbulent field may be removed by averaging several neighboring pixels to derive an averaged local ordered field curvature (e.g., applying the Pil15 or Pat17 technique mentioned in Section 2.2.1), but the KTH method may still fail when the magnetic field is dominated by the turbulent component. The KTH method also is not applicable in the extremely strong field case where the matter can move only along the field lines (i.e., ψ ∼ 0, α ∼ π/2, and R_B ∼∞). In some regions with strong rotation, the effect of rotation can lead to a change in the resulting direction of motion, which may be mitigated by adding a centrifugal term in the MHD equation. If the temperature variation is irregular and significant over a map, the intensity gradient should be replaced by the density gradient in the KTH analysis. The biggest uncertainty of the KTH method may come from the basic assumption that the intensity (or density) gradient traces the resulting motion of the MHD force equation. Koch et al. (2013) has investigated this assumption analytically and suggested that its validity relies on the difference between the velocity and density gradient directions. Although the numerical studies of the HRO method can be regarded as a partial test of the KTH method, the KTH method itself has not been fully compared with simulations yet. Further numerical investigations will be of significance to understand the uncertainty of this method in different physical conditions.

4.2 Observations

Most observational studies using the KTH method only partially applied this method and studied the distribution of the intermediate parameters. Here we summarise the observational studies of each parameter.

4.2.1 Magnetic field versus intensity gradient: δ

The angle δ = 90° − α is the relative orientation between the magnetic field and the intensity (or density) gradient. Basically, studying the distribution of δ is similar to studying the distribution of ϕ_B-N introduced by the HRO technique.

Similar to ϕ_B-N, the interpretation of the value and distribution of δ is still not yet well established. Koch et al. (2013) and Koch et al. (2014) suggested that a bimodal distributions of δ can be interpreted as a sign of collapse, where the angle δ measures how efficiently the magnetic field inhibits a collapse. Koch et al. (2013) and Koch et al. (2014) further proposed that the angle δ can used as a tracer of the evolution stage of star-forming regions. They suggested that δ is spatially randomly distributed due to the lack of a gravitational center in the early phase (Type I). In a later stage, elongated dust structures appear in star-forming regions. The magnetic field is parallel to the major axis of elongated dust structures that are created by large-scale flows and/or turbulence (Type IIB) or is perpendicular to the major axis of elongated dust structures when the gravity has just started to shape the field (Type IIA). Types IIB and IIA may further evolve into one and the same system (Type III) where a dominant gravity drags the magnetic field to form a radial, pinched, or hourglass shape. Koch et al. (2014) analysed a sample of 50 star-forming regions with CSO and SMA observations and found that the mean |δ| values for Types IIB, IIA, and III are 51° (CSO and SMA), 30° (CSO) to 34° (SMA), and 30° (SMA only), respectively. The categorization of different δ types is an empirical characterization. The correspondence between the δ types and different evolution stages of low-mass and high-mass star formation at different scales is still unclear and warrants further investigations.

Studies of individual star-forming regions have revealed different values and distributions of δ. Tang et al. (2013) found that the magnetic field is less correlated with the dust intensity gradient at larger scales revealed by JCMT or CSO, while the two angles are more correlated at smaller scales observed by SMA. This trend is consistent with the results from the observational HRO studies (Section 3.2). With SMA polarization observations toward the massive dense cores in the DR21 filament, Ching et al. (2017) found that the magnetic field and the intensity gradient are misaligned in rotation-like cores and are aligned in non-rotation-like cores, which suggests the magnetic field could be distorted by the rotation. Girart et al. (2013) found random distributions of δ with an average value of 40° in the massive dense core DR21(OH) with SMA observations. With ALMA polarization observations toward the massive clumps in the IRDC G28.34, Liu et al. (2020) found random distributions of δ and average δ values of 40° and 46° in the dense cores in two massive clumps MM1 and MM4, respectively. In another IRDC G14.225, Añez-López et al. (2020) found that the δ value in the Hub-N region is mostly small with their CSO observations. More observational studies of δ are still required to better understand its general distributions at different scales and evolution stages of star formation.

4.2.2 Magnetic field versus local gravity: ω

The angle ω measures the relative orientation between the magnetic field and the local gravity. The angle ω may be used as a measure for how effectively gravity can shape the magnetic field structure (Koch et al., 2012a), but further explanations of the distribution of ω are still yet to be established.

There are few observational studies of ω. Koch et al. (2012a) found an average ω of 13° in W51 e2, suggesting the magnetic field lines are mainly dragged by gravity in this core. With ALMA polarization observations, Koch et al. (2018) found average sin ω values of 0.4, 0.41, and 0.47 (i.e., ω ∼ 23°, 23°, and 27°) in three massive cores W51 e2, W51 e8, and W51N, respectively. The distributions of ω in three cores vary, but all show some magnetic channels with sin ω ∼ 0. Koch et al. (2018) suggested that the gas collapse can proceed in free-fall in these magnetic channels without magnetic resistance, which may have some interesting implication for the star formation rate. Liu et al. (2020) found average ω values of 34° and 36° in the dense cores in massive clumps G28-MM1 and G28-MM4, respectively, which are lower than the average δ values in the same region. This shows that the magnetic field is more aligned with the gravity direction than the intensity (or density) gradient, which suggests ω may be better than δ (or ϕ_B-N of the HRO technique) for the study of the correlation between the magnetic field and gravity in high density regions. An interesting future study point of ω may be studying its variation at different densities in a similar way to the studies of the angle ϕ_B-N.

4.2.3 Intensity gradient versus local gravity: ψ

The angle ψ measures the relative orientation between the intensity (or density) gradient and the local gravity. The angle ψ may indicate how effective gravity can shape the density structure, but further analytical implications of ψ are still yet to be investigated.

There are few observational studies of the angle ψ. Koch et al. (2012a) found an average ψ of 20° in W51 e2. Liu et al. (2020) found average ψ values of 30°, 22°, and 28° in the dense cores in massive clumps G28-MM1, G28-MM4, and G28-MM9, respectively. These observations suggest that the gravity direction is closely aligned with the intensity gradient in high density molecular cores.

4.2.4 Magnetic force versus gravitational force: Σ_B

The field significance parameter Σ_B measures the ratio between the magnetic force and the gravitational force if the gas pressure is negligible. The implication of Σ_B is clear, but the accuracy of the estimated force ratio is not tested by simulations yet.

Based on a sample of 50 sources, Koch et al. (2014) found that the different types of star-forming regions categorized by the angle δ (see Section 4.2.1) show clear segregation in Σ_B values. Type IIB sources where the magnetic field is aligned with the clump/core major axis have an average Σ_B of 1.29 (CSO) to 1.49 (SMA), which suggests that type IIB sources are supported by the magnetic field and do not collapse on average. Type IIA sources where the magnetic field is perpendicular with the clump/core major axis have an average Σ_B of 0.69 (CSO) to 0.74 (SMA), which suggests that type IIA sources are collapsing on average. Type III sources at a later stage have a smaller average Σ_B value of 0.59 (SMA only), suggesting an even more dynamical collapsing.

Individual studies have mostly found Σ_B ≲ 1 in star-forming dense clumps/cores (e.g., W51 e2, W51A, W51N, DR21(OH), clumps in G34.43, cores in G28.34, and clumps in G14.225, Koch et al., 2012a; Koch et al., 2012b; Tang et al., 2013; Girart et al., 2013; Tang et al., 2019; Liu et al., 2020; Añez-López et al., 2020). Specifically, Tang et al. (2013) found that the lower density structures in W51N revealed by the CSO and JCMT have Σ_B ∼ 0.71–1.17, while the higher density structures revealed by SMA have Σ_B ∼ 0.5. Koch et al. (2012a) and Koch et al. (2012b) found smaller Σ_B values in higher column density regions in W51 e2 and W51A. These findings indicate that the gravity plays an more dominant role in higher density regions, which agrees with the observational DCF studies (see Section 2.4.1.2). On the other hand, Liu et al. (2020) found that the Σ_B values are higher in more evolved dense cores in the IRDC G28.34, which might suggest a more dynamical star formation at earlier evolution stages.

4.2.5 Magnetic field strength

The KTH method can be used to estimate the magnetic field strength with Equation 23. The accuracy of the estimated field strength is still unclear as it has not been tested by simulations yet.

By far, there are only two observational studies that has applied the KTH method to estimate the field strength. Koch et al. (2012a) found an average field strength of 7.7 mG in W51 e2. Stephens et al. (2013) estimated a field strength of 3.4 mG with the KTH method in L1157-mm, which is ∼2.5 times higher than the DCF estimation of 1.4 mG in the same work.

4.2.6 Mass-to-flux-ratio to critical value

The mass-to-flux-ratio to critical value λ_KTH can be estimated from the field significance parameter Σ_B through Equation 25. Similarly, the accuracy of λ_KTH is yet to be investigated by simulations.

Koch et al. (2012b) found larger λ_KTH values in higher density regions in W51 e2, which agrees with the trend from DCF estimations (see Section 2.4.1.2). Other two observational studies in DR21(OH) and G14.225N have found that the λ_KTH value estimated with the KTH method approximately agrees with the value of DCF estimations (Girart et al., 2013; Añez-López et al., 2020).

5 Summary

The recent improvement of instrumental sensitivity and development of new techniques (e.g., the VGT) have led to an increasing number of observations that reveal the POS component of magnetic field orientations in star-forming molecular clouds. In this review, we discuss the developments and limitations of the DCF and KTH methods that quantify the dynamic importance of magnetic fields in star-forming molecular clouds based on the field orientations and the HRO analysis that characterize the statistical relation between field orientations and density structures. We also summarize the observational studies using these methods and discuss their implications on star formation.

The original DCF method is based on several assumptions: the total magnetic field is composed of a mean field component and a turbulent field component; the energy equipartition; isotropic turbulence; and the turbulent-to-mean or -total field strength ratio is traced by angular dispersions. The ordered field component is considered instead of the mean field component (e.g., the ADF method) if there are curved ordered magnetic field structures. We suggest that the ordered field and turbulent field of a particular region are local properties and are dependent on the scale range (i.e., the beam resolution to the maximum recoverable scale of interferometric observations or the region size) of the region of interest. There is still a debate on whether there is an equipartition between the turbulent magnetic field and the turbulent kinetic field or between the coupling-term magnetic field and the turbulent kinetic field in the sub-Alfvénic regime, while both equipartitions are not satisfied for super-Alfvénic turbulence. The energy non-equipartition can be the biggest uncertainty in the DCF method, which should be further investigated with simulations and observations. The uncertainty from anisotropic turbulence may be insignificant for the DCF estimations in self-gravitating regions. The turbulent-to-underlying or -total field strength ratio can be expressed as different forms of angular dispersions, but each has its limitations. The ADF method correctly accounts for the beam-smoothing effect, interferometric filtering effect, and the ordered field structure, but its applicability for quantifying the turbulent field and LOS signal integration may need further numerical investigations. The correction factor for the most widely used formula $B_{\text{pos}}^{\text{u}}\sim \sqrt{{\mu }_0\rho }\frac{\delta v_{\text{los}}}{\delta {\varphi }_{\text{obs}}}$ decreases with increasing density. The DMA presents an important improvement for the DCF method by analytically accounting for the mean field inclination angle and the turbulence anisotropy in the non-self-gravitating regime. A further extension of the DMA in the self-gravitating regime would be of significance.

A compilation of previous DCF estimations suggests a scenario that magnetically sub-critical low-density clouds gradually form super-critical high-density substructures. The critical column density is around 3.4 × 10²¹ cm⁻² on average, which needs to be better constrained and may differ in different clouds. The gravity may be more dominant in high-mass star formation than low-mass star formation. The average state of dense substructures within molecular clouds is approximately trans-Alfvénic if the energy equipartition assumption is satisfied, or super-Alfvénic if the energy equipartition assumption is unsatisfied for some of the sources.

Observational HRO studies mainly focus on the alignment between the magnetic field and density structure. Low-resolution HRO studies have found a general trend of transition from a preferentially parallel alignment at low column densities to a perpendicular alignment at higher column densities. This observational trend agrees with trans-to-sub-Alfvénic simulations, which indicates that the star-forming molecular clouds are trans-to-sub-Alfvénic. This trans-to-sub-Alfvénic state is consistent with the results derived from other techniques (e.g., the VGT, Hu et al., 2019b). The analytical explanation for the transition from parallel to perpendicular alignment is still unclear, but may be related to changes of the local Alfvénic Mach number, A₁ + A₂₃ term, mass-to-flux-ratio, and/or ∇v. The transition occurs at 10²¹–10²² cm⁻², which agrees with the critical column density derived from DCF estimations. But it is unclear whether the two transition column densities are related. High-resolution HRO studies have revealed a possible transition from perpendicular alignment back to random alignment at high column density sub-regions. The reason for this reverse transition is also unclear, but may be related to the impact of accretion gas flows, outflows, disk rotation, and/or the projection effect.

The advantage of the KTH method compared to the DCF method is that it does not require the information on the velocity dispersion. However, the uncertainty of the KTH method is still unknown since it has not been fully tested by simulations. Results from observational KTH studies on the relative alignment between the magnetic field and intensity (density) gradient within dense clumps/cores approximately agree with those of the observational HRO studies. The value and density-varying trend of the mass-to-flux-ratio and the magnetic field strength derived from the KTH method approximately agree with those derived from the DCF estimations.

Author contributions

JL, QZ, and KQ contributed to the outline of the review. JL led the writing of the manuscript. QZ and KQ read, commented on, and edited the manuscript.

Funding

JL acknowledges the support from the EAO Fellowship Program under the umbrella of the East Asia Core Observatories Association. KQ is supported by National Key R&D Program of China grant No. 2017YFA0402600. KQ acknowledges the support from National Natural Science Foundation of China (NSFC) through grant Nos. U1731237, 11590781, and 11629302.

Acknowledgments

We thank the referees for constructive comments that improved the clarity of this paper. JL thanks Prof. Martin Houde for insightful discussions about the correlation between turbulent magnetic and kinetic fields. JL thanks Heshou Zhang and Suoqing Ji for helpful discussions on the general concept of MHD turbulence.

Conflict of interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher’s note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Footnotes

¹In this paper, the mean field refers to the vector-averaged magnetic field (B^m) and the ordered field refers to the curved large-scale ordered magnetic field (B^o). We also use the term “underlying field” (B^u) to refer to either the mean field or ordered field since many previous studies did not explicitly differ the two. The ordered field and the mean field are equivalent if the large-scale ordered field lines are straight.

²The equations in SI units in this paper can be transformed to Gaussian units by replacing μ₀ with 4π.

³The inclination angle γ corresponds to the angle between the 3D mean field and the LOS throughout this paper.

⁴The field angle considers the direction of the magnetic field and is in the range of -180°–180° (Li et al., 2022), while the position angle ϕ only considers the orientation of the magnetic field and is in the range of -90°–90°. Due to the 180° ambiguity of dust polarization observations, only the magnetic field position angle is observable.

⁵The Hil09 variant of ADF directly estimates the ${(\frac{B_{\text{pos}\perp }^{\text{t}}}{B_{\text{pos}}^{\text{o}}})}_{\text{adf,int}}$. For the Hou09 and Hou16 variants, the ${(\frac{B_{\text{pos}\perp }^{\text{t}}}{B_{\text{pos}}^{\text{o}}})}_{\text{adf,int}}$ is derived by dividing the estimated ${(\frac{B_{\text{pos}\perp }^{\text{t}}}{B_{\text{pos}}^{\text{o}}})}_{\text{adf}}$ by a factor of $\sqrt{N_{\text{adf}}}$, where N_adf is the number of turbulent cells contained in the column of dust probed by the telescope beam.

⁶Note that there is a factor of $\sqrt{3}$ difference between the correlation length of the autocorrelation funciton in Houde et al. (2009) and the correlation length of a Kolmogorov-like turbulent power spectrum.

⁷Here we use f for the analytically derived correction factor to differ with the numerically derived correction factor Q_c.

⁸Note that there is a 90° difference between the ϕ_B-N defined in Soler et al. (2013) and in the subsequent HRO papers. In Soler et al. (2013), ϕ_B-N is defined as the angle between the column density gradient and the POS magnetic field orientation.

⁹A significant amount of the diffuse ISM is in the atomic phase, while molecular clouds are in the molecular phase. Here we use N_{H_{2}} to describe the column density in diffuse ISM and molecular clouds for uniformity.

¹⁰Planck Collaboration et al. (2016b) did not explicitly indicate the time of the snapshot they studied. We consulted their corresponding author to confirm this information.

¹¹It is unclear whether these numerical studies has excluded the non-turbulent motions (e.g., infall, outflow, and/or rotation.) in the calculation of ${\mathcal{M}}_A$, so the ${\mathcal{M}}_A>1$ at high densities found in these simulations may be arguable.

¹²The angle α is equivalent to the ϕ_B-N in Section 3 if the intensity structure perfectly traces the density structure.

References