Magnetic Field Orientation in Self-Gravitating Turbulent Molecular Clouds
Lucas Barreto-Mota, Elisabete M. de Gouveia Dal Pino, Blakesley Burkhart, Claudio Melioli, Reinaldo Santos-Lima, Luís H. S. Kadowaki
DDraft version January 12, 2021
Preprint typeset using L A TEX style emulateapj v. 12/16/11
Magnetic Field Orientation in Self-Gravitating Turbulent Molecular Clouds
Barreto-Mota, L. , de Gouveia Dal Pino, E. M. , Burkhart, B. , Melioli, C. , Santos-Lima, R. , and Kadowaki, L. H. S. Draft version January 12, 2021
AbstractStars form inside molecular cloud filaments from the competition of gravitational forces with turbu-lence and magnetic fields. The exact orientation of these filaments with the magnetic fields dependson the strength of these fields, the gravitational potential, and the line-of-sight (LOS) relative to themean field. To disentangle these effects we employ three-dimensional magnetohydrodynamical numer-ical simulations that explore a wide range of initial turbulent and magnetic states, i.e., sub-Alfv´enic tosuper-Alfv´enic turbulence, with and without gravity. We use histogram of relative orientation (HRO)and the associated projected Rayleigh statistics (PRS) to study the orientation of density and, inorder to compare with observations, the integrated density relative to the magnetic field. We findthat in sub-Alfv´enic systems the initial coherence of the magnetic is maintained inside the cloud andfilaments form perpendicular to the field. This trend is not observed in super-Alfv´enic models, wherethe lines are dragged by gravity and turbulence and filaments are mainly aligned to the field. ThePRS analysis of integrated maps shows that LOS effects are important only for sub-Alfv´enic clouds.When the LOS is perpendicular to the initial field orientation most of the filaments are perpendicularto the projected magnetic field. The inclusion of gravity increases the number of dense structuresperpendicular to the magnetic field, reflected as lower values of the PRS for denser regions, regardlessof whether the model is sub- or super-Alfv´enic. The comparison of our results with observed molecularclouds reveal that most are compatible with sub-Alfv´enic models.
Subject headings: galaxies: ISM — galaxies: star formation — ISM: kinematics and dynamics — stars:formation — turbulence INTRODUCTION
Molecular clouds (MCs) are the cradles of star forma-tion in our Galaxy. These are generally cold regions, withaverage temperatures often ranging between 10 and 50K,that exhibit filamentary structures produced by the in-terplay of supersonic magnetohydrodynamic (MHD) tur-bulence and self-gravity (Bergin & Tafalla 2007). Thedevelopment of new observational techniques and severaltheoretical advances have helped with the understandingof several aspects of MCs, but still several questions areleft unanswered (Burkhart et al. 2009, 2015a; Pattle &Fissel 2019; Girichidis et al. 2020).The presence of supersonic turbulence (with sonicMach number M s = v/c s >
1, where v is the tur-bulent velocity and c s the sound speed), initially in-ferred through line width observations (Larson 1981;Padoan et al. 1999), is one of the most important phys-ical ingredients in these environments and has beenextensively studied both theoretically (McKee & Os-triker 2007; Krumholz et al. 2007; Padoan & Nordlund2011; Hennebelle & Chabrier 2013; Federrath & Baner-jee 2015; Mocz et al. 2017) and observationally (Scalo& Elmegreen 2004; Bergin & Tafalla 2007; Andr´e et al.2010, 2014). Turbulence is also reflected in the struc- Department of Physics and Astronomy, Rutgers University,136 Frelinghuysen Rd, Piscataway, NJ 08854, USA Center for Computational Astrophysics, Flatiron Institute,162 Fifth Avenue, New York, NY 10010, USA Instituto de Astronomia, Geof´ısica e Ciˆencias Atmosf´ericasda USP, Dipartimento Astronomia di Bologna ture hierarchy of these clouds and consequently in thestatistical properties of the gas, such as the probabilitydistribution functions (PDFs) of the density and columndensity field and the power-spectrum of density and ve-locity (see, Landau & Lifshitz 1959; Lazarian & Pogosyan2008; Federrath et al. 2009; Collins et al. 2012; Federrath& Klessen 2013; McKee & Ostriker 2007; Hennebelle &Falgarone 2012). The motion of the fluid can be alsoaffected by the presence of magnetic fields since ionizedgas is also present and thus need to be accounted for tofully explain the dynamics of MCs (Padoan et al. 2000;V´azquez-Semadeni et al. 2011; Li et al. 2012; Burkhartet al. 2015b).The exact ratio between thermal pressure and mag-netic pressure (i.e., the plasma beta parameter) or theAlf´enic Mach number (i.e., the ratio between the velocityof the turbulence and the Alfv´en velocity, M A = v/v A )inside MCs are not well-known, but there are severaltechniques that provide some insight into the intensityand morphology of the magnetic field. For example, themagnetic field intensity along a line of sight (LOS) canbe estimated using the Zeeman effect, while its morphol-ogy in the plane-of-sky (POS) is more commonly ob-served using the polarized emission and absorption oflight from asymmetric dust grains (Heiles et al. 1993;Crutcher 2012). These grains are assumed to have theirsmaller axis aligned to the magnetic field, with the di-rection of the polarization of the thermal emission beingperpendicular to the magnetic field, and the polarizationfrom any background light that passes through the cloudand is absorbed by these aligned dust grains being paral-lel to the magnetic field. Observations are usually madein the visible and near-infrared spectrum, with more re- a r X i v : . [ a s t r o - ph . GA ] J a n Barreto-Mota et al.cent techniques also probing the polarization from thesubmilimiter spectra (Hull et al. 2013, 2014; Stephenset al. 2013, 2014; Davidson et al. 2014; Hull & Plambeck2015). The nature of this alignment is still a matterfor debate, but evidence points to the action of radiativetorques (RATS, see Hoang & Lazarian 2008; Hoang et al.2018, for further details). More recently, another methodhas been proposed with the idea of velocity anisotropiesin radio position-position-velocity channel maps, whichcan estimate the direction of the plane of the sky compo-nent of the magnetic field and also give a lower limit onthe intensity of the magnetic field (Burkhart et al. 2014;Esquivel et al. 2015; Lazarian et al. 2015; Kandel et al.2016).Recently, observations made by
Planck Telescope (Planck Collaboration et al. 2016),
Herschel andBLASTpol (Palmeirim et al. 2013; Soler et al. 2017),have revealed the relation between the filamentary struc-tures and the magnetic field of these regions. They havefound that denser structures usually appear perpendicu-lar to the magnetic field. Palmeirim et al. (2013), in par-ticular, have observed not only large filaments in Tau-rus MC perpendicular to the magnetic field as well asless dense striations of gas parallel to the magnetic fieldwhich seem to be flowing towards the filament. Simi-lar behavior has been observed by Planck Collaborationet al. (2016); (see also Soler 2019) in several other clouds,with less dense structures appearing mainly aligned tothe projected magnetic field in the plane of the sky.Previous studies have investigated the alignment be-tween magnetic fields and density structures using 3Dmagnetohydrodynamic (MHD) simulations. Soler et al.(2013) and Hull et al. (2017) analyzed the alignment be-tween structures and the magnetic field using statisticaltools like those we will employ in this work. They foundthat most of the dense filaments are nearly perpendicu-lar to the magnetic fields. This relation is also presentin the column density maps. Soler & Hennebelle (2017)have also investigated the same issue and concluded thatthe observed change in relative orientation between col-umn density structures and the projected magnetic fieldin the plane of sky, from mostly parallel at low columndensities to mostly perpendicular at the highest columndensities, would be the result of gravitational collapseand/or convergence of flows. Seifried et al. (2020), on theother hand, have recently found that a transition betweenstructures parallel to the magnetic field to structures per-pendicular to it not always can be seen, depending on theline of sight (LOS).At smaller scales (i.e., down to 1000 AU scales), Hullet al. (2017) and Mocz et al. (2017) also studied the mor-phology of the magnetic field around collapsed cores us-ing ALMA data and 3D MHD simulations. They foundthat only in very strongly magnetized systems there ispreservation of the field direction from cloud to diskscales and an hourglass-shaped field at scales smallerthan 1000 AU. G´omez et al. (2018) has also studied thecollapse inside filaments through 3D MHD simulations.They noted that the magnetic field around the filamentis primarily perpendicular to the structure and the col-lapse along the filament would later bend the magneticfield lines creating U-shapes.In this work, we extend upon these studies investi-gating how the alignment of density structures with the magnetic field is related to the presence of gravity andthe sonic and Alfv´enic Mach numbers associated to theturbulence and how projection effects may affect obser-vations for different lines of sight. For this aim, we per-formed and analysed several 3D MHD simulations withinitial conditions compatible with observed MCs, withand without the inclusion of self-gravity.This paper is structured as follows, in Section 2 we dis-cuss the simulated models, their initial conditions andthe statistical methods used to analyze our results. InSection 3 we present the results obtained from our mod-els. In Section 4 we compare these results with obser-vations published in the literature and other numericalstudies. In Section 5 we discuss and compare our resultswith precious works and, finally, in Section 6 we presentour conclusions. SIMULATION SETUP
In order to study the evolution of molecular cloudenvironments and how the initial stages of star forma-tion are affected by turbulence, magnetic fields and self-gravity, we consider three-dimensional MHD simulationsof two families of models. In the first set, the simula-tions only have the presence of MHD turbulence, withno self-gravity. These are the numerical simulations ofisothermal molecular clouds performed by Cho & Lazar-ian (2003) and Burkhart et al. (2009). The second setof models is composed by isothermal models with self-gravity. We run these new MHD simulations using amodified version of the code developed by (see Kowalet al. 2007; Le˜ao et al. 2009; Santos-Lima et al. 2010).Both sets of simulations have the similar initial condi-tions, as we will describe below.Our simulations solve the ideal MHD equations in theirconservative form: ∂ρ∂t + ∇ · ( ρ v ) = 0 (1) ∂∂t ( ρ v ) + ∇ · (cid:104) ρ vv + (cid:16) P + B π (cid:17) I − BB π (cid:105) = ρ g + F (2) ∂ B ∂t + ∇ · ( vB − Bv ) = 0 , ∇ · B = 0 (3)where ρ , v and B are density, velocity and magneticfield, respectively, P is the thermal pressure of the gas, g is gravity, I is the identity dyadic tensor and F is thesource term for the turbulence driving. We consider acold cloud with a temperature of 10K, which is of the or-der of the estimated temperature from cold dust emissioninside cold dark clouds.An ideal isothermal equation of state implies that thepressure can be written as P = ρc s , where c s = (cid:112) k b T /m is the isothermal sound speed of the gas, ¯ m is the averagemass of the gas given by ¯ m = µm H , with the meanatomic weight µ = 2 . ρ g = − ρ ∇ Ψ is the gravita-tional force due to self-gravity. This term is consideredonly in the second family of models of our study. Thegravitational potential Ψ obeys the Poisson equation: ∇ Ψ = 4 πGρ. (4)Finally, the source term F appearing on the RHS termagnetic Field Orientation in GMCs 3of Eq. 2 is responsible for the turbulent energy injection.Turbulence is driven solenoidally continuously at everytime-step. The forcing is calculated in Fourier spacearound a characteristic wavelength that defines the in-jection scale ( l inj = 1 / L for models with self-gravityand l inj = 1 / L for models without self-gravity, whereL is the size of the domain). The turbulent energy cas-cades down to the (numerical) dissipation scale withinone dynamical time, L/c s . Self-gravity is turned on afterone dynamical time, which is at least two turnover timesfor the models with lower Mach numbers.The first set of models, which neglects self-gravity, has ψ = 0 in the equations above. These simulations werebuilt using a third order accurate hybrid essentially non- oscillatory (ENO) scheme to solve the MHD equations(Cho & Lazarian 2003; Burkhart et al. 2009), while theset of models with self-gravity employed a total varia-tion diminishing (TVD) method (Kowal et al. 2009; Le˜aoet al. 2009; Santos-Lima et al. 2010). To solve the Pois-son equation, a multigrid method was used (see del Valleet al. 2015, for further details). Both codes are based onthe Godunov method and use a Runge-Kutta procedure(e.g. Londrillo & Del Zanna 2000; Del Zanna, L. et al.2003) for time integration. Boundary and initial conditions
We consider a cubic Cartesian domain with periodicboundaries. For the first family of models, the units oflength are given by the size of the injection scale ( L inj ).The rms velocity ( δV ) is kept close to unity so that thevelocity can be seen in units of δV , and B / (4 πρ ) / isthe Alfv´en velocity in the same units. The unit of time isgiven by the turnover time of the largest eddy, L inj /δV .Density is in units of the initial ambient density ρ . Theremaining quantity units are all derived from these ones.The simulations that consider self-gravity were per-formed using resolutions of 256 and 512 cells in thethree directions of a uniform grid. Our lower resolutionsimulations were mainly used for testing and calibrationof the initial conditions of the models that consider self-gravity. With that said, previous works have used similarresolutions of 256 to 512 for studying molecular cloudswith good statistical convergence of the results alreadyin 256 resolution, at least for non self-gravitating turbu-lent models (e.g. Kowal & Lazarian 2007; Burkhart et al.2009; Santos-Lima et al. 2010). The statistical behaviorof models with a resolution of 256 show very similar re-sults to their high resolution counterparts (see AppendixA).To evaluate the relative importance of turbulence com-pared to gravity, we will define the virial parameter, α vir ,following Bertoldi & McKee (1992): α vir ∼ E k E g = 5 v c s RGM , (5)where E k and E g are the kinetic and gravitational ener-gies of the system, respectively, v is the one-dimensionalrms velocity and M is the mass evaluated over a sphereof radius R .In the second set of simulations that consider gravity,units of length are given by the size of the domain ( L ),and units of velocity are given in terms of the isother-mal sound speed ( c s ), which implies time unit L/c s . The magnetic field is scaled such that B c.u. = B / (4 π ) / .Density is given in units of initial ambient density ρ .The initial setup is built to ensure the minimum condi-tion for collapse, such that the virial parameter (definedin Eq. 5) is α vir ∼ . M s (with M s = v/c s being sonicMach number and v being the characteristic velocity ofthe turbulence), the Alfv´enic Mach number M A and, inthe case of simulations with gravity, the virial parameter α vir (Eq. 5). These parameters are described in Table 1.The table also gives the corresponding initial thermal tomagnetic pressure ratio ( β ) for each model, the initialfree-fall time ( t ff = (3 π/ Gρ ) / ), calculated for theinitial ρ , for the simulations that consider self-gravity.We will use this time scale to compare the evolution ofdifferent self-gravitating models.Models without self-gravity are identified by the prefix Turb and models that consider self-gravity are identi-fied by the prefix
Grav . These simulations are drawnfrom the Catalog for Astrophysical Turbulence Simula-tions (CATS); (Burkhart et al. 2020). It is importantto note that there is a small difference in the parame-ters considered for the purely turbulent, non-gravitatingmodels and the ones that consider self-gravity. Some ofthe models without self-gravity have M A ∼ . M A ∼ . M s ∼ . M s ∼ .
8. However, the differences are so small thatthe comparison between them is not compromised.
General characteristics of the simulated models
In this subsection we present the general visual char-acteristics of the models used in this study.In Figure 1 we show the three-dimensional views of theself-gravitating models before (left) and after (right) self-gravity is turned on. On the top of Figure 1, we showthe supersonic, sub-Alfv´enic model
Grav Ms7.0 Ma0.6 and in the bottom we show the supersonic super-Alfv´enicmodel
Grav Ms7.0 Ma2.0 . When self-gravity is included(right-diagrams), fragmentation and filamentary struc-ture formation is enhanced and the collapse of the dens-est regions (clumps) of these filaments leads to star for-mation.A closer inspection of the models of Figure 1 showsthat the distribution of the magnetic field lines is deter-mined mainly by turbulence, and this effect is more pro-nounced in the case of super-Alfv´enic turbulence. Whatdetermines if the lines become parallel or perpendicularto a given density structure, or if they are twisted insidethe domain is whether the turbulence is sub-Alfv´enic orsuper-Alfv´enic. The twisting of the lines is more pro-nounced in the super-Alfv´enic case due to larger turbu-lent energy relative to the magnetic energy. This also af-fects the filament structure, which appears more chaoti-cally distributed with respect to magnetic field lines whenthe magnetic field is weaker. In the case of the sub-Alfv´enic model, we see that the lines are less distortedby turbulence and later on, when self-gravity becomesdominant, most of the filaments seem to be nearly per-pendicular to the magnetic field lines (right top panel).In this case, collapse is not inhibited along the lines. Barreto-Mota et al.
TABLE 1Initial conditions for all simulated models, with and without self-gravity
Simulation M s M A n ( cm − ) β = P th P mag t ff ( M yrs ) Resolution Turbulence GravityTurb Ms2.0 Ma0.7 2.0 0.7 117.65 0.302 ——– 256 yes noTurb Ms4.0 Ma0.7 4.0 0.7 444.81 0.080 ——– 256 yes noTurb Ms7.0 Ma0.7 7.0 0.7 1779.25 0.020 ——– 256 yes noTurb Ms2.0 Ma2.0 2.0 2.0 117.65 2.469 ——– 256 yes noTurb Ms4.0 Ma2.0 4.0 2.0 444.81 0.653 ——– 256 yes noTurb Ms7.0 Ma2.0 7.0 2.0 1779.25 0.163 ——– 256 yes noGrav Ms1.8 Ma0.6 1.8 0.6 117.65 0.302 7.23 512 yes yesGrav Ms4.0 Ma0.6 4.0 0.6 444.81 0.080 3.72 512 yes yesGrav Ms7.0 Ma0.6 7.0 0.6 1779.25 0.020 1.86 512 yes yesGrav Ms1.8 Ma2.0 1.8 2.0 117.65 2.469 7.23 512 yes yesGrav Ms4.0 Ma2.0 4.0 2.0 444.81 0.653 3.72 512 yes yesGrav Ms7.0 Ma2.0 7.0 2.0 1779.25 0.163 1.86 512 yes yes
In Figures 2 and 3 we present 2D (column-density)maps of super-Alfv´enic models with M A = 2 .
0, in twodifferent snapshots. Figure 2 compares models with dif-ferent sonic Mach numbers, while Figure 3, comparesmodels with different LOS. The models in Figure 2 con-sider M s = 1 . , . .
0, as indicated on the top ofeach column density map. In the left column of thisFigure, we show the initial snapshot t = 0 . t ff , cor-responding to the time when the turbulence has com-pletely evolved throughout the domain and has reacheda steady state regime, before self-gravity is turned-on.In the right column, we show the final snapshot for thesame models, when self-gravity becomes dominant. Wenotice an increase in the formation of filamentary struc-tures both with increasing M s and with the introductionof self-gravity.In Figure 3 all the models have initial M s = 7 .
0, andwe show the initial snapshot (with no self-gravity) in thetop row and the final snapshot (when self-gravity hasbecome dominant) in the bottom.In Figure 2 all column density maps were integratedalong a line of sight (LOS) perpendicular to the initialmagnetic field. Other LOS tests are presented in Figure3, but they all show very similar characteristics. Noneof the LOS have any distinctive characteristic, as onemight expect for super-Alfv´enic turbulence. The direc-tion of the projected magnetic field in the plane of thesky ( B ⊥ ) is also shown and has been produced using alinear integral convolution method (LIC, Cabral & Lee-dom 1993).Figure 4 presents the column-density maps at ini-tial snapshots (with fully developed turbulence and noself-gravity) of sub-Alfv´enic models with M A = 0 . M s = 1 . , .
0, and 7 .
0, respectively). From left to right,the Figure depicts maps integrated along different LOS,namely, perpendicular, at an angle of 45 ◦ , and parallelto the initial direction of B , respectively.Figure 5 shows the final snapshot (when self-gravityhas become dominant) for the same sub-Alfv´enic mod-els as in Figure 4, for comparison. Different from whatwe see in the super-Alfv´enic models, the LOS here isimportant, changing the distribution of the B ⊥ , and in- fluencing observed filaments.The comparison between sub-Alfv´enic and super-Alfv´enic models shows that both depict dense filamentsseparated by diffuse interstellar gas, but these under-dense structures when seen integrated along a LOS, ingeneral, seem to be larger (more coherent) in the super-Alfv´enic models with no self-gravity. This effect is morepronounced in the LOS perpendicular to the initial B (see Figure 2 left and top of Figure 3). Supersonic tur-bulence leads to shocks and compression of gas and mag-netic field lines, particularly in the early phase, beforeself-gravity sets in. In the super-Alfv´enic case, these ef-fects are more efficient in the building-up of large struc-tures because the magnetic field strength is smaller thanin the sub-Alfv´enic models. In the latter, stronger mag-netic pressure gradients offer larger resistance to the ac-cumulation of the overdense structures by shock com-pression. When self-gravity becomes important, frag-mentation and collapse will eventually dominate over thesupport provided by magnetic fields and turbulence inboth cases, but the general imprints left earlier in theformation of the large scale filaments by turbulence andmagnetic fields remain. General Statistics of the Simulations: Gravityversus Turbulence
Supersonic turbulence in MCs can both enhance andinhibit over-dense regions that may eventually achievesufficient conditions for collapse (Melioli et al. 2006; Le˜aoet al. 2009; Mocz & Burkhart 2019). It is possibleto study the interaction between turbulence and grav-ity through the use of statistical tools such as the onepoint and two point statistics, e.g., probability distribu-tion functions (PDF) and the power spectrum (Collinset al. 2012), respectively.In this subsection, we investigate the one and two pointstatistics of our two sets of simulations. In this waywe are able to benchmark them relative to theoreticalexpectations regarding the interaction and transition ofturbulent supported regions to self-gravitating collapsingregions.
Density PDFs agnetic Field Orientation in GMCs 5
Fig. 1.—
Comparison of the 3D distribution of density, represented by the iso-contours in light blue and purple, and magnetic fields,represented by the orange lines, for the initial ( t = 0 . t ff ; without self-gravity) and final snapshots ( t = 0 . t ff ; with gravity) for models Ms7.0 Ma0.6 grav (top) and
Ms7.0 Ma2.0 grav (bottom). The orange tubes represent the lines of
Due to interacting independent shock events, theisothermal density PDF in turbulent regions is well rep-resented by a lognormal distribution (Vazquez-Semadeni1994; Padoan et al. 1997; Scalo et al. 1998; Burkhart2018): p ( s ) = 1 (cid:112) πσ s exp (cid:18) − ( s − s ) σ s (cid:19) , (6)where s ≡ ln( ρ/ρ ), σ s is the standard deviation of thelognormal, s gives the value of s for the mean density and can be related to the width as: s = − σ s . Theturbulent sonic Mach number is related to width of thelognormal as: σ s = ln[1 + b M s ] , (7)where b is a dimensionless turbulent forcing parameter(Federrath et al. 2008) related to the solenoidal and com-pressive modes of the turbulence. For purely solenoidal For solenoidal turbulence driving, ∇ · δv = 0. For compressiveturbulence driving, ∇ × δv = 0 (Federrath et al. 2008). Barreto-Mota et al.
Fig. 2.—
Column density maps with LIC method applied to B ⊥ for super-Alfv´enic models with M A = 2 . M s = 1 .
8, and 4 . t = 0 . t ff ), representative of models with fully devel-oped turbulence without self-gravity. The right column presentsthe final snapshot of each model, at which self-gravity has becomeimportant. In all panels we show the column density distributionintegrated along X (the direction perpendicular to the initial field). turbulence driving, we have b = 1 / b = 1 . p P L ( s ) ∝ exp( − αs ) for s > s t , where s t is the transitional nor-malized density value between the diffuse gas given bythe lognormal distribution and the power-law tail. Ifthe piecewise PDF is continuous and differentiable than(Burkhart et al. 2017): s t = 12 (2 | α | − σ s , (8)where α is the power-law index.The transition density between lognormal and power-law is related closely to the density at which gravity takesover the dynamics of the cloud, i.e., the so-called criti-cal density for collapse (Burkhart & Mocz 2019). Weconsider that above a critical density ρ c all matter con-tributes to star formation. To estimate ρ c we consider themodel presented in Padoan & Nordlund (2011), whichtakes into account the contribution of the magnetic fieldfor the determination of the critical density for star for-mation and can also be used for sub-Alfv´enic turbulentconditions. The critical density is defined as: ρ c ρ = 0 . ζ − α vir M s (1 + 0 . β − / ) / (1 + β − ) / , (9)where ζ ≤ ζL is the turbulence integral scale ( L being the size of the system), and β is the ratio betweenthermal pressure and magnetic pressure. In the limit β → ∞ , we recover the hydrodynamical case. Figure 6 compares the evolution of the PDF of thedensity (Eq. 6) for two models analyzed in the pre-vious figures. The sub-Alfv´enic model is presented inthe left and the super-Alfv´enic in the right panel. Theblue dark curve represents t = 0 . t ff . Initially, thesuper-Alfv´enic model shows a wider spread of densityvalues when compared to the sub-Alfv´enic model (as inV´azquez-Semadeni & Garc´ıa 2001; Burkhart et al. 2009).The black dotted line is the fitted lognormal PDF (Eq.6) for t = 0 . t ff and the red dashed-dotted line is thepower-law fitted for the same time. The fitted index,and the region where the fit was considered, are indi-cated in the plot. The magenta vertical dashed line isthe estimated s t (Eq. 8) for the fitted power-law indexand the green vertical dashed line is the critical densityfor the magnetized case (Eq. 9). We note that both co-incide in the case of the super-Alfv´enic model. As thematerial collapses, the power-law becomes shallower andeventually approaches ρ c (Burkhart 2018).The sub-Alfv´enic model shows a small tail that devi-ates from the lognormal distribution, which is not presentin the super-Alfv´enic case. Since at this stage, gravity isnot acting in the system, this deviation from the lognor-mal is most likely caused by the presence of a strong mag-netic field (Burkhart & Lazarian 2012). At later times,the power-law reflects the action of gravity, as discussedpreviously in the literature (Ballesteros-Paredes et al.2011; Collins et al. 2012; Girichidis et al. 2014; Myers2015; Burkhart et al. 2017; Mocz et al. 2017; Burkhart2018).Since our code does not support a treatment of adap-tive mesh refinement or sink particles, we consider thatour results are valid until the power-law index reaches avalue α ∼ .
5. This has been chosen in accordance to pre-vious studies that indicate power-law tails from observedclouds with an index up to this value (see e.g., Table 1from Burkhart 2018, and references therein). In fact, theevolution of the cloud should result in a power-law withan index that converges to −
1. However, effects due tothe lifetime of the MC or to the LOS may yield steepervalues for the observed power-law index (Girichidis et al.2014; Guszejnov et al. 2018).
Power Spectrum
The power spectrum can provide additional informa-tion about the development of turbulence present in thecloud. An incompressible fluid with fully developed hy-drodynamical turbulence follows a Kolmogorov powerspectrum ( P D ( k ) dk ∝ k − / dk ). However the slope ofthe spectrum may change in the presence of magneticfields, shocks and gravity.The collapse of structures due to the action of grav-ity produces very shallow slopes that may even becomepositive valued, similar to a delta-function (Federrath &Klessen 2013; Burkhart et al. 2015a).Figure 7 shows the evolution of the 1D power spec-trum of density evolution in time for the same twomodels analyzed in Figure 6. The dotted red line is areference to the expected power-law from a Kolmogorov The power spectrum was calculated along the x-axis (perpen-dicular to the initially homogeneous magnetic field). The powerspectrum obtained along the other axes showed a very similar be-havior. agnetic Field Orientation in GMCs 7
Fig. 3.—
Column density maps with LIC method applied to B ⊥ for super-Alfv´enic models with M A = 2 . M s = 7 .
0, for three different LOS. Top row shows the column density map at the initial snapshot ( t = 0 . t ff ) (with fully developedturbulence and no self-gravity), for each LOS. Bottom panel shows the final snapshot (when self-gravity becomes important) for the sameLOS. From left to right the column density distribution is integrated along X (the direction perpendicular to the initial field), 45 ◦ withregard to the initial field, and Z (parallel to the initial field) directions. cascade ( P D ( k ) dk ∝ k − / dk ). Initially, in both sim-ulations it is possible to identify an inertial region inthe spectrum that roughly follows the same Kolmogorovslope. As matter accumulates around overdense regions,the power spectrum at higher wavenumbers k (or smallerlength scales) starts to flatten. This is a confirmation ofthe action of gravity at later times of the evolution of thecloud. QUANTIFYING THE RELATIVE ORIENTATIONBETWEEN DENSITY STRUCTURES AND THEMAGNETIC FIELD
In the previous sections we have demonstrated that oursimulations show statistical and visual behavior compati-ble with previous works which investigated the transitionof supported turbulent gas to collapsing star forming gas.The primary aim of this work is to quantify the role thatgravity and supersonic turbulence play in shaping therelationship between the magnetic field and collapsingdense structures. We focus on quantifying the impor-tance of gravity vs. turbulence in the relative orientationbetween ISM density/column density structures and theembedded magnetic field.The distribution or histogram of relative orientation(HRO) between the magnetic field and density is de-scribed in detail in Soler et al. (2013). The basics ofthe statistics are as follows: Gradients of density will be perpendicular to isodensity contours. Comparing thedensity gradient with the magnetic field orientation, it ispossible to evaluate how the filamentary structures arealigned with the magnetic fields.The relative orientation can be defined as:tan φ = B × ∇ ρ B · ∇ ρ . (10)With φ as defined above, we can evaluate the histogramof cosφ within different density bins. A peak in thishistogram around cosφ = 0 means that the field is per-pendicular to the density gradient, in other words, thefield is aligned with the structures. Similarly, cosφ = ± φ + 90 ◦ in this case. In three dimensions, two randomvectors have a higher probability of being perpendicularto each other than being parallel, that is why we chooseto use cos φ in 3D. This is not the case in two dimensionsand we can use simply the information from φ . This isalso what is commonly used in observations, e.g. Planck(Planck Collaboration et al. 2016; Soler et al. 2017). Barreto-Mota et al. Fig. 4.—
Column density maps with LIC method applied to B ⊥ for sub-Alfv´enic models with M A = 0 . M s , at the initialsnapshot (with fully developed turbulence and no self-gravity). From top to bottom panels we show models with M s = 1 . , . . ◦ with regard to the initial field and Z (parallel to the initial field) directions. For observational data, the angles are calculated be-tween the gradient of column density vector and the es-timated B ⊥ from polarization data obtained from theStokes parameters, i.e., a set of values that characterizethe polarization of electromagnetic waves.In order to compute these parameters for polarized ra-diation from the simulations, we will use the same as-sumptions as Falceta-Gon¸calves et al. (2008); (see alsoPoidevin et al. 2015). We assume that only thermal emis-sion is emitted by grains that are perfectly aligned withthe magnetic field. The dust abundance (which is notexplicitly accounted for in our simulations) is consideredto be proportional to the gas density, as well as the inten- sity of its emission. Finally, we consider that all grainsemit at a single temperature. With these assumptions,for each cell of the computational domain we calculate: q = ρ cos 2 ψ sin i,u = ρ sin 2 ψ sin i, (11)where ρ is the local density, ψ is the local angle of align-ment, determined by the projection of the local magneticfield onto the plane of the sky, and i is the angle betweenthe magnetic field and the line of sight (LOS). Integrat-ing q , u and ρ along a chosen line of sight results on theStokes parameters Q , U :agnetic Field Orientation in GMCs 9 Fig. 5.—
Column density maps with LIC method applied to B ⊥ for the same sub-Alfv´enic models presented in Figure 4, but in anevolved time when self-gravity has become dominant. The snapshot time for each map is indicated on the top of each diagram. From top tobottom, we show M s = 1 . , . .
0, respectively. From left to right, we show the LOS integration along X (the direction perpendicularto the initial field), 45 ◦ with regard to the initial field and Z (parallel to the initial field) directions. Q = (cid:90) q dl,U = (cid:90) u dl, (12)and the column density: N H = (cid:90) ρ dl. (13)This way, the intensity of the polarization vector E and its direction will be calculated as: E = (cid:112) Q + U , Φ B = 1 / U/Q ) . (14)We can use Equation 10, where ∇ ρ can also be re-placed by ∇ N H and B by B ⊥ if we use column densityand polarization maps, to build a histogram to evaluatethe angular distribution.To analyze the behavior of filaments inside molecularclouds we need to evaluate a wide range of density val-ues that spans 2 to 3 orders of magnitude. To checkhow the alignment between density structures and themagnetic field occurs at different scales, we analyze thedensity information in several bins. To guarantee com-0 Barreto-Mota et al. Fig. 6.—
Time evolution of the density PDF for the super-Alfv´enic model
Grav Ms7.0 Ma0.6 (left) and the sub-Alfv´enic model
Grav Ms7.0 Ma2.0 (right). The green vertical dashed line indicates the critical density for star formation (Eq. 9). The magenta ver-tical dashed line represents the transition density as defined in Eq. 8. The black dashed line and the red dashed-dotted line are the fittedlognormal and power-law to t = 0 . t ff , respectively. Fig. 7.—
Power spectrum of 3D density for models with M s = 7 . M A = 0 . M A = 2 .
0. The red dashed line represents the Kolmogorov power-law ( k − / ) and the green dashed line represents thecase for Burgers Turbulence ( k − ), for reference. parable statistics for each density bin, the density rangeis divided into bins with the same number of grid cellsand then the HRO is calculated for each bin.In order to study the relative orientation betweenthe (column) density gradient and the (projected) mag-netic field, an additional statistical method can be used,namely, the projected Rayleigh statistics or PRS. This can be calculated as (Jow et al. 2018): Z x = Σ ni cos θ i (cid:112) n/ , (15)where θ i ∈ [ − π, π ] is the set of angles between the twovector quantities that we want to characterize and n isagnetic Field Orientation in GMCs 11the total number of angles in our set. Positive valuesof Z x are indicative of strong parallel alignment betweenthe two vectors, while negative values indicate a strongperpendicular alignment between them.Jow et al. (2018) argue that, in the limit of n → ∞ ,the PRS approximates the standard normal distribution.Therefore, for a general distribution of angles, the vari-ance of Z x can be estimated as σ Z x = 2 Σ ni (cos θ i ) − ( Z x ) n . (16)The error of each measurement of Z x will be given bythe equation above.We differentiate between the PRS applied to the den-sity distribution of our 3D simulated cube and the PRSapplied to the integrated simulated column density mapsalong a given line of sight. When analyzing a 3D distri-bution, we will calculate the PRS values for the anglesbetween the gradient of density ( ∇ ρ ) and the magneticfield ( B ). In this case, Eq.15 will be referred as Z D .For the integrated density maps along a given line ofsight (column density), what is being evaluated is theangle between the polarization pseudo-vector ( E , insteadof the magnetic field projected on the plane of the sky, B ⊥ ) and the gradient of column density ( ∇ N H ). In thiscase Eq.15 will be referred as Z D .It is important to clarify the meaning of positive andnegative values of Z D and Z D . For the integrated col-umn density maps, Z D estimates the average angle be-tween the gradient of column density, ∇ N H , and thedirection of the polarization pseudo-vector, E . Given anisocontour of column density, since ∇ N H is perpendic-ular to it and E is perpendicular to the projected mag-netic field in the plane of the sky, B ⊥ , the angle betweenthem, φ , will be also between B ⊥ and the direction ofthe density structure (iso-contour). With this in mind,if E is parallel to the gradient of column density ∇ N H ,then Z D >
0, and if perpendicular, Z D < Z D , we are directly comparing the mag-netic field direction ( B ) with the density gradient ( ∇ ρ ),so positive values of Z D indicate parallel alignment be-tween the two vectors while negative values indicate thatthe two vectors are perpendicular. Results of the PRS and HRO
In this work, all the PRS calculations considered 20bins, both for Z D and Z D . Tests with additional binsdid not add any relevant information. As we will see,the information provided by both criteria Z D and Z D must be seen as complementary to each other to pro-vide a whole picture of the relative distributions of thestructures and their magnetic fields. Models without self-gravity
Figure 8 shows the PRS analysis for all the sub-Alfv´enic models and for the three different lines of sight.First, when the LOS is perpendicular to the initial uni-form field, Z D has a broader variation in time for higherdensities.Positive values of Z D indicate that E , the directionof the polarization vector, is parallel to ∇ N H , which in-dicates that the projected magnetic field in the plane ofsky ( B ⊥ ) is parallel to the iso-contours of N H . When the LOS is parallel to the initial magnetic field, the PRSanalysis returns only positive values for most densitieswith very little variation. This happens because B ⊥ inthis case results from motions perpendicular to the maincomponent of the field and this results in a random fielddistribution as seen in Figures 4 and 5. This is the casefor the maps on the right column diagrams from Fig-ure 5, which show that the coherence length is smallerin this LOS. Since we are projecting only the plane ofthe sky component of the field, we do not see the maincomponent, only perturbations perpendicular to it. Fig. 8.—
PRS time evolution for all sub-Alfv´enic models (with M A = 0 .
7) without self-gravity. From left to right the PRS isapplied to the LOS along X (the direction perpendicular to theinitial field), 45 ◦ with regard to the initial field and Z (parallelto the initial field) directions. From top to bottom initial M s =2 . , . , .
0, respectively.
In the other two LOS of Figure 8, as the sonic Machnumber is increased (keeping the same M A = 0 . less aligned with the local projected magnetic field to theplane of sky. This seems to be counter-intuitive to whatone should expect, since this effect seems to be larger forlarger turbulent motions (larger M s ) relative to the mag-netic field strength, where compression effects should beeven stronger. However, looking at the column densitymaps of Figures 2 through 5, we note that the increaseof turbulence (increase of M s ) causes more fragmenta-tion and the formation of more numerous smaller anddenser structures. This effect is more pronounced forthe LOS along X (where the projected magnetic field tothe sky has a larger component aligned to the originalmagnetic field) and less pronounced as we go to the LOSalong Z (where the projected magnetic field to the skyhas a smaller component aligned to the original magneticfield). In other words, in these sub-Alfv´enic models, onlythe densest and smallest structures that develop from2 Barreto-Mota et al.increased fragmentation in the more supersonic (larger M s ) at latter stages of evolution, tend to align with theintrinsic magnetic fields, and this effect is observable onlyfor LOS perpendicular or with angles around 45 ◦ to theoriginal field.Figure 9 shows similar PRS analysis as in Figure 8, butfor the super-Alfv´enic models. When M A = 2 . Z in the sub-Alfv´enic case. Therefore, in general, the mag-netic field presents itself parallel to the filaments in thesesuper-Alfv´enic, supersonic flows. Fig. 9.—
PRS time evolution for all super-Alfv´enic models (with M A = 2 .
0) without self-gravity. From left to right the PRS isapplied along LOS X (the direction perpendicular to the initialfield), 45 ◦ with regard to the initial field and Z (parallel to theinitial field) directions. From top to bottom initial M s = 2 . , . .
0, respectively
Figure 10 depicts the Z D analysis for the models with-out self-gravity. This analysis method indicates that theoverall gradient of density is mostly perpendicular to B (and hence, structures are actually mostly aligned withthe intrinsic magnetic fields), since Z D is negative forall densities. This applies to both set of simulations,sub-Alfv´enic and super-Alfv´enic. The sub-Afv´enic mod-els show slightly higher values of Z D towards denserregions compared to their super-Alfv´enic counterparts.These results indicate that the compression of the linesby the gradient of pressure is the dominant factor due tothe supersonic turbulence. When M A (cid:29) Fig. 10.— Z D analysis (Equation 15) applied to the 3D distri-bution of density and magnetic fields for all simulated models with-out self-gravity (see Table 1) in three different snapshots indicatedin the inset. On the left are presented the sub-Alfv´enic (initial M A = 0 . M A = 2 . M s = 2 . , . .
0, respec-tively. sure) and perpendicular to B , resulting the alignment ofthe magnetic field with the density filament. Now, as thesystem gets more sub-Alfv´enic, B will be more intenseand offer a greater resistance to compression.This Z D analysis seems to be, at first sight, a little incontradiction with the Z D analysis performed before forthe column density distributions, at least for the denseststructures. However, one has to have in mind that the Z D analysis is subject to projection effects, while the Z D analysis provides intrinsic values. Furthermore, thefact that no positive values of Z D appear does not meanthat there are no regions where ∇ ρ is parallel to B (i.e.,where dense filaments are normal to the intrinsic mag-netic fields). The first step to calculate the PRS (both fordensity and for column density distributions) is to dividethe density distribution into bins (density intervals) withthe same number of cells (see Soler et al. 2013), and then Z D (for density) and Z D (for column density) are cal-culated inside each bin. Since at very high densities wehave fewer structures, the last bin may include structuresthat have very different values of density and therefore,may have very different alignment with the field whichmay be hindered by averaging. As we will see below,the HRO analysis can help to distinguish the presence ofperpendicular filaments in the densest regions.agnetic Field Orientation in GMCs 13The different colors of the lines in Figure 10 indicate3 different snapshots at which each simulation was an-alyzed (t=500c.u., 600c.u., and 700c.u.). For both sub-Alfv´enic and super-Alfv´enic sets, as time goes by, thedensity variations are only caused by compression andrarefaction due to turbulent motions. Z D values do notchange much along time, but they do change with dif-ferent M A and M s values. For M A = 0 .
7, strongermagnetic fields become dynamically more important andforce the motion of the turbulent flow along the magneticfield lines, thus increasing the values of Z D at higherdensities, specially as the sonic Mach number increasescausing the formation of smaller denser structures. Thisdoes not occur at same degree in the super-Alfv´enic case.While the density distribution shows only a change inthe values of Z D as the Alfv´enic Mach number changes(Figure 10), for the column density along a given LOS,this behavior is realized only partially in Z D (in Figures8 and 9).To exemplify the distribution of angles at different den-sity bins, Figure 11 shows the histograms of relative ori-entations (HRO); (Eq. 10, Section 3) for two differentmodels, Turb Ms7.0 Ma0.7 (left) and
Turb Ms7.0 Ma2.0 (right). With φ being the angle between the density gra-dient and the local magnetic field (Eq. 10), cos( φ ) = 0means that B is perpendicular to the density gradientwhile cos( φ ) = ± φ shows a clear peak around zero for every bin, exceptfor the densest one, where the histogram is almost flat,meaning that there are more or less the same numberof structures parallel and perpendicular to the field linesat this density bin (this is compatible with our previ-ous analysis of Figures 8 and 4). At the same time,when we look to the super-Alfv´enic model, there is ahigher count of cos φ around zero, i.e. with the struc-tures mostly aligned to B at all densities, even for thedensest regions (which is also compatible with the previ-ous analysis).The change in alignment at different densities, spe-cially for the sub-Alfv´enic models, reflects in the PRSanalysis, as there are less negative values for Z D as wego to denser regions (Figure 10). However, these regionsare not numerous enough nor big enough to bring Z D to positive values. Models with self-gravity
In models with no self-gravity the turbulent motionis the main agent modifying the magnetic field distribu-tion. In this section we now investigate the effects ofself-gravity on the orientation of the magnetic field anddensity.Figure 12 shows the calculated values of Z D for allmodels. As the self-gravitating regions accrete, the den-sity gradient becomes less perpendicular to the magneticfield at denser regions (i.e., B tends to become moreperpendicular to these collapsing regions). This resultis similar to what was seen in the previous section (Fig-ure 10) when no self-gravity was present. The addition isthat, higher densities are achieved as time passes. As thefluid streams more easily along the magnetic field lines(since in the normal direction magnetic pressure gradi- ents will inhibit the motion and provide support againstgravity), dense structures will accumulate by gravity ac-tion mainly perpendicularly to the field direction. Thiscan happen for both, the sub and super-Alfv´enic cases,since at smaller scales magnetic fields become more andmore important, as they are brought along with thecollapsing regions, but is more pronounced in the sub-Alfv´enic cases.The slopes of the curves presented in Figure 12 alsobehave in a similar way as in Figure 10, changing for dif-ferent M s in the sub-Alfv´enic case, while in the super-Alfv´enic case, the slope does not change much. Onceagain, no positive values of Z D are seen, but the changeof the slope as time passes indicates that regions where ∇ ρ is parallel to B at higher densities become more im-portant at later times. This is of course, due to the actionof gravity creating collapsed regions to where the flowof matter converges. Compared to the models with noself-gravity where filaments are formed by compressionforces only, in the models with self-gravity we see thatlower density regions are still dominated by the interplaybetween turbulent motions and the magnetic field whilehigher density regions become dominated by an interplaybetween the action of gravity and magnetic fields.Similarly to the left panel of Figure 11, the right panelshows the HRO curves for two self-gravitating models, Grav Ms7.0 Ma0.6 and
Grav Ms7.0 Ma2.0 . This wasevaluated over the last output of these models, when t = 0 . t ff (see the right side of Figure 1 which showsthe filamentary 3D distributions for these models). Com-pared to the models that do not include self-gravity, thereis an enhancement in the number of regions perpendic-ular to the magnetic field at higher densities. With theaction of gravity, even the super-Alfv´enic model shows achange in the number of counts of cos φ = ± . Z D isstill negative to all values of density, but the influence ofgravity interacting with denser regions is clear.The critical density (vertical green dashed line in Fig-ure 12) seems to be related to the densest bins, for modelswith sonic Mach numbers M s = 4 . .
0. There is a jump in the values of density between the penultimateand last points in the diagrams that comes from the den-sity range considered in the bins (as can be seen in Figure11, the 20 th bin has a wider range of densities comparedto the 19 th ). In the super-Alfv´enic case, the penultimatepoint is very close to ρ c for all times considered, while inthe sub-Alfv´enic case the penultimate point approaches ρ c as the system evolves.Tracing the critical density of a system using the PRSis an interesting possibility, and the exact relation be-tween the two can be further explored following the evo-lution of the alignment between structures and the mag-netic field at smaller scales and at latter times of thecollapse. However due to the lack of an adaptive meshwith increasing resolution in the densest regions in ourmodels, this is out of the scope of this work.Still, the results discussed up to this point are very sim-ilar to what was shown in the previous section for modelswithout self-gravity. We note that self-gravity does in-crease the number of denser regions where B is parallelto the gradient of density in the sub-Alfv´enic models,but it has little or no effect when it comes to the simula-tions with M A = 2 .
0. At smaller scales this is probably4 Barreto-Mota et al.
Fig. 11.—
HRO Histograms of cos( φ ) (Eq. 10, Section 3) for different bins of density for the sub-Alfv´enic (left) and super-Alfv´enic (right)models without self-gravity. Only bins 1 st , 10 th , 17 th ,18 th ,19 th and 20 th are shown. not true, since we expect that magnetic fields should bebrought along with the fluid during collapse, at somepoint these cores must become sub-Alfv´enic and onceagain B would influence how the gas collapses. How-ever, the simulations do not have enough resolution tofollow the process up to this point.In order to compare the simulations with observationswe need once again to integrate along a defined LOS. Wewill follow same directions used in section 3.1.1.Figures 13 and 14 show the PRS analysis ( Z D ) for theintegrated density (column density) distribution alongthe three different LOS. As in Figures 8 (for initial M A = 0 .
7) and 9 (for initial M A = 2 . Z D for a LOS makingan angle of 45 ◦ with respect to the initial field, and theright column, is for a LOS parallel to the initial magneticfield. From top to bottom, each line has, respectively, M s ∼ . , . . α = 1 . Z D (black curve) is similar to the simulations withoutself-gravity. However, as the densest regions collapse, theeffect seen in Z D (Figure 12) is more pronounced for Z D .From the sub-Alfv´enic models that do not consider self-gravity (see Figure 8) we saw that for higher sonic Machnumbers, and when the LOS is not parallel to initial mag-netic field, Z D decreases as the density increases. Thatis exactly what is seen in the models from Figure 13 whent = 0.0 t ff . However, with gravity acting over the sys-tem all models evolve to a similar distribution of Z D ,with lower densities having positive ( B aligned with thefilaments) values and higher densities showing negativevalues ( B perpendicular to the filaments). This turns outto be the case even when the LOS is parallel to the initialfield (where we see less of the original magnetic field ori-entation and more of the random component, see the lastcolumn of Figure 13). In particular, for M s = 1 .
8, thistrend is not initially present, but the action of gravity re-agnetic Field Orientation in GMCs 15
Fig. 12.— Z D analysis (Equation 15) applied to the 3D distri-bution of density and magnetic fields for all simulated models withself-gravity (see Table 1), for different snapshots depicted in theinset. On the left we present the sub-Alfv´enic (initial M A = 0 . M A = 2 .
0) models.From top to bottom initial M s = 2 . , . .
0, respectively. Thegreen dashed line indicates the critical density for star-formation(Eq. 9). sults in negative values of Z D at later times for all LOS.In summary, self-gravity does affect the Z D distributionof dense regions of sub-Alfv´enic models more clearly thanin models with no gravity, as we should expect, and mayhelp to distinguish between different observed IS regions(see Section 4).For the super-Alfv´enic case ( M A = 2 .
0; see Figure14) the scenario is different, the polarization vector E appears mainly aligned to the column density gradient ∇ N H (i.e. Z D > Z D valuesare achieved. COMPARISON WITH OBSERVATIONS
To compare the diagnosis of the simulated models pre-sented in the previous Section with observations we useclouds observed by
Planck Satellite (Planck Collabora-tion et al. 2016) and BLASTPol (Soler et al. 2017). Theseobjects have been also analyzed by Jow et al. (2018).The
Planck Satellite observed thermal emission anddust polarization in 7 bands between 30 and 353 GHz.
Fig. 13.—
PRS time evolution for all sub-Alfv´enic (with M A =0 .
6) models with self-gravity. From left to right the PRS is ap-plied along X (the direction perpendicular to the initial field), 45 ◦ with regard to the initial field and Z (parallel to the initial field)directions. From top to bottom initial M s = 2 . , . . Fig. 14.—
PRS time evolution for all super-Alfv´enic (with M A = 2 .
0) models with self-gravity. From left to right the PRSis applied along X (the direction perpendicular to the initial field),45 ◦ with regard to the initial field and Z (parallel to the initialfield) directions. From top to bottom initial M s = 2 . , . .
0, respectively.
In the case of Planck Collaboration et al. (2016), theobservations were made using the High Frequency In-strument at 353 GHz. The molecular clouds reportedhave distances estimated between ∼
150 pc and ∼ ◦ , these dis-tances imply that the clouds have sizes up to ∼
100 pc.The sonic Mach number estimated for these regions mayvary. For instance, Polaris MC has M s varying between ∼ − M s ∼ ∼ −
30K (Kirk et al.2013). Our models are compatible with these values.Vela C observations were obtained using BLAST-Polto estimate the magnetic field direction (Galitzki et al.2014; Fissel et al. 2016; Gandilo et al. 2016), and
Hes-chel , to derive the column density maps (Hill et al. 2011).BLAST-Pol used three wavelengths centered at 250, 350and 500 µm . For the column density, Heschel , SPIREand PACS data were used, with observations made at160 (PACS), 250, 350 and 500 µm . Previous studiesplace the cloud around a 700 pc distance, with a totalmass of more than 10 M (cid:12) . Despite being a massivecloud, it is still at an early stage of evolution, and someauthors claim that only one or two O-type stars havebeen formed (Soler et al. 2017). The estimated temper-ature ranges between ∼ −
30K inside the cloud (Hillet al. 2011).
PRS and angular distribution analysis
Figure 15 shows the comparison of the PRS presentedin Section 3.1.2 with the results obtained for Vela C and afew of the clouds observed by Planck Collaboration et al.(2016).The PRS analysis is sensitive to the number of pointsin the sample. This means that directly comparing thevalues of Z D with the ones obtained in observations re-quires some caution and we discuss this later in this sec-tion. Nevertheless, our results offer some insight into themagnetic field and what kind of structure distribution wecan expect inside the observed MCs.First, we discuss the general distribution of the pro-jected magnetic field on the sky ( B ⊥ ) for most clouds. Ingeneral, B ⊥ appears parallel to a single direction insidethese molecular clods, with a smaller fraction being ran-domly distributed. Chamaleon-Musca and Aquila, forinstance, are good examples of this behavior (see Fig-ure 3, left, from Planck Collaboration et al. 2016). Asdiscussed in the previous section, this is a characteristicobserved mainly in our sub-Alfv´enic models.From the integrated maps obtained from our super-sonic, sub-Alfv´enic models (Figures 4 and 5), B ⊥ be-comes more chaotic as the LOS gets closer to the direc-tion of the mean magnetic field. The LOS then catchesthe effect of regions where locally the field can appear asa random twist, while still showing a general coherencefor most of the cloud. The presence of such a charac-teristic in the distribution of the observed B ⊥ is andindication that these clouds are sub-Alfv´enic.Evidence in favor of a sub-Alfv´enic description of theturbulence in these clouds is the PRS analysis presentedin Figure 15. Jow et al. (2018) have revised the datapresented in the works mentioned above in order to ap-ply the PRS method described in Section 3. This revealsthat the observed relation between density gradients andthe projected magnetic field in the sky is similar to thesub-Alfv´enic models we simulated, going from positive to negative as density grows, i.e. less dense regions appearmore aligned to the magnetic field, while dense regionsappear more perpendicular to it . For most clouds, theturbulent models approximately produce a behavior of Z D compatible with the observations. The different re-gions of Vela C, in particular, show higher column densityvalues and most of the structures appear perpendicularto the projected magnetic field, since they have a nega-tive tail of Z D at higher densities. This is only achievedin the sub-Alfv´enic models that consider self-gravity (seeFigure 15). Aquila, on the other hand, can be describedby a sub-Alfv´enic model without self-gravity. This con-clusion is supported by the PRS data and the informa-tion from the observed column density. Aquila does notshow large gradients of column density, which is char-acteristic of our turbulent models without self-gravity.Moreover, the cloud is known for having star formationactivity inside isolated pockets (Prato et al. 2008; Eiroaet al. 2008), and since the PRS data was calculated us-ing a bigger portion of the cloud, it makes sense thatthe observed PRS follows a distribution similar to ourturbulent models.Complementary information is provided in Figures 3and 4 from Planck Collaboration et al. (2016), thatpresent the observed integrated column densities of theclouds Aquila, Chamaeleon-Musca, Taurus and Ophiu-cus on the left side, and the histograms of relative ori-entations between the projected magnetic fields and thedensity gradient of the structures in these clouds. Com-paring with the HROs of our models in Figure 11, we notea similar behavior with the sub-Alfv´enic models, i.e., aswe go from less dense to denser structures, the relativeorientation between the magnetic field and the filamentsgoes from aligned to perpendicular, particularly in thecase of Chamaeleon-Musca and, at some extent, Taurus.Note that Figure 11 shows the HRO applied to the 3Dstructures and not to column density maps, still the be-havior is similar.One effect that can also influence the comparison ofthe PRS from observations with our simulated models isrelated to the resolution and field of view that is pos-sible to achieve with the telescopes. It is important tohighlight that due to self-gravity, the units considered inthe code are scale dependent. The size of the clouds inour simulations is 10 pc, while the clouds observed canextend up to hundreds of parsecs. This means that ourresults are more representative of sub-structures insidethe clouds and not so much of the global formation ofthe cloud and their surroundings. The PRS results re-ported in Jow et al. (2018) are for entire regions observedby the Planck satellite . Even though, our study presentsresults that are qualitatively comparable to the observa-tions and hence, one can always argue that due to theself-similarity nature of the turbulent clouds, the generalbehavior at very large scales does not differ much of thatin the intermediate scales inside the clouds, at least inscales where self-gravity is not dominant yet (Bialy &Burkhart 2020). To account for the scaling effect above,we may take the self-gravity, sub-Alfv´enic simulationswith M s = 1 .
8. The Poisson equation considers the nor-malized potential: Note that this trend appears consistently only in sub-Alfv´enicmodels. agnetic Field Orientation in GMCs 17
Fig. 15.—
Top: Z D (Eq. 15) calculated from observations (Planck Collaboration et al. 2016; Soler et al. 2017). Bottom: Z D (Eq. 15)calculated from our simulated models considering all LOS for models with M s = 4 and 7. The PRS when only turbulence is present is onthe left. The PRS calculated for the final snapshot of models that consider self-gravity is shown on the right. The Alfv´enic Mach numberof the models is indicated above each plot. ∇ Ψ (cid:48) = ∇ (cid:16) Ψ4 πG (cid:17) = ρ. (17)The normalization implies a gravitational constant incode units as follows: G c.u. = G N ρ (cid:18) Lc s (cid:19) , (18)where G N is the gravitational constant, ρ is the initialdensity of the simulation, and L the size of the domain. With this in mind, it is possible to re-scale our modelsas long as we keep G c.u. the same, i.e. the ratio in theright-hand side of the equation, ρ ( L/c s ) , must be keptconstant (see alternative ways of scaling gravity keepingconstant the virial parameter in e.g., Mckee et al. 2010).We consider a region of 40 pc, which is approximatelythe estimated extension of the clouds observed by PlanckCollaboration et al. (2016). Considering the same tem-perature for the re-scaling (10K), the average density inthis larger system is around 6 cm − to keep the ratio8 Barreto-Mota et al. ρ ( L/c s ) constant. The final result is shown in Figure16.On the left side of Figure 16, the original column den-sity and PRS that were extracted from the simulationare shown for the re-scaled system. The integration wasmade along a direction 45 ◦ inclined with respect to theinitial magnetic field. On the bottom left diagram ofFigure 16, the red line is the obtained PRS from obser-vations of Chamaeleon-Musca (same as in the top dia-gram of Figure 15), and the black line is the PRS cal-culated for model Ms1.8 Ma 2.0 grav . The right-handside of the figure shows the same column density map ofthe left, but convolved with a 2D Gaussian kernel (seeSoler et al. 2013, for further details). The image hasbeen smoothed roughly to the same spatial resolution ofthe observations, as it was done by Planck Collaborationet al. (2016). In our simulation, the smoothing lengthcorresponds to about 7 cells. We note that the PRS forthe smoothed map is quite different at higher densities,even if the column density is not much different fromthe original one, thus providing different information. Infact, we see that the PRS calculated for the smoothedmap is more similar to the observed one (red curve). In-trinsically (i.e., examining the high resolution map of theleft), these structures are aligned to the magnetic fieldfor all column densities (explaining the positive PRS),but with the loss of resolution at the smaller scales ofthe observations (in the map of the right), the denseststructures actually appear perpendicular to the magneticfield. Both maps share similarities with Chamaeleon-Musca, but the alignment indicated by the PRS in thesmoothed map, while being compatible with the observa-tions, is actually not representative of the real behaviorof the projected magnetic field onto the sky (as suggestedby the simulation of the left panels).
A closer look inside the Molecular Clouds
Observations made by Palmeirim et al. (2013) using
Herschel
Telescope have revealed several smaller struc-tures around the filament B211/3 in Taurus molecularcloud. The dense filament (B211/3) appears perpendicu-lar to the magnetic field around it, while less dense struc-tures (the striations observed by Palmeirim et al. 2013)are parallel to the projected magnetic field onto the sky.According to the scale indicated in the Figure, the size ofthe region is about ∼ pc × pc . The separate regions ofVela C as defined in Soler et al. (2017) also have similarsizes.To evaluate the behavior of these smaller, denser re-gions, both of Taurus and Vela C, Figure 17 shows thetime evolution of the column density of two regions ofsimilar size to these clouds extracted from our self-gravitymodels with M s = 7 .
0. The density integration wasalong a LOS making an angle of 45 ◦ with respect to theoriginal magnetic field. On the left, we have M A = 0 . M A = 2 . t = 0 . t ff ), it is possibleto see filaments both parallel and perpendicular to themagnetic field in the sub-Alfv´enic model. At this time,as only turbulence and magnetic fields are present, the compression motions tend to align the filaments with themagnetic lines (an effect which is more pronounced inthe super-Alfv´enic model on the right side of the figure),while the stronger magnetic fields imposed by the sub-Alfv´enic regime tend to oppose resistance to alignment,through their tension and pressure gradient forces. Astime passes, matter flows along the lines and results indenser filaments perpendicular to the field. In the super-Alfv´enic case, the magnetic field is dragged with the flowresulting in projected magnetic fields aligned to the fila-ments. Along time, it is also possible to see dense regionswhere the field is perpendicular to the structures, whichis evidenced in Figure 18 for t = 0 . t ff , with Z D < Z D obtained in our sub-Alfv´enic models seems a better representation of what isobserved in the sky. Note that the values of Z D forthe super-Alfv´enic case when t = 0 . t ff and several re-gions have already collapsed, also show negative valuesfor higher column densities, which is expected as dis-cussed in Section 3.1, but the distribution for smallerdensities and at previous times is not compatible with theobservation. On the other hand, the sub-Alfv´enic modelshows very similar behavior when compared to the obser-vations. In particular, the PRS analysis obtained fromthe sub-Alfv´enic model in Figure 15 resembles what isseen in Vela C.Our models were made considering only solenoidaldriving turbulence. However as discussed in Section 1,compressive modes may change the distribution of thefilaments and this may affect the resulting PRS analysis.We will study this in a forthcoming work. DISCUSSION
In this section we compare our results with previousstudies and summarize our findings.
Comparison of our results with previous studies inthe literature
Soler et al. (2013) have first analyzed the alignmentof structures with the magnetic fields in a molecularcloud using statistical tools like those employed in thiswork. However, in their study turbulence was not con-stantly driven in the simulated system, and thus wasallowed to decay with time. Also, they did not con-sider sub-Alfv´enic models and investigated only a singlesonic Mach number. Their highest magnetized modelhad M A = 3 .
16 and M s = 10, so that they could notinvestigate most relevant dynamical effects of the mag-netic fields in the evolution of star forming systems, as inthe present work. Compared to our most similar model, Ms7.0 Ma2.0 grav , their results present significant dif-ferences. They only consider a LOS perpendicular tothe initial magnetic field and find a distribution of Z D that is closer to our distributions for sub-Alfv´enic mod-els. Also, the general coherence of B ⊥ in their columnagnetic Field Orientation in GMCs 19 Fig. 16.—
Left: Column density map along a LOS 45 ◦ inclined with respect to the initial magnetic field for our model Ms1.8 Ma2.0 grav (top) and the PRS calculated for the respective map (bottom). Right: Same map, but convolved with a Gaussian kernel in a similarprocesses to the one made for observations from Planck Collaboration et al. (2016) (top) and the PRS of the respective map (bottom);(seealso Soler et al. 2013, for further details). density maps is only observed in our sub-Alfv´enic mod-els. This difference is most likely due to the fact thatturbulence was not continuously driven. As turbulencedecays their system approaches the sub-Alfv´enic regime,with gravity and magnetic pressure becoming the mainforces acting on the fluid. Hence, the collapse primarilyoccurs along the field lines, resulting in dense structuresperpendicular to the magnetic field. As turbulence isdominated by the magnetic field it is unable to bend thelines, thus explaining why their results are more compa-rable with our sub-Alfv´enic models.In a more recent study, using the same set of simula-tions as in the work above, Soler & Hennebelle (2017)derived an expression for the time evolution of the an-gle between the density gradient and the magnetic field in the turbulent MCs and found that these two quan-tities are preferentially either parallel or perpendicularto each other. In our simulations, we identify a simi-lar trend only for the evolved sub-Alfv´enic models withself-gravity (see Figure 13). In addition, these authorshave concluded that the observed change in the relativeorientation between column density structures and theprojected magnetic field in the plane of sky, from mostlyparallel at low column densities to mostly perpendicularat the highest column densities, would be the result ofgravitational collapse and/or convergence of flows. Thistrend is also identified in our models, specially in thesub-Alfv´enic ones with self-gravity (see Figures 13 and18).Our results are also in accordance to those reported0 Barreto-Mota et al.
Fig. 17.—
Time evolution of zoomed-in regions extracted from models
Ms7.0 Ma0.6 grav (left) and
Ms7.0 Ma2.0 grav (right). Bothmaps were integrated along a LOS 45 ◦ inclined with respect to the initial magnetic field. The time considered is indicated above each map. by Seifried et al. (2020) who investigated only cloudswith no self-gravity. As in our models without gravity,they do not observe perpendicular structures to the meanfield, particularly in their lower magnetization models.As we have seen, a sub-Alfv´enic cloud observed in a LOSparallel to the magnetic field and a super-Alfv´enic cloudunder the action of gravity can yield similar results withregard to the general alignment between magnetic fieldand filaments.Hull et al. (2017) have also performed 3D MHD simula-tions in order to study the alignment at smaller scales inorder to compare with observations made by ALMA. As initial conditions, they considered a single sonic Machnumber ( M s = 10) and different Alfv´enic Mach num-ber cases, including a trans-Alfv´enic and a sub-Alfv´enicone. The general setup of their models is similar to theones presented in this paper, with column density mapsalso comparable to ours, but their analysis focus on theformation of collapsed cores. They have studied the mor-phology of the magnetic field around these cores for dif-ferent scenarios in order to compare to ALMA obser-vations. In our work however, we have focused on theformation of filaments and the interaction of these struc-tures with the magnetic fields. Both works can thus beagnetic Field Orientation in GMCs 21 Fig. 18.—
Time evolution of the PRS analysis ( Z D , Eq. 15) for the same maps presented in Figure 17. Solid lines represent the PRSalong time calculated from our models, sub-Alfv´enic on the left-hand side and super-Alfv´enic on the right hand side. Dashed lines representthe observed PRS of some regions presented in the top diagram of Figure 15. seen as complementary to each other, with the resultspresented here giving some insight about the birthplaceof cores that will ultimately form stars.G´omez et al. (2018) studied the structure of magneticfields inside self-gravitating filaments in turbulent envi-ronments. They note that the magnetic field around thefilament is primarily perpendicular to the structure andthe collapse along the filament would later bend the mag-netic field lines creating U-shapes. However, they arguethat the lack of resolution, as well as the decrease of po-larization in observations would not allow to detect thesefeatures. This is consist with the result we found in Fig-ure 16, where the decrease of resolution, simulated bythe smoothing of the image, revealed a complete changeof the behavior in the PRS.A final remark regarding the result from Figure 16 isin order. Though applied to Chamaeleon-Musca cloud,it has served to illustrate how re-scaling may influencethe results obtained from observations. In particular,the similarities seen in the general behavior of B ⊥ andthe formation of larger collapsed regions perpendicularto it, as well as the comparison of the PRS showing areasonable match in Figure 16, indicate that this issueshould be further explored in future works. CONCLUSIONS
In this work, we have explored, by means of 3D MHDsimulations of isothermal MCs, how the relative align-ment between density structures and the magnetic fieldsis affected by different regimes of MHD turbulence andby self-gravity. We have also examined how projection ef-fects of these structures on the plane of sky may be com-pared with polarization observations of distinct clouds.We considered models with initially homogeneous mag-netic field and density, and different sonic Mach numbers( M s ∼ . , . .
0) and Alfv´enic Mach numbers( M A ∼ .
7, and M A ∼ . pc × pc × pc vol-ume. To ensure the gravitational collapse of our cloudsin the presence of self-gravity, we considered an initial ratio between turbulent and gravitational energy density α vir ∼ . Z D or Z D , respectively.Computing the histograms of the relative orientationsbetween the magnetic fields and the density gradientsand applying the PRS analysis (Jow et al. 2018) en-ables better understanding of the interaction of den-sity structures with the magnetic field inside molecularclouds. Considering the results from the analysis of thecolumn density maps with models that do not considerself-gravity (Section 3.1.1) we have found that: • The LOS in the the sub-Alfv´enic model caseschanges the PRS distribution ( Z D ; Equation 15),and this is also affected by the magnitude of thesonic Mach number. For smaller sonic Mach num-bers ( M s ∼ .
0) most of the polarization vector E appears parallel to the column density gradient( ∇ N H ), Z D > Z D can be seen fordenser regions (specially if the LOS is not parallelto the initial magnetic field). This indicates thatthese denser smaller regions have a greater contri-2 Barreto-Mota et al.bution of structures perpendicular to B ⊥ , due tothe dominance of the magnetic forces that preventgas motion across the lines, facilitating their flowalong them. • For super-Alfv´enic models, the less intense mag-netic field does not show the same coherence alonglarge length scales, since the lines are more easilytwisted by the turbulent motions, and the columndensity maps show similar characteristics at differ-ent LOS. The PRS analysis also reflects this, with Z D presenting only positive values (often rangingbetween ∼ −
20) for all LOS. Different sonic Machnumbers also show no clear effect in the alignmentof the structures with the field.When gravity is included these scenarios change.For the models with self-gravity (Section 3.1.2) we havefound that: • There is an enhancement of dense structures per-pendicular to the projected magnetic field to thesky, even when looking along the mean field. Whenthe system is sub-Alfv´enic, smaller values of Z D are present for all LOS at higher column densityvalues. At later times, all models show at somedegree a change from positive to negative values in Z D which implies the presence of dense, collapsingstructures with magnetic fields normal to them, asone should expect, since the collapse is easier alongthe magnetic fields. • For the super-Alfv´enic models, Z D is positivefor most LOS, indicating that most structures arealigned with the projected magnetic field to the sky.This effect is less prominent for M s = 7 .
0, wheremore fragmentation again propitiates the forma-tion of denser collapsed cores inside filaments, andsince magnetic fields are aligned to them, thesecores appear as perpendicular structures to thefield and therefore yield lower values of Z D .We conclude that self-gravity can create structures per-pendicular to the magnetic field, even for mean super-Alfv´enic Mach numbers. As evidenced also by the PDFand power spectrum, dense regions of MCs are clearlydominated by gravity, while less dense regions are mainlyaffected by turbulence, as lower densities still sit in thelognormal branch of the PDF (see Figure 6). Effectsof projection due to the LOS may change the observedalignment for less dense regions (see Figure 14), but stilloverdense regions show smaller values of Z D (this effectcan be seen by comparing the bottom panel of Figure 14with the right panel of Figure 18).The observed behavior of the PRS in column densitymaps is a result of the intrinsic distribution of filamentsand magnetic fields inside the molecular clouds. Withregard to the 3D analysis of the density structures, wehave found that: • For the sub-Alfv´enic models, the inclusion of grav-ity helps the creation of structures perpendicularto the strong magnetic field, but supersonic turbu-lence favors the formation of less dense filamentsparallel to the field (see Figure 1). This is more eas-ily realized through the histograms of cos( φ ) (top diagrams of Figure 11), with cos( φ ) = ± Z D (Figures 4 and 5), thePRS analysis of the density shows higher values of Z D with higher sonic Mach numbers, specially athigher densities. • For the super-Alfv´enic models, magnetic field linesappear mostly aligned to filaments when only tur-bulence is considered, with the sonic Mach havinglittle or no effect. But models with gravity havethe densest structures perpendicular to the mag-netic field at later times.Finally, in Section 4 we have compared the results de-scribed above with observations made by
Planck , Her-schel and BLASTPol. The comparison indicates that,qualitatively, our sub-Alfv´enic models can better repro-duce the characteristics of observed clouds. Not only thebehavior of the observed Z D , but also the general co-herence of the magnetic field projected on the plane ofthe sky ( B ⊥ ), is compatible with our sub-Alfv´enic mod-els for most clouds. There are clouds where twists of B ⊥ DATA AVAILABILITY
The data of this article will be available upon requestto the corresponding author Note that positive values of Z D should be reflected as negativevalues of Z D , and vice-versa (see Section 3). agnetic Field Orientation in GMCs 23 REFERENCESAndr´e, P., Di Francesco, J., Ward-Thompson, D., et al. 2014, inProtostars and Planets VI, ed. H. Beuther, R. S. Klessen, C. P.Dullemond, & T. Henning, 27Andr´e, P., Men’shchikov, A., Bontemps, S., et al. 2010, A&A,518, L102Ballesteros-Paredes, J., Hartmann, L. W., V´azquez-Semadeni, E.,Heitsch, F., & Zamora-Avil´es, M. A. 2011, MNRAS, 411, 65Bergin, E. A., & Tafalla, M. 2007, ARA&A, 45, 339Bertoldi, F., & McKee, C. F. 1992, ApJ, 395, 140Bialy, S., & Burkhart, B. 2020, ApJ, 894, L2Burkhart, B. 2018, ApJ, 863, 118Burkhart, B., Collins, D. C., & Lazarian, A. 2015a, ApJ, 808, 48Burkhart, B., Falceta-Gon¸calves, D., Kowal, G., & Lazarian, A.2009, ApJ, 693, 250Burkhart, B., & Lazarian, A. 2012, ApJ, 755, L19Burkhart, B., Lazarian, A., Balsara, D., Meyer, C., & Cho, J.2015b, ApJ, 805, 118Burkhart, B., Lazarian, A., Le˜ao, I. C., de Medeiros, J. R., &Esquivel, A. 2014, ApJ, 790, 130Burkhart, B., & Mocz, P. 2019, ApJ, 879, 129Burkhart, B., Stalpes, K., & Collins, D. C. 2017, ApJ, 834, L1Burkhart, B., Appel, S., Bialy, S., et al. 2020, arXiv e-prints,arXiv:2010.11227Cabral, B., & Leedom, L. C. 1993, Imaging vector fields using lineintegral convolution, Tech. rep., Lawrence Livermore NationalLab., CA (United States)Cho, J., & Lazarian, A. 2003, MNRAS, 345, 325Collins, D. C., Kritsuk, A. G., Padoan, P., et al. 2012, ApJ, 750,13Crutcher, R. M. 2012, ARA&A, 50, 29Davidson, J. A., Li, Z. Y., Hull, C. L. H., et al. 2014, ApJ, 797, 74del Valle, M. V., Romero, G. E., & Santos-Lima, R. 2015,MNRAS, 448, 207Del Zanna, L., Bucciantini, N., & Londrillo, P. 2003, A&A, 400,397Eiroa, C., Djupvik, A. A., & Casali, M. M. 2008, The SerpensMolecular Cloud, ed. B. Reipurth, Vol. 5, 693Esquivel, A., Lazarian, A., & Pogosyan, D. 2015, ApJ, 814, 77Falceta-Gon¸calves, D., Lazarian, A., & Kowal, G. 2008, ApJ, 679,537Federrath, C., & Banerjee, S. 2015, MNRAS, 448, 3297Federrath, C., & Klessen, R. S. 2013, ApJ, 763, 51Federrath, C., Klessen, R. S., & Schmidt, W. 2008, ApJ, 688, L79—. 2009, ApJ, 692, 364Fissel, L. M., Ade, P. A. R., Angil`e, F. E., et al. 2016, ApJ, 824,134Galitzki, N., Ade, P. A. R., Angil`e, F. E., et al. 2014, inGround-based and Airborne Telescopes V, ed. L. M. Stepp,R. Gilmozzi, & H. J. Hall, Vol. 9145, International Society forOptics and Photonics (SPIE), 257 – 267Gandilo, N. N., Ade, P. A. R., Angil`e, F. E., et al. 2016, ApJ,824, 84Girichidis, P., Konstandin, L., Whitworth, A. P., & Klessen, R. S.2014, ApJ, 781, 91Girichidis, P., Offner, S. S. R., Kritsuk, A. G., et al. 2020,Space Sci. Rev., 216, 68G´omez, G. C., V´azquez-Semadeni, E., & Zamora-Avil´es, M. 2018,MNRAS, 480, 2939Guszejnov, D., Hopkins, P. F., & Grudi´c, M. Y. 2018, MNRAS,477, 5139Heiles, C., Goodman, A. A., McKee, C. F., & Zweibel, E. G.1993, in Protostars and Planets III, ed. E. H. Levy & J. I.Lunine, 279Hennebelle, P., & Chabrier, G. 2013, ApJ, 770, 150Hennebelle, P., & Falgarone, E. 2012, A&A Rev., 20, 55Hill, T., Motte, F., Didelon, P., et al. 2011, A&A, 533, A94Hoang, T., Cho, J., & Lazarian, A. 2018, ApJ, 852, 129Hoang, T., & Lazarian, A. 2008, MNRAS, 388, 117Hull, C. L. H., & Plambeck, R. L. 2015, Journal of AstronomicalInstrumentation, 4, 1550005Hull, C. L. H., Plambeck, R. L., Bolatto, A. D., et al. 2013, ApJ,768, 159Hull, C. L. H., Plambeck, R. L., Kwon, W., et al. 2014, ApJS,213, 13 Hull, C. L. H., Mocz, P., Burkhart, B., et al. 2017, ApJ, 842, L9Jow, D. L., Hill, R., Scott, D., et al. 2018, MNRAS, 474, 1018Kainulainen, J., & Federrath, C. 2017, A&A, 608, L3Kandel, D., Lazarian, A., & Pogosyan, D. 2016, MNRAS, 461,1227Kirk, J. M., Ward-Thompson, D., Palmeirim, P., et al. 2013,MNRAS, 432, 1424Kowal, G., & Lazarian, A. 2007, ApJ, 666, L69Kowal, G., Lazarian, A., & Beresnyak, A. 2007, ApJ, 658, 423Kowal, G., Lazarian, A., Vishniac, E. T., & Otmianowska-Mazur,K. 2009, ApJ, 700, 63Krumholz, M. R., Klein, R. I., & McKee, C. F. 2007, ApJ, 656,959Landau, L. D., & Lifshitz, E. M. 1959, Fluid mechanicsLarson, R. B. 1981, MNRAS, 194, 809Lazarian, A., Andersson, B. G., & Hoang, T. 2015, Grainalignment: Role of radiative torques and paramagneticrelaxation, 81Lazarian, A., & Pogosyan, D. 2008, ApJ, 686, 350Le˜ao, M. R. M., de Gouveia Dal Pino, E. M., Falceta-Gon¸calves,D., Melioli, C., & Geraissate, F. G. 2009, MNRAS, 394, 157Li, P. S., Myers, A., & McKee, C. F. 2012, ApJ, 760, 33Londrillo, P., & Del Zanna, L. 2000, ApJ, 530, 508McKee, C. F., & Ostriker, E. C. 2007, Annual Review ofAstronomy and Astrophysics, 45, 565Melioli, C., de Gouveia Dal Pino, E. M., de La Reza, R., & Raga,A. 2006, MNRAS, 373, 811Mocz, P., & Burkhart, B. 2019, ApJ, 884, L35Mocz, P., Burkhart, B., Hernquist, L., McKee, C. F., & Springel,V. 2017, ApJ, 838, 40Myers, P. C. 2015, ApJ, 806, 226Padoan, P., Bally, J., Billawala, Y., Juvela, M., & Nordlund, ˚A.1999, ApJ, 525, 318Padoan, P., & Nordlund, ˚A. 2011, ApJ, 730, 40Padoan, P., Nordlund, A., & Jones, B. J. T. 1997, in AmericanInstitute of Physics Conference Series, Vol. 393, AmericanInstitute of Physics Conference Series, ed. S. S. Holt & L. G.Mundy, 101–104Padoan, P., Zweibel, E., & Nordlund, ˚A. 2000, ApJ, 540, 332Palmeirim, P., Andr´e, P., Kirk, J., et al. 2013, A&A, 550, A38Pattle, K., & Fissel, L. 2019, Frontiers in Astronomy and SpaceSciences, 6, 15Planck Collaboration, Ade, P. A. R., & Aghanim. 2016, A&A,586, A138Poidevin, F., Falceta-Gon¸calves, D., Kowal, G., Pino, E. D. G. D.,& Magalh˜aes, A.-M. 2015, Highlights of Astronomy, 16, 391Prato, L., Rice, E. L., & Dame, T. M. 2008, Where are all theYoung Stars in Aquila?, ed. B. Reipurth, Vol. 4, 18Santos-Lima, R., Lazarian, A., de Gouveia Dal Pino, E. M., &Cho, J. 2010, ApJ, 714, 442Scalo, J., & Elmegreen, B. G. 2004, ARA&A, 42, 275Scalo, J., V´azquez-Semadeni, E., Chappell, D., & Passot, T. 1998,ApJ, 504, 835Schneider, N., Andr´e, P., K¨onyves, V., et al. 2013, ApJ, 766, L17Seifried, D., Walch, S., Weis, M., et al. 2020, arXiv e-prints,arXiv:2003.00017Soler, J. D. 2019, A&A, 629, A96Soler, J. D., Ade, P. A. R., Angil`e, F. E., et al. 2017, A&A, 603,A64Soler, J. D., & Hennebelle, P. 2017, A&A, 607, A2Soler, J. D., Hennebelle, P., Martin, P. G., et al. 2013, TheAstrophysical Journal, 774, 128Stephens, I. W., Looney, L. W., Kwon, W., et al. 2013, ApJ, 769,L15—. 2014, Nature, 514, 597Vazquez-Semadeni, E. 1994, ApJ, 423, 681V´azquez-Semadeni, E., Banerjee, R., G´omez, G. C., et al. 2011,MNRAS, 414, 2511V´azquez-Semadeni, E., & Garc´ıa, N. 2001, ApJ, 557, 727
APPENDIX
APPENDIX A: RESOLUTION EFFECTSFig. 19.—
Column density maps with LIC method applied to B ⊥ for sub-Alfv´enic models with M A = 0 . M s = 7 . , for three different LOS. Top row shows the column density map at the initial snapshot ( t = 0 . t ff ) (withfully developed turbulence and no self-gravity), for each LOS. Bottom panel shows the final snapshot (when self-gravity becomes important)for the same LOS. From left to right, the column density distribution is integrated along X (the direction perpendicular to the initial field),45 ◦ with regard to the initial field, and Z (parallel to the initial field) directions. Figures 19 and 20 depict column density maps and PRS time evolution, respectively, for the same sub-Alfv´enicmodels with M A = 0 . M s = 7 . .A quick comparison with the higher resolution counterparts (in Figures 4 and 5), shows that the column density mapsare quite similar. The PRS in Figure 20 shows some differences relative to the higher resolution model in Figure13, particularly at later times when gravity has caused the collapse of several structures and the differences due toresolution become more obvious at these smaller scales, but the general behavior discussed in Sections 2.2 and 3 isalready present at these smaller resolution models.agnetic Field Orientation in GMCs 25 Fig. 20.—
PRS time evolution for the same sub-Alfv´enic model presented in Figure 19. From left to right the PRS is applied to the LOSalong X (the direction perpendicular to the initial field), 45 ◦◦