Probing 3D Magnetic Fields Using Thermal Dust Polarization and Grain Alignment Theory

Magnetic fields are ubiquitous in the Universe and are thought to play an important role in various astrophysical processes. Polarization of thermal emission from dust grains aligned with the magnetic field is widely used to measure the 2D magnetic field projected onto the plane of the sky, but its component along the line of sight is not yet constrained. Here, we introduce a new method to infer 3D magnetic fields using thermal dust polarization and grain alignment physics. We first develop a physical model of thermal dust polarization using the modern grain alignment theory based on the magnetically enhanced radiative torque alignment theory. We then test this model with synthetic observations of magnetohydrodynamic simulations of a filamentary cloud with our updated POLARIS code. Combining the tested physical polarization model with synthetic polarization, we show that the B-field inclination angles can be accurately constrained by the polarization degree from synthetic observations. Compared to the true 3D magnetic fields, our method based on grain alignment physics is more accurate than the previous methods that assume uniform grain alignment. This new technique paves the way for tracing 3D B-fields using thermal dust polarization and grain alignment theory and for constraining dust properties and grain alignment physics.


INTRODUCTION
Magnetic fields (B-fields) are ubiquitous in the universe and thought to play an important role in astrophysical processes, including the evolution of the interstellar medium (ISM), formation and evolution of molecular clouds and filaments, star and planet formation (see Crutcher 2010; Hennebelle & Inutsuka 2019 for reviews), and cosmic ray transport.Polarization of starlight (Hall 1949;Hiltner 1949) and polarized thermal emission (Hildebrand 1988;Planck Collaboration et al. 2015b) induced by aligned dust grains (hereafter dust polarization) is the most popular tool to probe the projected B-fields on the plane of the sky (POS), i.e., two-dimensional (2D) B-fields (Hull & Zhang 2019;Pattle et al. 2022).The polarization angle dispersion is routinely used to measure the strength of the projected B-fields, B POS , using the Davis-Chandrashekhar-Fermi (DCF) technique (Davis 1951;Chandrasekhar & Fermi 1953).
Significant advances in polarimetric techniques from ground-based single-dish telescopes (JCMT, Wardthiemhoang@kasi.re.kr Thompson et al. 2017) and interferometers, including SMA and ALMA (Hull & Zhang 2019), space telescope (Planck, Planck Collaboration et al. 2015b), and airborne telescope (SOFIA/HAWC+, Dowell et al. 2010) provide a vast amount of dust polarization data, which allow us to measure the strength of 2D B-fields from a wide range of astrophysical scales, from the full galaxies to molecular clouds and filaments to prestellar cores, protostellar envelopes and disks (Pattle et al. 2022;Tsukamoto et al. 2022).To achieve more accurate measurements of B-field strengths using dust polarization, significant efforts have been invested in improving the DCF technique over the last decades (e.g., Ostriker et al. 2001;Falceta-Gonçalves et al. 2008;Hough et al. 2008;Cho & Yoo 2016;Lazarian et al. 2022.Yet, the measurements based on the DCF technique are still for the 2D B-fields only.Although the Zeeman splitting effect has been used to accurately measure the LOS B-fields using several spectral lines (H, OH, CN, CCS), Zeeman measurements are still seldom due to the weak signal and long-observing time (Crutcher 2010).
Astrophysical magnetic fields are three-dimensional, and thus accurate measurements of 3D B-fields are necessary to understand the dynamical role of B-fields in the ISM evolution and star formation (Hull & Zhang 2019;Pattle et al. 2022), especially in the formation and evolution of interstellar filaments (Hennebelle & Inutsuka 2019) and prestellar core collapse.However, to date, accurate measurements of interstellar 3D B-fields are not yet available.Over the last few years, various techniques have been proposed to constrain 3D B-fields by combining the different tracers, e.g., the combination of dust polarization with Faraday rotation (Tahani et al. 2018) (see Tahani 2022 for a review), the combination of dust polarization with velocity gradient (Lazarian et al. 2022;Hu & Lazarian 2023a,b).
Dust polarization itself has the potential of probing 3D B-fields because the polarization angle reveals the POS component and the polarization degree constrains the inclination angle of B-fields.Indeed, Lee & Draine (1985) first presented a model for starlight polarization induced by partially aligned grains, which depends on dust properties (shape, size, and composition), degree of grain alignment with the B-field, and the inclination angle of the mean B-fields and the fluctuations of the local field with respect to the mean B-fields.Therefore, in principle, accurate measurements of the dust polarization degree can constrain the B-field's inclination angle and the tangling when dust properties and grain alignment are well constrained.Using the observed polarization degree of silicate feature at 9.7 µm, Lee & Draine (1985) suggested that the B-fields along the LOS toward the Backlin-Neugebauer (BN) object are inclined by an angle sin 2 γ max < 0.5 and is highly ordered.Since then, there has been little progress in this direction because there was a lack of a quantitative theory of grain alignment that can accurately predict the grain alignment degree.
The last two decades have witnessed significant progress in grain alignment physics (e.g., Hoang 2022;Hoang et al. 2022, see Tram & Hoang 2022 for the latest reviews).The leading theory of grain alignment is based on radiative torques (RATs, Dolginov & Mitrofanov 1976;Draine 1996;Lazarian & Hoang 2007b;Hoang & Lazarian 2008).The RAT theory is qualitatively tested with various observations (Andersson et al. 2015;Lazarian et al. 2015).The modern understanding of grain alignment based on RATs established that grain alignment degree depends on the fraction of grains aligned at an attractor point with angular momentum greater than the thermal value (i.e., high-J attractors), denoted by f high−J .For the RAT theory, the fraction f high−J depends on the grain shape, the radiation field, the relative angle between the radiation direction and the magnetic field (Lazarian & Hoang 2007a;Hoang & Lazarian 2009).It can span between 0.2 − 0.5 (Herranen et al. 2021).The RAT alignment theory was implemented in a po-larized radiative transfer code, POLARIS (Reissl et al. 2016), where the value of f high−J is a fixed input parameter.The POLARIS code can predict thermal dust polarization maps and polarization spectra for different astrophysical environments.
As a first step toward an ab initio model of dust polarization, Lee et al. (2020) presented a physical model of thermal dust polarization using the RAT alignment theory, assuming the fixed f high−J and uniform B-fields.However, for grains with embedded iron inclusions, the value f high−J is increased due to the joint effect of both RATs and magnetic relaxation (Davis & Greenstein 1951;Jones & Spitzer 1967), aka magnetically enhanced RAT alignment or MRAT (Hoang & Lazarian 2008, 2016a).The enhanced magnetic relaxation is caused by the enhanced magnetic susceptibility of grains due to the inclusion of iron clusters into dust grains, which is most likely because Fe is among the abundant elements in the universe, and up to 95% of Fe is missing from the gas (Jenkins 2009;Dwek 2016).For grains with embedded iron inclusions, numerical simulations reveal the increase of grain alignment efficiency with magnetic relaxation rate (Hoang & Lazarian 2016a;Lazarian & Hoang 2019).The MRAT alignment mechanism can now predict the degree of grain alignment as a function of the local gas density, radiation field, grain sizes, and dust magnetic properties (Hoang & Lazarian 2008, 2016a;Lazarian & Hoang 2019).Therefore, it is necessary to develop a physical model for thermal dust polarization based on the MRAT mechanism, which can be used with the observed polarization degree to infer the B-field inclination angle and then 3D B-fields.
Recently, Chen et al. (2019) proposed a technique to infer 3D B-fields in molecular clouds using polarization of thermal dust emission.However, the authors assumed a simple model of thermal dust polarization in which the grain alignment efficiency and the B-fields do not change within the cloud.King et al. (2019) studied the effect of grain alignment on synthetic polarization and found it important, but their approximated alignment degree is based on the early model of RAT alignment by Cho & Lazarian (2005) that assumed that all grains larger than a critical size for RAT alignment are perfectly aligned (i.e., ideal RAT theory).
Moreover, the disregard of the B-field fluctuations along the LOS in Chen et al. (2019) is only valid for the case of strongly magnetized medium (i.e., sub-Alfvén turbulence, M A = δv/v A = δB/B < 1) where δv = v l for the turbulence velocity.However, polarimetric observations toward molecular clouds and dense clouds/cores reveal that the Alfvén turbulence is usually sub-Alfvénic , Pattle et al. 2021;Hoang et al. 2022;Ngoc et al. 2023;Karoly et al. 2023;Tram et al. 2023).However, the turbulence becomes super-Alfvénic with M A > 1 in high-density regions of n H > 10 7 cm −3 (see Fig 2f in Pattle et al. 2022).Therefore, in the majority of star-forming regions (from clouds to cores to protostellar envelope of n H < 10 7 cm −3 , observations reveal sub-Alfvénic turbulence, i.e., the B-field fluctuations are less important due to weak Alfvénic turbulence and strong B-fields, and the polarization would depend crucially only on the inclination angle and the grain alignment efficiency.However, a different compilation of B-field strengths by Liu et al. (2022) suggests trans-Alfvénic turbulence.Hu & Lazarian (2022, 2023a) improved Chen's method by considering the fluctuation of the B-field along the LOS predicted by anisotropic MHD turbulence.Using the perfect alignment and MHD simulations, the authors showed that Chen's method is valid for sub-Alfvénic turbulence of M A ≪ 1.Moreover, for 0.4 < M A < 1, the deviation of Chen's result and the actual value is above 10 degrees.However, as Chen et al., Hu & Lazarian (2022, 2023a) assumed that the grain alignment degree is constant within the molecular cloud.Therefore, the variation of grain alignment must be taken into account to reliably constrain the 3D B-field based on the polarization degree.
Here, we introduce a general technique by including the alignment efficiency using the modern theory of the MRAT grain alignment.For simplicity, we still assume the alignment efficiency is constant along the LOS, but we take into account the variation of the grain alignment efficiency across the cloud surface.By doing so, we can account for both the effect of and the grain alignment variation and magnetic turbulence across the cloud.Ultimately, one must also consider the variation of the grain alignment efficiency along the LOS, but it is challenging because of the fluctuations of both grain alignment and the magnetic turbulence.
The structure of our paper is as follows.In Section 2 we introduce a physical model of thermal dust polarization based on the MRAT mechanism and present a new method for constraining 3D B-fields using dust polarization.In Section 3 we will test the dust polarization model using synthetic polarization observations of MHD simulations.In Section 4 we will present our numerical results, compare them with the polarization from the analytical model, and apply our method for inferring the inclination angles using the synthetic polarization degree map.The effects of magnetic turbulence on the dust polarization model and inclination angles are discussed in Section 5.An extended discussion of our results and implications is presented in Section 6, and a summary is shown in Section 7.

A Physical Model of Thermal Dust Polarization
We first construct a physical model of thermal dust polarization using the theory of grain alignment based on radiative torques and magnetic relaxation (Hoang & Lazarian 2016a).

Grain Extinction and Polarization Cross-Section
Consider an oblate spheroidal grain shape with the symmetry axis â1 , which is also the grain axis of maximum inertia moment.The effective size of the grain of volume V gr is a = (3V gr /(4π)) 1/3 .For simplicity, we consider the composite dust model in which silicate and carbonaceous material are mixed in a single dust population, which is considered the leading dust model of the ISM (Draine & Hensley 2021).
The extinction cross-sections of the spheroidal grain for the electric field of the incident light parallel/perpendicular to the symmetry axis are denoted by , respectively.The polarization cross-section of the oblate grain is given by C pol = C ⊥ − C ∥ .The total extinction cross-section, obtained by averaging the extinction cross-section over an ensemble of grains with random orientation, is approximately given by C ext = (2C ⊥ + C ∥ )/3 (see e.g., Hoang et al. 2013;Draine & Hensley 2021).
Let dn/da be the grain size distribution.Then, n d (a) = n −1 H dn/da with n H being the hydrogen number density is the density of grains with a size in the range of a, a + da per H. Here, we assume the powerlaw size distribution of n −1 H dn/da = Ca −3.5 with C the normalization constant (Mathis et al. 1977).
The dust extinction and polarization coefficients are calculated as where a min and a max are the minimum and maximum grain sizes, respectively.

Thermal Dust Polarization from Aligned Grains
The modern grain alignment theory establishes that interstellar grains are aligned with their axis of maximum inertia moment, also the grain's shortest axis (â 1 ) along the angular momentum J due to various internal relaxation processes, including Barnett relaxation (Purcell 1979), inelastic relaxation (Lazarian & Efroimsky 1999;Efroimsky & Lazarian 2000) and nuclear relaxation (Lazarian & Draine 1999).The grain angular momentum J rapidly precesses around the ambient B-field B due to fast Larmor precession which determines the axis of grain alignment in the majority of astrophysical environments (see Hoang 2022 for details).Let ξ be the alignment angle between J and B (see Figure 1).The net alignment degree of the grain's axis of maximum inertia moment with B is described by f align (a).
Let ẑ be the line of sight (propagation direction of light).The observer's coordinate system is described by xŷẑ with xŷ representing the POS.For convenience, we assume that the mean B-field B 0 lies in the plane xẑ and makes an inclination angle γ with respect to ẑ, and the component projected onto the POS, B pos , is directed along the x axis (see Figure 1).Let C x , C y be the extinction cross-sections for the incident electric field E along x and ŷ, respectively.
The polarization and extinction cross-section averaged over the grain's Larmor precession for the light propagating along the ẑ are given by ( Lee & Draine 1985;Hoang et al. 2013;Planck Collaboration et al. 2015a) For the sake of simplicity, we first assume that grain alignment, dust properties, and B-fields do not change along the LOS.As a result, within the modified picketfence approximation (MPFA) model of polarization (Dyck & Beichman 1974;Draine & Hensley 2021), the intensity of polarized and total thermal emission for a composite dust model of mixed silicate and carbonaceous material is where T d is the mean dust temperature 1 , and N H is the hydrogen gas column density along the LOS (see also Appendix A).
In the electric dipole limit of λ/2πa ≫ 1, the cross-section is proportional to the grain volume, i.e., C pol,ext ∝ V gr (a) = 4πa 3 /3.Let C pol (a) = η pol V gr (a), 1 Here, for the physical model, we ignore the dependence of the grain temperature on the grain size and assume all grain sizes have the same mean temperature.then, the polarized intensity can be written as dan d (a)(4πa 3 /3)f align sin 2 γ, (7) where the mass-weighted alignment efficiency is given by and V d = (4πa 3 /3)n d (a)da is the total dust volume per H respectively, and η pol = σ pol /V d .Similarly, the emission intensity can be rewritten as From Equations ( 8) and (10), we calculate the degree of thermal dust polarization: where the extinction and polarization coefficients are given by Equations (2), and p i = σ pol /σ ext is the intrinsic polarization degree which only depends on the grain properties (shape and size distribution), i.e., larger p i for larger elongation (e.g., Lee & Draine 1985;Draine & Hensley 2021).Equation ( 11) reveals that the dust polarization degree depends on three key parameters, including intrinsic polarization (p i , i.e., the grain shape), grain alignment degree (⟨f align ⟩), and the projection of the mean B-field onto the POS (sin 2 γ).

Grain Alignment Theory
The key parameter required for the polarization model (Eq.11) is the grain alignment function, f align (a).Here, we derive f align using the modern theory of grain alignment based on radiative torques (RATs) and magnetic relaxation.
Traditionally, the alignment efficiency for an ensemble of grains with size a is described by the Rayleigh reduction factor, R (Greenberg 1968).The Rayleigh reduction factor describes the average alignment degree of the axis of major inertia of grains (â 1 ) with its angular momentum (J , i.e., internal alignment) and of the angular momentum with the ambient B-field (i.e., external alignment), given by where β is the angle between â1 and J and ξ is the angle between J and B, and ⟨...⟩ describes the averaging over an ensemble of grains (see, e.g., Hoang & Lazarian 2014, 2016a,b).Here, R = 0 for randomly oriented grains, and R = 1 for perfect internal and external alignment.Therefore, we can define the grain alignment function as follows (Hoang & Lazarian 2014): • Grain alignment by RATs only In general, grain alignment by RATs can occur at low-J attractors and high-J attractors (Weingartner & Draine 2003;Hoang & Lazarian 2008, 2016a).Let f high−J be the fraction of grains that can be aligned by RATs at high-J attractors.For grain alignment by only RATs (e.g., grains of ordinary paramagnetic material, Davis & Greenstein 1951), high-J attractors are only present for a limited range of the radiation direction that depends on the grain shape (Draine & Weingartner 1997;Hoang & Lazarian 2008;Herranen et al. 2019).
Numerical simulations in Hoang & Lazarian (2008, 2016a) show that if the RAT alignment has a high-J attractor point, then, large grains can be perfectly aligned because grains at low-J attractors would be randomized by gas collisions and eventually transported to more stable high-J attractors by RATs.On the other hand, grain shapes with low-J attractors would have negligible alignment due to gas randomization.For small grains, numerical simulations show that the alignment degree is rather small even in the presence of iron inclusions because grains rotate subthermally (Hoang & Lazarian 2016a).
The first parameter of the RAT alignment theory is the minimum grain size required for grain alignment, denoted by a align .Let n H and T gas be the local gas density and temperature.Let γ rad , u rad , and λ be the anisotropy degree, radiation energy density, and the mean wavelength of the radiation field.The minimum alignment size by RATs is given by Hoang et al. (2021): where ρ = ρ/(3 g cm −3 ) is the normalized mass density of grain material, T gas,1 = T gas /10 K, n 3 = n H /(10 3 cm −3 ), γ −1 = γ rad /0.1, U = u rad /u ISRF with u ISRF = 8.6 × 10 −13 erg cm −3 being the radiation energy density of the interstellar radiation field (ISRF) in the solar neighborhood (Mathis et al. 1983), and F IR is a dimensionless parameter that describes the grain rotational damping by infrared emission.For the ISM and dense clouds, F IR ≪ 1, and can be omitted in Equation ( 14).Therefore, the Rayleigh reduction factor can be parameterized by the fraction of grains aligned at low-J and high-J attractors (Hoang & Lazarian 2014): where Q X = ⟨(3 cos 2 β − 1)/2⟩ describes the internal alignment of grain axes with J (Roberge & Lazarian 1999), and the external alignment of J with B is assumed to be perfect (Hoang & Lazarian 2008, 2016a).At high-J attractors, Q high−J X can reach 100% due to suprathermal rotation, but at low-J attractors with thermal rotation, the value of Q low−J X is smaller which follows the local thermal equilibrium (Boltzmann) distribution of the angle β with roughly 30% (Hoang & Lazarian 2016a,b).2Therefore, the second term is subdominant, and the alignment efficiency R only depends on f high−J .
• Grain alignment by MRAT Composite grains considered in this paper contain silicate, carbonaceous material, and clusters of metallic iron or iron oxide.Let N cl be the number of iron atoms per cluster, and ϕ sp be the volume filling factor.The existence of metallic iron/iron oxide clusters makes composite grains become superparamagnetic material that has magnetic susceptibility increased by a factor of N cl from ordinary paramagnetic material (Davis & Greenstein 1951;Jones & Spitzer 1967).The enhanced magnetic susceptibility significantly increases the rate of magnetic relaxation (denoted by τ −1 mag,sp ) and the degree of grain alignment (Lazarian & Hoang 2008;Hoang & Lazarian 2016a;Lazarian & Hoang 2019).
The strength of magnetic relaxation for rotating composite grains is defined by the ratio of the magnetic relaxation rate relative to the gas randomization rate (denoted by τ −1 gas ) where T d,1 = T d /10 K is the normalized dust temperature, p = p/5.5where pµ B is the mean magnetic moment per iron atom and µ B = eℏ/2m e c is the Bohr magneton, B 2 = B/10 2 µG is the normalized B-field strength, and ϕ sp,−2 = ϕ sp /10 −2 .Above, k sp (Ω) is the function of the grain angular velocity Ω which describes the suppression of the magnetic susceptibility at high angular velocity (see, e.g., Hoang 2022).Detailed numerical calculations in Hoang & Lazarian (2016a) show the increase of f high−J and R with δ mag .To model the increase of f high−J with δ m , as in Giang et al. (2023), we introduced the following parametric model: Above, for alignment grain driven by RATs only, i.e., for δ mag < 1, we assumed f high−J = 0.25.This choice is based on a statistical study of grain alignment by RATs for an ensemble of grain shapes in Herranen et al. (2019).
For a given grain size distribution and magnetic properties, using a align , δ mag and R, one can calculate ⟨f align ⟩ by plugging f align (a) into Equation ( 9).It can be seen that ⟨f align ⟩ cannot reach 1 because there exists a fraction of small dust of size below a align which are weakly not aligned.Therefore, even in the ideal RAT theory, the maximum value of ⟨f align ⟩ < 1.

Constraining the B-field inclination angle
Equation ( 11) shows the dependence of the dust polarization degree on three parameters, including the intrinsic polarization degree, p i , the mass-averaged alignment efficiency, ⟨f align ⟩, and the inclination angle γ.For a given grain shape, the intrinsic polarization is known (Draine & Hensley 2021).The MRAT alignment theory yields ⟨f align ⟩.Therefore, if the thermal dust polarization along a LOS is measured to a degree p obs , we can infer the mean inclination angle of the B-field along this LOS by setting the polarization degree predicted from Equation ( 11) to p obs , as follows: Equation ( 18) provides the direct link between the inclination angle as a function of the observed polarization degree, intrinsic polarization, and grain alignment.Therefore, when the latter parameters are constrained, one can infer the inclination angle.

The idealized case of uniform grain alignment
Chen et al. ( 2019) considered the idealized case of uniform grain alignment and assumed the constant polarization, p 0 = p i ⟨f align ⟩.In this special case, the polarization degree (Eq.11) and the inclination angle (Eq.18) return to the results in Chen et al. (2019).The latter is given by where the subscript Chen denotes the results from Chen et al. ( 2019).Chen et al. ( 2019) obtained the polarization p 0 using the maximum polarization observed toward the cloud, p max , which requires the B-field lying in the POS (i.e., sin 2 γ = 1) and perfect grain alignment of f align = 1.Using these criteria for Equation (11) one obtains which yields, The determination of the polarization parameter, p 0 , from Equation ( 21), in practicality, is uncertain because the B-fields within the observed region may not lie in the POS or the grain alignment is not perfect, which does not satisfy the criteria for the maximum polarization.In our proposed method, the intrinsic polarization, p i , is based on the input grain elongation, which can be constrained with polarization observations (see e.g., Draine & Hensley 2021).
Comparing Equation ( 18) with ( 19), one can see that the difference between the inclination angle for the realistic case with alignment is different from that of perfect alignment by where we have assumed 1+p obs ≈ 1 and 1+2p 0 ⟨f align ⟩ ≈ 1.
Therefore, the alignment degree affects significantly the inclination angle inferred from the perfect alignment (PA) model.Table 1 shows the importance of grain alignment.The variation of the inclination angles in the presence of grain alignment, assuming a fixed maximum polarization and observed polarization degree of p obs = p max /2.The perfect alignment yields the inclination angle of 41.25 • , but the imperfect alignment yields a wide range of the inclination angle.If the alignment efficiency is decreased to 60%, the inclination angle is decreased by a factor of 2. The inclination angle cannot be retrieved if the alignment degree is reduced to less than 50%.
Equation ( 18) has the physical solution of γ only when sin 2 γ ≲ 1.However, in reality, one may encounter singular situations of sin 2 γ > 1 because the alignment efficiency f align ≪ 1, as shown in Table 1.

SYNTHETIC POLARIZATION OBSERVATIONS OF MHD SIMULATIONS
The physical model of thermal dust polarization described by Equation ( 11) is obtained assuming the homogeneous gas, dust, and B-field properties along the LOS.In reality, these properties change along the LOS due to the turbulent nature of the ISM and grain evolution.Therefore, to infer the inclination angle using the observed polarization degree, we first need to test whether the physical polarization model can adequately reproduce the observed polarization.To do so, we will post-process MHD simulations with our POLARIS code to obtain the synthetic dust polarization and compare it with the results obtained from the model.We use a snapshot taken at time t = 0.783 kyr of MHD simulations that follows the evolution of a filamentary molecular cloud under the effects of self-gravitating and B-fields by Ntormousi & Hennebelle (2019).The filamentary cloud was set up with a length of 66 pc and an effective thickness of 33 pc onto a three-dimensional data cube of 256 3 pixels.The cloud is supersonic with M s ∼ 5−10, and the gas motion is mostly driven by selfgravity and turbulent fragmentation.The turbulence is initially added with the initial ratio of the free-fall time and the turbulence crossing time q turb = t ff /t turb set to unity, and the energy spectrum followed by Kolmogorov's power law with the slope of -5/3.The initial mean B-field has a strength of B = 5µG and is perpendicular to the filament axis.The simulation box size of 66 pc was chosen to cover the whole morphology of the cloud.

MHD simulation datacube
Figure 2 shows the maps of the gas column density with the segments showing the B-field orientation (left panel) and the Alfvén Mach number averaged along the LOS (right panel) from the MHD datacube.The gas column density changes significantly within the filament, from N H ∼ 10 21 − 10 22 cm −2 in the outer region to N H ∼ 4 − 5 × 10 23 cm −2 in the central region.The filament is mostly sub-Alfvénic with ⟨M A ⟩ ∼ 0.2 − 0.8.

Dust and Grain Alignment Models
Here, we consider the composite model of interstellar dust with 67.5% of silicate and 32.5% of carbonaceous grains.The choice of composite dust model is advantageous because the polarization degree at long wavelengths in Equation ( 11) is independent of dust temperature.It is also motivated by the latest studies that advocate for a single dust component of the ISM (Draine & Hensley 2021;Hensley & Draine 2023).The grain size distribution is described by the power law with dn ∝ Ca −3.5 da in the range of grain sizes from a min = 5 nm to a max = 0.25 µm (Mathis et al. 1977).The grain shape is assumed to be oblate spheroid with the axial ratio of s = 2.The choice of such an elongation is for the sake of convenience because of the existing pre-calculated cross-sections in the POLARIS code.However, it only affects the intrinsic polarization p i and does not affect the resulting conclusions (see Appendix C for the case of smaller elongation with s = 1.4).
To quantify the effect of grain alignment on synthetic dust polarization and inferred inclination angle, we consider the different grain alignment models.Table 2 shows the list of considered grain alignment models, including perfect alignment (PA), ideal RAT, and MRAT with different magnetic properties.For the radiation field, we consider grains to be irradiated by only the ISRF, which is described by the average radiation field in the solar neighborhood (Mathis et al. 1977).For grain magnetic properties, we consider both paramagnetic grains (PM) with the fraction of iron of f p = 0.1, and superparamagnetic (SPM) grains with increasing levels of iron inclusion N cl from 50 to 10 3 .Following the MRAT mechanism (Hoang & Lazarian 2016a), SPM grains are predicted to achieve a higher alignment degree than PM grains due to higher magnetic susceptibil- SPM, N cl = 50, 10 2 , 10 3 ity, with a higher fraction of grains being aligned with B-fields (see Equation 13 -17).

Synthetic Polarization Modeling with POLARIS
The RAT alignment physics was implemented in the POLARIS code by Reissl et al. (2016), whereas the MRAT mechanism was recently incorporated in the updated version by Giang et al. (2023).Here, we use the latest version of POLARIS code by Giang et al. (2023) for our numerical studies.For a detailed description of the implementation of RAT alignment physics and polarized radiative transfer, please see Reissl et al. (2016).For a detailed description of the alignment model parameters, see Giang et al. (2023).

NUMERICAL RESULTS
Here we first show the numerical results from synthetic polarization observations and compare them with the polarization calculated by the analytical model.We consider the different viewing angles, denoted by γ view , which is defined be the angle between the mean B-field B 0 and the LOS (see Figure 1).We will then apply our technique to infer the inclination angle of the B-field using the synthetic data and obtain the full 3D B-field strength.

Grain alignment size and mass-weighted grain alignment efficiency
For numerical calculations of the mass-weighted alignment degree with POLARIS, denoted by ⟨f align ⟩, we calculate f align for each cell using the alignment size a align and δ mag and integrate it along the LOS.
Figure 3 shows the map of the mass-weighted grain alignment efficiency for the different alignment models, including perfect alignment, ideal RAT alignment, and MRAT alignment (see Table 2 for details), assuming the mean B-field in the POS (γ view = 90 • ).For the model of perfect alignment, all grain sizes are perfectly aligned with B-fields with a constant ⟨f align ⟩ = 1 in the entire filament.In the case of ideal RAT, by contrast, the alignment degree is strongly affected by gas randomization, particularly in the inner regions with high gas density (see, e.g., (Hoang et al. 2021)).Then, the value of ⟨f align ⟩ decreases from 0.9 to 0.1 along the filamentary cloud.
For the models of MRAT alignment, the grain alignment efficiency significantly depends on the magnetic properties of grains.When grains are paramagnetic, only 25% of their population can be at high-J attractors by RATs, while the remaining grains at low-J attractors are aligned with a lower alignment degree of 30% (see Section 2.1.3).This produces a lower alignment efficiency with ⟨f align ⟩ < 0.5.For SPM grains, the increasing levels of embedded iron enhance the magnetic relaxation and the fraction of grains at high-J attractors (see Equation 17).Consequently, these grains can be efficiently aligned by MRAT, with a maximum alignment efficiency up to ⟨f align ⟩ = 0.9, the same as the maximum value achieved for the Ideal RAT model.
To compare with synthetic results, we also calculate the mass-weighted alignment degree directly through Equation ( 9), denoted by ⟨f ana align ⟩, using the values of a align and δ mag from Equations ( 14) and ( 16), respectively, with the input parameters of the mean gas density ⟨n H ⟩ and B from MHD simulations.Figure 4 shows the comparison of ⟨f align ⟩ obtained from synthetic observations and ⟨f ana align ⟩ obtained from our analytical calculations for the different models of grain alignment.Regardless of grain magnetic properties, the results show a good agreement between synthetic alignment from PO-LARIS and analytical calculations.However, as seen in Figure 4, ⟨f align ⟩ is slightly higher than ⟨f ana align ⟩.This is because, in the synthetic modeling by POLARIS, we take the variation of the anisotropic degree of interstellar radiation fields from 0.8 in the outer regions to 0.2 in the inner regions of the cloud rather than the standard ISM anisotropic degree of 0.1 in the entire cloud as we adopt in analytical calculations.Therefore, the analytical formula can be used to calculate the mass-weighted alignment efficiency of grains.

Synthetic Polarization Versus Analytical Polarization Model
Figure 5 shows the maps of synthetic dust polarization at the wavelength of 870 µm obtained from the PO-LARIS code for different grain alignment models, assuming the mean B-field in the POS (γ view = 90 • ).
For the perfect alignment model, the polarization degree is highest and can reach a maximum value of P syn,max ≈ 57.8% at the outer region.However, it decreases considerably in the central region (blue area) due to the depolarization effect caused by magnetic turbulence.For the model of Ideal RAT, the polarization degree is lower than the PA model due to the lack of alignment of small grains, which becomes more important toward dense regions.As a result, the polarization degree decreases from the maximum value of 50.8% in the outer and diffuse region to ∼ 5% at the central region of highest density.
For the model of PM grains, the overall polarization degree reduces with p syn,max ≈ 26.9% due to the lower alignment efficiency.The polarization degree is significantly enhanced for grains with high iron levels and can recover to the Ideal RAT case when N cl reaches 10 3 with P syn,max ≈ 50.2%.In particular, compared to the PA model, the alignment models based on RATs predict a more prominent polarization hole in the center of the filament (blue areas), which is caused by the loss of RAT alignment (i.e., increase in a align ) due to the increase of gas randomization and the strong attenuation of the ISRF.
Figure 6 shows the resulting maps of polarization degree, but for different viewing angles, γ view = 0 • (i.e., the mean B-field parallel to the LOS), 30 • , 60 • and γ view = 90 • (i.e., the mean B-field perpendicular to the LOS).One can see that the polarization degree is significantly reduced from the optimal values at γ view = 90 • to the smallest values at γ view = 0 • due to the projection effect.Note that the polarization degree is still considerable of P syn,max < 16% when the mean B-field is along the LOS (γ view = 0 • ) because the local B-fields are not along the LOS due to magnetic turbulence.
To test whether our physical polarization model can adequately describe the synthetic polarization, we calculate the polarization degree, P ana , from Equation (11), by using p i (determined by the maximum polarization for the PA model with γ view = 90 • , i.e., optimal con-Figure 3. Map of mass-weighted grain alignment efficiency from POLARIS for the different alignment models.The black contours show the grain alignment efficiency of 0.5, 0.7, and 0.85, respectively.The alignment efficiency is uniform with ⟨f align ⟩ = 1 for the perfect alignment case.However, the grain alignment is strongly affected by the gas randomization, resulting in the variation of alignment efficiency within the cloud with ⟨f align ⟩ < 1 in the Ideal RAT case.The alignment efficiency is much lower to ⟨f align ⟩ < 0.5 when grains are paramagnetic, and can be enhanced for SPM grains having high levels of iron inclusions.ditions for polarization), the mean alignment efficiency ⟨f align ⟩ obtained from our numerical simulations and the mean inclination angle of B-fields ⟨γ B ⟩ from MHD simulations.
The left panel of Figure 7 compares the synthetic polarization degree obtained from POLARIS with the polarization degree calculated by the analytical formula given by Equation ( 11) for various models of grain alignment and γ view = 90 • .The right panel shows the same results but for different viewing angles of γ view = 0 • , 30 • , 60 • and 90 • and the model of Ideal RAT alignment.In all cases, a good correlation between P syn and P ana is observed, implying that the analytical formula could be used to describe the synthetic polarization.However, the results obtained from synthetic observations tend to be lower than those from analytical calculations.This is caused by the fact that the synthetic polarization can account for the effect of B-field fluctuations relative to the mean B-field along the LOS, whereas the analytical model assumes a uniform B-field along the LOS.This issue will be studied in Section 5.As shown in Figure 7, the analytical dust polarization model (Eq.11) can describe the synthetic polarization.We now apply our method in Section 2 to infer the inclination angles and compare them with the actual angles from MHD simulations.

Inferring inclination B-field angles from synthetic polarization
Figure 8 shows the map of inclination angles obtained from our synthetic observations (γ align syn ) for the different models of grain alignment, assuming the viewing angle γ view = 90 • .Generally, the inclination angle is about 60-80 degrees and tends to decrease to 40-50 degrees toward the dense regions due to the effects of magnetic fluctuations on dust polarization degree (see Figure 7).When applying our method with the variation of grain alignment, the inclination angle does not change considerably.However, it will return NAN values in the densest regions, which is a result of the reduced grain alignment efficiency by enhanced gas randomization and attenuation of ISRF (as predicted in Table 1).The situation is worse for PM grains with small ⟨f align ⟩ < 0.5, and the value of γ align syn cannot be extracted in some regions of the cloud.This can be resolved when grains have iron inclusions with due to their efficient alignment by MRAT, and the inclination angles can be retrieved from the outer to the inner regions of the cloud.
Using the synthetic polarization data, we also apply the Chen technique and infer the mean inclination angle for each grain alignment model, denoted by γ Chen syn , which is shown in Figure 9.In comparison with our method, the inclination angles from Chen's method vary significantly with the magnetic properties of grains, with the increased inclination angle for SPM grains with high levels of embedded iron.
To quantify the quality of our method, we compare the inclination angles of B-fields extracted from synthetic dust polarization using both our method and Chen's method with the true inclination angles from MHD simulations.For each cell of the simulation box, the true inclination angles of B-fields are defined as and the mean inclination angle is calculated by The inferred inclination angles can also be compared with the density-weighted of the mean inclination angle along the LOS, denoted by ⟨γ ρ ⟩, which is determined by where γ B (i,j,k) are the local inclination angles derived from Equation ( 23).
Figure 10 shows the distribution of the inclination angles obtained from our method (solid color lines) and Chen's method (dashed color lines), compared with the true angle ⟨γ B ⟩ (solid black line) and the densityweighted angle γ ρ (dash-dotted black line).The most probable values of the inclination angles (i.e., the inclination angle at the peak distribution) from each technique are summarized in Table 3.As one can see, the peak inclination angle from our method tends to be close to the case of perfect alignment (yellow line) and only different by 5 − 10 • , whereas that from Chen's method is lower and differs by 10 − 20 • from the true value.In particular, the mean inclination angle from our method is insensitive to the grain alignment model and grain magnetic properties.This is expected because we consider the alignment efficiency in deriving the inclination angle.
Figure 11 shows the map of inclination angles for the different viewing angles.One can see that, despite the variation in viewing angles, the inclination angles can still be extracted accurately.However, the cloud-scale inclination angles contribute significantly to the depolarization when γ view is small (i.e., the B-fields are parallel to the LOS), resulting in fewer NAN values when applying our method in deriving the inclination angles (see also in Table 1).The color solid lines show the results from our method, and the color dashed lines for Chen's method.The peak values from our method are less different than the true angle by 5-10 degrees and do not change regardless of grain magnetic properties.
Table 3. Summary of the most probable values (i.e., the peak inclination in the histograms) of the inclination angles derived from the true B-fields in MHD simulation, the mean inclination angles from the Chen technique and our method.Figure 12 shows the histogram of the offsets of the inclination angles obtained from our method and Chen's method for the different viewing angles (same as the right panel of Figure 10).Table 4 summarizes all probable values of the inferred inclination angles for the different viewing angles.It is clear that in spite of the change in viewing angles, our method still provides more accurate inclination angles than Chen's method for various models of grain alignment, with a difference of less than ∼ 5 − 10 degrees.The inferred inclination angles from Chen's method (dashed lines), on the other hand, strongly depend on grain alignment models and viewing angles, with the difference from the true inclination angle of 10 − 20 degrees.The Chen's method becomes less accurate for larger γ view as a result of the increasing depolarization effect by grain alignment loss.

Full B-field strength from synthetic polarization
With our inferred inclination angle, we can calculate the strength of the 3D B-field as follows where B POS = B 2 x + B 2 y is the B-field strength projected in the plane-of-sky, which is taken from the MHD simulation.
Figure 13 compares the full 3D B-field strength from Equation ( 26) with the actual strength ⟨B⟩ = 7.5µG from the MHD simulation.One can see that, by taking the inclination angles from our technique, the calculated 3D field strength differs from the true B-field strength by ∼ 0.5 − 1 µG, which is only about 10 − 20% of the actual value, and the inferred value is independent of the grain alignment models.Meanwhile, the lower inclination angles derived from Chen's technique (Figure 10) lead to higher values of 3D B-field strength, with the difference of ∼ 2 − 3 µG, i.e., about 40 − 60% from the true value (see right panel).

Physical dust polarization model with magnetic turbulence
In Section 2, we focus on the effects of grain alignment on the polarization degree by assuming the uniform inclination angle of the mean B-field along the LOS.In reality, the magnetic turbulence induces fluctuations in the local B-field and causes an additional depolarization effect.In the presence of magnetic turbulence, the polarization cross-section in the observer's xŷ coordinate system (see Figure 1) can be written as where ∆θ is the deviation angle between the local B-field (B) and the mean field (B 0 ), and describes the effect of magnetic turbulence on the polarization degree (Lee & Draine 1985).
The degree of thermal dust polarization (Eq.11) now becomes Setting p(λ) to p obs , we obtain the inclination angle For small fluctuations of ∆θ < 1 (i.e., sub-Alfvénic turbulence), one can make the approximation where δθ = ⟨(∆θ) 2 ⟩ 1/2 denotes the dispersion of the deviation angle.Equation ( 32) reveals that the increase in the B-field angle dispersion δθ (stronger magnetic turbulence) reduces the value of F , which decreases the polarization degree.Therefore, if the turbulence is isotropic, the Bfield angle dispersion δθ corresponds to δϕ, but their correlation depends on the viewing angles, as shown in Figure 26 (see Appendix B for details).

Comparison to synthetic dust polarization
For our synthetic observations of MHD simulations, the deviation angle between the local and the mean Bfield is calculated in each cell of the simulation box as where B, B 0 are the local and mean B-fields from MHD datacube.From that, we integrate the cos 2 (∆θ) over the LOS and calculate the dispersion angle δθ, the factor F turb and the analytical polarization degree as shown in Equation ( 30).
Figure 14 shows the resulting maps of the dispersion angle δθ in the entire filamentary cloud, for the different viewing angles.For a given viewing angle, the magnetic fluctuations are weaker in the outer region and tend to increase toward denser regions with δθ > 50 • .Moreover, the fluctuations tend to increase when the viewing angle decreases.
Figure 15 shows the maps of F turb calculated using the angle dispersion maps in Figure 14.The value of F turb decreases from the outer to the inner region due to stronger turbulence.Moreover, the depolarization is stronger when the mean field is closer to the LOS, also due to stronger fluctuations.
Figure 16 (left panel) shows the comparison of the synthetic polarization to the analytical model similar to the previous Figure 7 but in the presence of random Bfields described in Equation ( 30).In all cases, when the additional term F turb is included, the analytical model becomes in agreement with the synthetic polarization than the model without the turbulence term shown in Figure 7.

Inferred inclination angles and 3D B-field strengths
Following this step, we include the factor F turb for improving the calculation of inclination angles through dust polarization degree, as described in Equation ( 31).We take the small test demonstrating the dependence of sin 2 γ on the grain alignment efficiency ⟨f align ⟩ and the turbulence factor F turb as shown in Figure 17, assuming p obs = p max /2.The values of sin 2 γ will be larger than 1 when ⟨f align ⟩ < 0.5 (i.e., significant loss of grain alignment) or F turb < 0.5 (i.e., high magnetic turbulence), and can return NAN values when extracting the B-field inclination angles.
Figure 18 shows the results of the inferred inclination angles when F turb is considered, denoted by γ align,F syn .Compared to Figure 11, the maps show more NAN      31) with respect to the grain alignment efficiency ⟨f align ⟩ and the turbulence factor F turb = 1, 0.8, 0.5, considering p obs = pmax/2 = 0.15.The inferred inclination angles cannot be extracted when ⟨f align ⟩ < 0.5 or F turb < 0.5.
data points where the thermal dust polarization is dominantly affected by both the magnetic turbulence and the grain alignment loss.As a result, the remaining data points demonstrate the regions in which the magnetic turbulence is weaker, and thus depolarization is mainly by the mean B-field inclination angles and, thus, potentially provide better values of inclination angle.Table 5 summarizes the most probable values of inferred inclination angle, considering the effects of random fields.In comparison to Table 4, the inclination angles inferred from our method have improvements and are much closer to the actual values from MHD simulations (γ ∧ ρ ).The angles inferred using Chen's method do not change because of their assumption of uniform B-fields.
Figure 19 shows the difference in the inclination angles inferred from our polarization model when accounting for the effect of turbulence with respect to the actual B-field inclination angles.Compared to Figure 12, the inferred inclination angle is much better in agreement with the true inclination angle by only 2-6 degrees in all cases of grain alignment models and various viewing angles.
From the improvement of the inferred inclination angles, we calculate again the strength of 3D B-fields, as previously shown in Equation ( 26). Figure 20 illustrates the histogram of the 3D B-field strength and its offset when the effect of magnetic turbulence is included.Compared to previous results in Figure 13, the derived strength from our method is significantly improved and shows more precisely the mean B-field strength of ⟨B⟩ ∼ 5 − 7.5 µG.The B-field strength inferred from Chen's method remains unchanged.

A physical model of thermal dust polarization based on the MRAT alignment theory
An accurate, physics-based model of dust polarization is necessary for interpreting polarimetric data and constraining dust physics and B-fields.In this paper, we introduced a physical model of thermal dust polarization, which is the function of four key parameters, including the intrinsic polarization degree (p i ), mass-weighted alignment efficiency (⟨f align ⟩), the mean B-field inclination angle (sin 2 γ), and a term describing the depolarization effect caused by magnetic turbulence, F turb (δθ) (see Eq. 11 for the model without magnetic turbulence and 30 for the model with magnetic turbulence).The intrinsic polarization p i is determined by the grain shape elongation.In contrast to previous models that constrain grain alignment using polarization data (Draine & Hensley 2021;Hensley & Draine 2023) or treat it as an input parameter, here, the value of ⟨f align ⟩ is calculated from the MRAT alignment theory, which depends on the local gas density, radiation field, grain magnetic susceptibility, and B-field strength (Hoang & Lazarian 2016a;Lazarian & Hoang 2021).The B-field inclination angle describes the projection effect of the mean B-field onto the POS.Finally, the term F turb depends on the nature of magnetic turbulence.
We have tested our physical model of dust polarization degree, both without and with magnetic turbulence, using synthetic polarization observations of MHD simulations for a filamentary cloud with our upgraded POLARIS code (see Section 3).We found that the physical polarization model without magnetic turbulence can reproduce the synthetic polarization (see Figure 7).However, the physical polarization model with magnetic turbulence reproduces better the synthetic polarization than the model without turbulence (see Figure 16).Our tested physical polarization model has broad implications for constraining B-fields, dust properties, and grain alignment.

Constraining 3D B-fields with dust polarization:
comparison of our method to previous studies 3D B-fields are crucial for understanding the dynamical role of B-fields in the ISM evolution and star formation (Hull & Zhang 2019;Pattle et al. 2022), especially in the formation and evolution of interstellar filaments (Hennebelle & Inutsuka 2019).However, to date, accurate measurements of 3D B-fields are not yet available.Chen et al. (2019) proposed a new method to infer    the inclination angle of the B-field using the analytical model of the dust polarization degree.However, the authors adopted an idealized model of dust polarization by assuming that both the B-field and grain alignment efficiency do not change along the LOS.HL23 improved Chen's method by considering the fluctuations of the Bfield along the LOS.However, both in Chen and LH23, the grain alignment is assumed to be uniform within MCs, which is unrealistic due to the varying local gas density and radiation fields in MCs.
In this paper, we introduced a general method for constraining the 3D B-field via dust polarization data and using the physical model of dust polarization by incorporating the variation of the grain alignment across the cloud using the modern MRAT alignment theory (Hoang & Lazarian 2016a;Lazarian & Hoang 2021).We applied our technique to synthetic polarization data and found that the inferred inclination angles are much closer to the true values of MHD data than Chen's method (see Fig. 12 and 19).Our new technique is particularly useful for tracing 3D B-fields in star-forming regions when the gas density and local radiation field change rapidly.In our follow-up studies, we will apply this method to real observational data of dust polarization to starforming regions by SOFIA/HAWC+ and JCMT/POL2 and present the results elsewhere.
We also note that both Chen et al. ( 2019) and Hu & Lazarian (2023a) derived the polarization p 0 using the maximum polarization observed from the cloud.This is not well constrained because the value of p 0 is only determined by the optical conditions for the maximum polarization, p max , including the perpendicular B-field and perfect alignment, which is not necessarily satisfied within a cloud.Moreover, the location with observed p max may have an inclined B-field and imperfect alignment.Therefore, there exists the uncertainty in p max and then p 0 .In our present study, the value of p i is directly determined by the grain shape/elongation, which is within a factor 2 different for the grain axial ratio from 1.4 − 2 (Draine & Hensley 2021).
In addition, the combination of our method using dust polarization with other B-field measurement techniques can further support probing both B-field and dust properties in star-forming regions.For instance, the Faraday rotation method (see, e.g., Tahani et al. 2018;Tahani 2022) can trace the orientation of B-fields pointing toward or away from us.By combining with our method using dust polarization, we can reconstruct both the morphology and directions of 3D B-fields, particularly in molecular clouds.Also, we note that the inclination angles of the mean B-fields can be potentially inferred by measuring both POS B-fields using dust polarization and LOS B-fields using Zeeman splitting in spectral lines (Crutcher 2010) with tan γ = B POS /B LOS .By comparing with the inclination angles derived from our technique using dust polarization degree, we can reversely retrieve the information on dust properties such as grain sizes, shapes, and magnetic properties of PM and SPM grains with iron inclusions.This is a new approach to constrain dust properties in star-forming regions, which will be presented in future studies.6.3.Can p × S describe the average grain alignment efficiency?
The product of the polarization degree (p) and the polarization angle dispersion (S), i.e., p × S, is usually used to probe the average efficiency of grain alignment using observational polarization data (Planck Collaboration et al. 2015b). Following Planck Collaboration et al. (2020), we calculate the polarization angle dispersion using the second-order structure function of the polarization angles: where the summation is taken over N (δ) pixels within the annulus centered at r having the inner radius δ/2 and outer radius 3δ/2.
In principle, S describes the level of B-field fluctuations, and larger/smaller S corresponds to stronger/weaker turbulence, which decreases/increases the observed polarization degree.Therefore, the product of p × S can characterize the average alignment efficiency (Planck Collaboration et al. 2015b;Planck Collaboration et al. 2020).To test this, we use the polarization degree and mean alignment efficiency from our synthetic observations.We calculate S by choosing the square box of δ = 3 × 3 pixels.Figure 21 shows p × S vs. ⟨f align ⟩ for different alignment models (left panel) and various viewing angles (right panel).It shows that p×S has an overall correlation with the mean alignment efficiency, ⟨f align ⟩.However, the correlation has a slope of α ∼ 0.4 for ⟨f align ⟩ < 0.7 in denser regions of the cloud.In lower-density regions with ⟨f align ⟩ > 0.7, the decreased level of magnetic turbulence is more prominent (Figure 14) than the increased grain alignment efficiency (Figure 3), resulting in a downward fluctuation in the product p × S.This value decreases significantly for lower viewing angles γ view due to the effect of inclination angles of the mean B-fields.This reveals that the product p × S may qualitatively describe the alignment efficiency, but it is not recommended to take this product as an accurate indication of grain alignment efficiency.For probing grain alignment with observations, it is suggested to use the mass-weighted alignment efficiency.

Dependence of the intrinsic polarization on grain shapes
The intrinsic polarization p i = σ pol /σ ext at farinfrared/sub-mm wavelengths much larger than the grain size is determined by the grain elongation (Lee & Draine 1985).For our numerical study, we assumed the grain elongation of s = 2 that has the intrinsic polarization of p i ≈ 0.58.Draine & Hensley (2021) computed the values of p i at λ = 850 µm for the Astrodust model and found that p i decreases by a factor of ∼ 2 from p i ≈ 0.58 at s = 2 to p i ≈ 0.26 for s = 1.4 (see Appendix C).Based on their starlight polarization, their constraint for the grain elongation is s ≳ 1.35.Moreover, from the RAT theory, we show that the maximum value of ⟨f align ⟩ ≈ 0.9 and it cannot reach 1 because the population of small grains of size a < a align ∼ 0.03 − 0.05 µm are weakly aligned.
We also note that our method of constraining 3D Bfields using dust polarization degree is based on the MRAT alignment theory for composite grains containing embedded iron clusters.The composite grain model is the most likely model of dust in molecular clouds and dense regions of star-forming regions due to the dust evolution via grain-grain collisions.The composite model is also the leading candidate dust in the diffuse ISM (Draine & Hensley 2021).Therefore, our method can be applied to various astrophysical environments, from diffuse to very dense regions.
Lastly, from observations of dust polarization, one can constrain the intrinsic polarization p i based on the maximum observed polarization degree and the B-field inclination angle and grain alignment efficiency.Nevertheless, due to grain evolution, the grain shape and elongation may change with local environments.However, we expect that the change in grain shapes and then p i is less sensitive than other parameters that affect dust polarization, such as the grain size, alignment degree, and B-fields.

Effects of beam size on inferred inclination angles from dust polarization degree
Above, we have ignored the effect of beam size on the synthetic dust polarization on inferred B-field inclination angles.To examine the effect of beam sizes, we apply a smooth Gaussian beam averaging and convolve with the synthetic polarization data, for three beam sizes of 5, 10, and 15 arcmin.This corresponds to the physical resolutions of 0.15, 0.3, and 0.45 pc, for the case of Orion A located at ∼ 420 pc (Großschedl et al. 2018).We include the beam convolution effect in the calculations of mean alignment efficiency ⟨f align ⟩, the depolarization factor by magnetic turbulence F turb and the synthetic polarization degree P syn , and recalculate the inferred inclination angles γ syn from Equation 31.
Figure 22 shows the beam convolution effect on the mean alignment efficiency ⟨f align ⟩ (top), the factor F turb (middle) and the synthetic polarization degree P syn (bottom), assuming the Ideal RAT alignment model and γ view = 90 • .Due to the beam convolution, all features smaller than the beam size are partially averaged out (see, e.g., King et al. 2018), which results in the decrease in the mean alignment efficiency ⟨f align ⟩, but the increase in F turb .However, the synthetic polarization degree, P syn , decreases with increasing the beam size, which implies that the depolarization caused by the Bfield tangling within the beam is stronger for the larger beam.
Figure 23 compares the synthetic polarization with the analytical polarization, same as in Figure 16, but when the beam convolution effect is taken into account.The Ideal RAT alignment and γ view = 90 • are considered.Due to the depolarization caused by the B-field tangling within the beam, the synthetic polarization degree is essentially lower than the analytical values of P ana obtained by Equation (30) that did not consider the effect of B-field tangling within the beam (see also Figure 7).
Figure 24 shows the maps of inferred inclination angles γ syn when the beam convolution is considered.The distributions of inferred inclination angles for different beam sizes are presented in Figure 25.One can see that for the larger beam size, the inferred inclination angles become lower than the real values and less accurate due to the disregard of the B-field tangling within the beam in the analytical polarization model. 3 Note that the choices of arcminute beam are comparable to the Planck beams (Planck Collaboration et al. 2020), which is applicable in inferring 3D B-fields in molecular clouds on a large scale of > 10 pc.In our follow-up study, we will further examine our technique in tracing 3D B-fields in smaller scales of filaments and dense cores, which is potentially applied in 3 Our technique could be improved by introducing one more parameter, e.g., F beam , that describes additional depolarization caused by B-field tangling within the beam to the analytical model (Eq.30).The value of F beam is smaller for larger beams due to stronger B-field tangling within the beam, which increases the value sin 2 γ given by Equation ( 31) and brings the inferred B-field inclination angles closer to the real values.This issue will be quantified in our follow-up studies.The solid lines show the running mean of the product p × S, while the dashed lines show the power-law fit to that running mean.In general, the product p × S has a good correlation with ⟨f align ⟩ with a slope of 0.4.However, the product of p × S strongly fluctuates for higher f align > 0.7, which is a result of the reduction of magnetic turbulence toward low-density regions and decreases significantly for lower γview.
recent single-dish observations from JCMT/POL-2 and SOFIA/HAWC+ with higher spatial resolution.
7. SUMMARY 1.We introduced a physical model of thermal dust polarization using the grain alignment theory based on RATs and magnetic relaxation.At far-IR/sub(-mm) wavelengths, the degree of thermal dust polarization is a function of the intrinsic polarization, the mass-weighted alignment degree, and the B-field inclination.
2. To test our physical polarization model, we performed synthetic polarization observations of MHD simulations for a filamentary cloud using the upgraded POLARIS code with the MRAT alignment theory.We found that the physical polarization model can accurately describe the synthetic dust polarization.
3. From synthetic observations, we found that the mass-weighted alignment efficiency by the MRAT mechanism varies significantly from the outer layer of the filament toward the central region due to the increase in the gas density and decrease of the interstellar radiation field.This reproduces the polarization hole in the synthetic data.
4. Using our tested physical polarization model and synthetic polarization data, we inferred the incli-nation angle and the strength of 3D B-fields and compared them with the actual values from MHD simulations.Compared to previous methods that disregarded grain alignment, our method with the MRAT alignment provides a more accurate inference of the inclination angle and 3D B-fields.
5. Our physical polarization model with magnetic turbulence can reproduce better the synthetic polarization and enables more accurate inference of the 3D B-fields than the polarization model without magnetic turbulence.
6. Our new method can be applied to real observational data to constrain 3D B-fields using dust polarization data.
7. When the 3D B-fields are available, our physical polarization model can help constrain the dust properties, such as the grain elongation and magnetic properties, and grain alignment physics.Figure 26 shows variation of the magnetic angle dispersion δθ (left panel) and F turb (right panel) with the polarization angle dispersion σ ϕ calculated from the un-sharp masking method for different viewing angles, assuming the squared box size of 3 × 3 pixels.The power-law fit is performed to the running means of δθ and F turb .Overall, the slope tends to be steeper for larger viewing angles.For small σ ϕ , the slope is rather shallow of |η| < 0.1, which implies the insignificant effect of turbulence on the dust polarization degree.The slope become steeper for σ ϕ > 5 • with |η| > 0.1, but it varies with viewing angles.For lower viewing angles (i.e., the mean field is along the LOS), the impact of magnetic turbulence is more prominent in the LOS rather than in the POS; consequently, the F turb is less correlated with the dispersion σ ϕ in the POS with the shallow slope η 2 ∼ −0.03.The effect of turbulence is more significant in the POS with increasing viewing angles γ view > 30 • , which spans η 2 ∼ −0.1 to −0.3.This implies that the factor F turb can potentially be retrieved from the polarization angle dispersion in the POS polarization map, but the effect of the inclination angle itself provides uncertainty.
Note that the angle dispersion is a function of M A with different slopes depending on the mode of MHD turbulence (Lazarian et al. 2022).In our future papers, we will apply this method to constrain 3D B-fields using dust polarization maps.

C. SYNTHETIC POLARIZATION OBSERVATIONS FOR THE ASTRODUST MODEL AND INFERRED INCLINATION ANGLES
So far, we have emphasized the capabilities of inferring inclination angles from dust polarization using our technique and applied it to the case of the composite model of interstellar dust with the highly elongated oblate shape of s = 2. Here, we assume the Astrodust model by Draine & Hensley (2021) with a less elongated oblate shape of s = 1.4 and perform synthetic polarization of MHD simulations as in Section 3. We apply our technique to the resulting synthetic data to infer the inclination angles.
Figure 27 illustrates the synthetic results of polarization degree for the Astrodust model.Since the intrinsic polarization decreases as grains become less elongated, the overall polarization degree for all alignment models reduces significantly compared to the previous results in Figure 5 by a factor of ∼ 2. For instance, p syn,max = 50.8%decreases to p syn,max = 25.5% for the Ideal RAT case -roughly close to the maximum polarization degree found in the diffuse ISM by Planck Collaboration et al. (2020).
Figure 28 shows the maps of the inferred inclinations using our technique when the effects of grain alignment and magnetic turbulence are considered for the Astrodust model and the mean fields perpendicular to the POS (i.e., γ view = 90 • ).Despite the decreased intrinsic polarization degree due to lower grain elongation, the values of the inferred inclination angle seem to be unchanged compared to the previous calculations illustrated in Figure 18.This implies that our technique can be applied to arbitrary grain shapes.

Figure 1 .
Figure 1.Illustration of grain alignment with the mean Bfield B, the line of sight, ẑ, and the sky plane xŷ.The mean B-field lies in the plane xẑ and makes an angle γ with the line of sight.

Figure 2 .
Figure 2. Map of the gas column density (left panel) and the Alfvén Mach number (right panel) from MHD simulations of the filamentary cloud.The cloud is mostly sub-Alfvénic , with the mean geometry of the B-field along the x-axis of the plane-of-sky, illustrated by the white segments (left panel).The gas density varies considerably within the cloud, which is expected to affect the grain alignment efficiency and the inferred inclination angles through dust polarization.

Figure 4 .
Figure 4. Comparison of mass-weighted alignment efficiency from synthetic observations (⟨f align ⟩) to the results from analytical formula (⟨f ana align ⟩).Data points show the results from numerical calculations, and solid lines show the running mean values.The dashed line shows the one-to-one correlation.A good correlation between ⟨f align ⟩ and ⟨f ana align ⟩ despite different models of grain alignment.

Figure 5 .
Figure5.Map of the synthetic polarization degree from POLARIS for the different grain alignment models and the mean B-field in the POS (viewing angle γview = 90 • ).The polarization degree from the RAT and MRAT mechanisms is maximum in the outer region and significantly decreases toward the central region of high gas density due to the loss of alignment efficiency by increasing gas collision.The polarization degree is lowest for PM grains and increases with the level of iron inclusions embedded in SPM grains.

Figure 6 .
Figure 6.Same as Figure 5, but for the different viewing angles between the mean B-field and the LOS, γview.The polarization degree decreases with decreasing the inclination angle γview due to the projection effect.

Figure 7 .
Figure 7. Left panel: Polarization degree from synthetic observations (Psyn) vs. numerical calculations with the analytical formula, Pana given by Equation (11), for the different grain alignment models, assuming the viewing angle γview = 90 • .Data points show the results from numerical calculations, and solid lines show the running mean values.Right panel: For different viewing angles and Ideal RAT alignment.The dashed line shows the power law of the slope of 1 when Psyn = Pana.Both results show a good correlation, but Psyn tends to be lower than Pana.

Figure 8 .
Figure 8. Maps of the inclination angle inferred from our technique for different grain alignment models and the viewing angle γ view = 90 • .The inclination angle is around 60-80 degrees and does not change significantly for the different grain alignment models.Nevertheless, the inclination angles cannot be retrieved toward denser regions due to the significant reduction of grain alignment.

Figure 9 .
Figure9.Same as Figure9, but being retrieved from Chen's technique with the assumption of uniform grain alignment.The local inclination angles are lower than those from our method and increase with increasing levels of iron inclusions inside SPM grains.

Figure 10 .
Figure 10.Left panel: Histogram of the inclination angle inferred from our technique and Chen's method for different grain alignment models.Right panel: the offset of the inclination angle to the density-weighted of the mean inclination angle γsyn −γρ.The color solid lines show the results from our method, and the color dashed lines for Chen's method.The peak values from our method are less different than the true angle by 5-10 degrees and do not change regardless of grain magnetic properties.

Figure 11 .
Figure 11.Resulting maps of inclination angles derived from our technique similar to Figure 8, assuming the variation of viewing angles between the mean B-field and the LOS.For smaller γview, the significant depolarization is mainly caused by the local inclination angles of B-fields, leading to fewer NAN values.

Figure 12 .
Figure 12.Histogram of the difference of the inferred inclination angles from Chen's and our technique with the actual B-field inclination angles, considering different viewing angles between the mean B-field and the LOS.Even though γview is changing, our method still derives more accurate inclination angles than Chen's technique.

Figure 13 .
Figure 13.Histogram of the full B-field strength obtained from our method compared to previous methods without uniform grain alignment.The viewing angle of γview = 90 0 is assumed.

Figure 14 .
Figure 14.Maps of the fluctuations of the local to the mean B-field, δθ, calculated for the different viewing angles.Fluctuations are stronger when the mean B-field becomes closer along the LOS (smaller viewing angles).

Figure 15 .
Figure 15.Maps of the term F turb describing the effect of magnetic turbulence on the polarization degree for the different viewing angles.F turb tends to decrease toward the inner region and with the viewing angles due to stronger magnetic fluctuations.

Figure 16 .
Figure16.Same as Figure7but with the presence of B-field turbulence in the polarization model described by Equation (30).The analytical model becomes in better agreement with synthetic results when the turbulence is included.

Figure 18 .
Figure 18.Similar to Figure 11 but when the effects of random B-fields along the LOS are considered described by Equation (31).The inferred inclination angles then depend on both random B-fields and the alignment loss, which can result in more NAN values.

Figure 19 .
Figure 19.Similar to Figure12but when accounting for the effect of magnetic turbulence.The inferred inclination angles from our method agree better with the true B-field inclination angles from MHD simulations, with the angle difference less than 5 • .The results from Chen's method remain unchanged due to their assumption of uniform B-field and alignment.

Figure 20 .
Figure20.Same as Figure13but when the effect of magnetic turbulence is included in the polarization model and calculations of inclination angles.The 3D B-field strength inferred from our method is much closer to the actual strength of the mean B-field.

Figure 21 .
Figure 21.Relationship between p × S and the mass-weighted alignment efficiency, ⟨f align ⟩ for the different alignment models (left panel) and viewing angles γview (right panel).The solid lines show the running mean of the product p × S, while the dashed lines show the power-law fit to that running mean.In general, the product p × S has a good correlation with ⟨f align ⟩ with a slope of 0.4.However, the product of p × S strongly fluctuates for higher f align > 0.7, which is a result of the reduction of magnetic turbulence toward low-density regions and decreases significantly for lower γview.

Figure 23 .
Figure 23.Same as Figure 16 but when the beam convolution is taken into account, assuming the Ideal RAT alignment model and γview = 90 • .The synthetic polarization Psyn decreases with increasing the beam size and is lower than the analytical model Pana due to the fluctuations of B-fields within the beam.

Figure 24 .
Figure24.The maps of inferred inclination angles derived from Equation 31 when the beam convolution effect is considered.Since magnetic turbulence components within the beam are partially destroyed by beam convolution, the effect of magnetic turbulence becomes less significant, which can result in fewer NAN values with increasing beam sizes.

Figure 25 .
Figure 25.Histogram of the inferred inclination angles in comparison to the true inclination angles from the MHD data, considering the beam convolution effect.The inferred inclination angles become lower and less accurate when observed with larger beams.

Figure 26 .
Figure26.Left panel: Variation of the magnetic fluctuation angle δθ with the polarization angle dispersion σ ϕ using the unsharp-masking method for the different viewing angles.Right panel: same as the left one, but for the term F turb .Color dashed lines show the results of power-law fit.A correlation between δθ and σ ϕ or between F turb and σ ϕ is observed for all cases, except when the mean field is along the LOS.

Figure 27 .
Figure 27.The resulting polarization degree for various grain alignment models similar to Figure 5, but assuming the Astrodust model with the grain elongation of s = 1.4 (Draine & Hensley 2021).The viewing angle γview = 90 • is considered.As grains become less elongated, the overall polarization degree reduces by a factor of 2.

Figure 28 .
Figure 28.The maps of inferred inclination angles derived from our technique, but when the Astrodust model and lesselongated grains are taken into consideration.The inferred inclination angles do not change in comparison to the previous results in Figure 18.

Table 2 .
Grain alignment models considered.

Table 5 .
Same as Table4but when the effect of magnetic turbulence is taken into account.