An observational method to measure the relative fractions of solenoidal and compressible modes in interstellar clouds

Brunt, C. M.; Federrath, C.

doi:10.1093/mnras/stu888

Abstract

We introduce a new method for observationally estimating the fraction of momentum density (ρv) power contained in solenoidal modes (for which ∇ · ρv = 0) in molecular clouds. The method is successfully tested with numerical simulations of supersonic turbulence that produce the full range of possible solenoidal/compressible fractions. At present, the method assumes statistical isotropy, and does not account for anisotropies caused by (e.g.) magnetic fields. We also introduce a framework for statistically describing density–velocity correlations in turbulent clouds.

turbulence, methods: statistical, ISM: clouds, ISM: kinematics and dynamics

1 INTRODUCTION

As the principal sites of star formation in the local universe, molecular clouds demand much observational and theoretical attention. Their structure is extremely complex, driven by the interaction of supersonic turbulence, gravity, and magnetic fields (e.g. Elmegreen & Scalo 2004; Mac Low & Klessen 2004; McKee & Ostriker 2007; Chapman et al. 2011; Heyer & Brunt 2012). Theoretical descriptions of molecular clouds must begin with the structure of the relevant physical fields (density, velocity, etc.) in three dimensions (3D), while our practical information is necessarily restricted to what can be measured from the projection of these fields on to the observational axes – two spatial and, for spectral line data, one (line-of-sight) velocity component. It is a crucially important, yet challenging, problem to relate the projected fields to intrinsic properties of the 3D physical fields.

Only limited information about the 3D density field, ρ(x, y, z), can be derived via analysis of the projected 2D column density field, N(x, y), which is observationally obtained by extinction measurements, dust emission, or integrated spectral line intensities. For the latter, we may be aided by (or perhaps limited by) density-selectivity of particular molecular transitions, yet molecular line emission is the only source of information for studying the dynamics of molecular clouds.

Arguably the most useful property of N(x, y) is that its Fourier transform, |$\tilde{N}(k_{x}, k_{y})$|⁠, is directly proportional to a 2D slice through the Fourier transform of the density field, |$\tilde{\rho }(k_{x}, k_{y}, k_{z})$|⁠, where the line-of-sight wavevector, k_z = 0; i.e. |$\tilde{N}(k_{x}, k_{y}) \propto \tilde{\rho }(k_{x}, k_{y}, k_{z} = 0)$|⁠. This allows, under the assumption of isotropy, the density power spectrum to be derived (e.g. Stutzki et al. 1998), the 3D density variance to be inferred (Fischera & Dopita 2004; Brunt, Federrath & Price 2010a, hereafter BFP), and an estimate of the 3D density probability distribution function (PDF) to be constructed (Brunt, Federrath & Price 2010b). The essence of the BFP method (relating the 2D normalized column density variance to the 3D normalized density variance) lies in determining the fraction of variance contained in a single 2D slice of the 3D power spectrum.

While BFP focused primarily on a method to relate density and column density statistics, they also presented a brief outline of the same method applied to velocity fields, though noting that the natural observational density-weighting of the velocity field would potentially cause problems. In this paper, we present an extended development of the outline method presented in BFP, subject to two key modifications. First, the density-weighting of the velocity field means that the physical field suitable for BFP-like analysis is the ‘momentum density’ field, |${\rho }{\boldsymbol {v}}$|⁠, rather than the velocity field alone. Secondly, the realization that only transverse (solenoidal) modes (for which |$\nabla \cdot {\rho }{\boldsymbol {v}} = 0$|⁠) are projected into 2D allows us to extend the BFP method to estimate the fraction of momentum density power that is held in solenoidal modes, if given an estimate of the total momentum density power (through a spectral line imaging observation). The ideas underpinning the extended method are presented below, along with a demonstration of its applicability using numerical simulations of turbulent clouds. The method as presented assumes statistical isotropy, so should not be applied to clouds for which significant anisotropy is observed or suspected due to (e.g.) the presence of a strong magnetic field at low Mach numbers (BFP).

The layout of this paper is as follows. In Section 2, we present the method for observationally estimating the fraction of momentum power in transverse modes, after first introducing the background considerations necessary for its formulation. In Section 3, we describe the numerical simulations to be used for testing purposes. In Section 4, we test the method by applying it to the numerical simulations, followed by a discussion (Section 5) and summary (Section 6). In the appendix, we examine statistical aspects of density–velocity correlation, which has a small effect on the method.

2 SOLENOIDAL AND COMPRESSIBLE MODES

In this section, we discuss geometrical properties of divergence-free (‘transverse’ or ‘solenoidal’) and curl-free (‘longitudinal’ or ‘compressible’) modes in a general vector field (Section 2.1), before considering the consequences of projection of such a field from 3D to 2D (Section 2.2). In the following, we will use transverse/solenoidal and longitudinal/compressible to refer, respectively, to divergence-free and curl-free components. Subsequently, we identify the momentum field as the most relevant physical field of interest for quantitative analysis and develop a method by which the fraction of power in transverse momentum modes may be estimated observationally (Section 2.3).

2.1 General considerations

The Helmholtz Decomposition Theorem (Helmholtz 1858) states that, as a function of position |$\boldsymbol {x} = (x, y, z)$|⁠, an arbitrary 3D vector field, |$\boldsymbol {F}(\boldsymbol {x})$|⁠, can be represented as the sum of a purely transverse field, |$\boldsymbol {F}_{\perp }(\boldsymbol {x})$|⁠, and a purely longitudinal field, |$\boldsymbol {F}_{||}(\boldsymbol {x})$|⁠:

\begin{equation} \boldsymbol {F}(\boldsymbol {x}) = \boldsymbol {F}_{\perp }(\boldsymbol {x}) + \boldsymbol {F}_{||}(\boldsymbol {x}), \end{equation}

(1)

where

\begin{equation} \boldsymbol {\nabla } \cdot \boldsymbol {F}_{\perp } = 0, \end{equation}

(2)

\begin{equation} \boldsymbol {\nabla } \times \boldsymbol {F}_{||} = 0. \end{equation}

(3)

In Fourier space, the equivalent relations to equations (2) and (3) are (using a tilde to represent the Fourier-transformed fields)

\begin{equation} \tilde{\boldsymbol {F}}(\boldsymbol {k}) = \tilde{\boldsymbol {F}}_{\perp }(\boldsymbol {k}) + \tilde{\boldsymbol {F}}_{||}(\boldsymbol {k}), \end{equation}

(4)

\begin{equation} \boldsymbol {k} \cdot \tilde{\boldsymbol {F}}_{\perp } = 0, \end{equation}

(5)

\begin{equation} \boldsymbol {k} \times \tilde{\boldsymbol {F}}_{||} = 0, \end{equation}

(6)

where |$\boldsymbol {k} = (k_{x}, k_{y}, k_{z})$| are wavevectors, and the Fourier-transformed field is defined, over a cubical spatial region of length L on each axis, by

\begin{equation} \tilde{\boldsymbol {F}}({\boldsymbol {k}}) = \displaystyle \int _{-L/2}^{L/2}\displaystyle \int _{-L/2}^{L/2}\displaystyle \int _{-L/2}^{L/2} \; {\rm{d}}^{3}{\boldsymbol {x}} \; {\boldsymbol {F}}({\boldsymbol {x}}) \; {\rm{e}}^{ -2{\pi } {\rm{i}} {\boldsymbol {k}} \cdot {\boldsymbol {x}} / L}. \end{equation}

(7)

In the following, we will assume that a reference frame can be chosen in which |$\langle {\boldsymbol {F}} \rangle = 0$|⁠, so that the spatial average, |$\langle {\boldsymbol {F}}^{2} \rangle$|⁠, is equal to the field variance, |$\sigma ^{2}_{\boldsymbol {F}}$|⁠. It is possible, making use of Parseval's Theorem, to calculate the field variances via Fourier space using

\begin{equation} \sigma _{{\boldsymbol {F}}}^{2} = {\frac{1}{L^{6}}} \displaystyle \sum _{k_{x}=-\infty }^{\infty } \sum _{k_{y}=-\infty }^{\infty } \sum _{k_{z}=-\infty }^{\infty } \; \tilde{{\boldsymbol {F}}} \cdot \tilde{{\boldsymbol {F}}}^{*}, \end{equation}

(8)

\begin{equation} \sigma _{{\boldsymbol {F}}_{\perp }}^{2} = {\frac{1}{L^{6}}} \displaystyle \sum _{k_{x}=-\infty }^{\infty } \sum _{k_{y}=-\infty }^{\infty } \sum _{k_{z}=-\infty }^{\infty } \; \tilde{{\boldsymbol {F}}}_{\perp } \cdot \tilde{{\boldsymbol {F}}}_{\perp }^{*}, \end{equation}

(9)

\begin{equation} \sigma _{{\boldsymbol {F}}_{||}}^{2} = {\frac{1}{L^{6}}} \displaystyle \sum _{k_{x}=-\infty }^{\infty } \sum _{k_{y}=-\infty }^{\infty } \sum _{k_{z}=-\infty }^{\infty } \; \tilde{{\boldsymbol {F}}}_{||} \cdot \tilde{{\boldsymbol {F}}}_{||}^{*}. \end{equation}

(10)

Note that the local orthogonality of |$\tilde{\boldsymbol {F}}_{\perp }$| and |$\tilde{\boldsymbol {F}}_{||}$| (ensuring that |$\sigma ^{2}_{{\boldsymbol {F}}} = \sigma ^{2}_{{\boldsymbol {F}}_{\perp }} + \sigma ^{2}_{{\boldsymbol {F}}_{||}}$|⁠) does not mean that the direct space fields |${\boldsymbol {F}}_{\perp }$| and |${\boldsymbol {F}}_{||}$| are also locally orthogonal. That is, at any single field point, we have

\begin{equation} F^{2} = {\boldsymbol {F}} \cdot {\boldsymbol {F}} = F^{2}_{\perp } + F^{2}_{||} + 2 {\boldsymbol {F}}_{\perp } \cdot {\boldsymbol {F}}_{||} \ne F^{2}_{\perp } + F^{2}_{||}, \end{equation}

(11)

since in general |${\boldsymbol {F}}_{\perp } \cdot {\boldsymbol {F}}_{||} \ne 0$|locally. However, this dot product vanishes when averaged over the entire space containing |${\boldsymbol {F}}$|⁠. (This is required by Parseval's Theorem, expressed in equations (8)–(10). Note that |$\tilde{{\boldsymbol {F}}} \cdot \tilde{{\boldsymbol {F}}}^{*} = \tilde{{\boldsymbol {F}}}_{\perp } \cdot \tilde{{\boldsymbol {F}}}_{\perp }^{*} + \tilde{{\boldsymbol {F}}}_{||} \cdot \tilde{{\boldsymbol {F}}}_{||}^{*}$| is satisfied at each point in Fourier space, so that |$\sigma _{{\boldsymbol {F}}}^{2} = \langle {\boldsymbol {F}} \cdot {\boldsymbol {F}} \rangle = \sigma _{{\boldsymbol {F}}_{\perp }}^{2} + \sigma _{{\boldsymbol {F}}_{||}}^{2} = \langle F^{2}_{\perp } \rangle + \langle F^{2}_{||} \rangle$|⁠, provided the averages are computed in the |$\langle {\boldsymbol {F}} \rangle = 0$| frame, as we have assumed.)

The Helmholtz Decomposition (equation 1) is unique, up to a vector constant, provided that the field |$\boldsymbol {F}$| falls to zero on its outer boundary. More generally, |$\boldsymbol {F}$| could in principle contain a contribution from a curl-free, divergence-free component, |$\boldsymbol {F}_{L}$|⁠, given by the gradient of a scalar harmonic field ϕ_L that satisfies the Laplace equation (∇²ϕ_L = 0). Specifically, the curl-free nature of |$\boldsymbol {F}_{L}$| follows from its definition as the gradient of a scalar, while its divergence-free nature requires that ϕ_L satisfies the Laplace equation. As a harmonic field, ϕ_L must obey the mean value theorem, i.e. that its value at any point |${\boldsymbol {x}}$| is equal to its average on any spherical surface of arbitrary radius surrounding |${\boldsymbol {x}}$|⁠. This means that ϕ_L can contain no local maxima or minima, and therefore that any maxima/minima must occur on its outer boundary. The properties of ϕ_L are therefore determined entirely by boundary conditions: |$\boldsymbol {F}_{L} = {\nabla }\phi _{L}$| will quantify large-scale, smooth gradients in |$\boldsymbol {F}$| that cannot be assigned to either |${\boldsymbol {F}}_{\perp }$| or |${\boldsymbol {F}}_{||}$|⁠.

For our study here, we use numerically simulated fields (with the momentum density |${\rho }{\boldsymbol {v}}$| playing the role of |$\boldsymbol {F}$|⁠) that obey periodic boundary conditions. For these fields, the multiplicity of choices for the ‘boundary’ and the condition of no local maxima/minima ensure that ϕ_L is a constant and therefore that |$\boldsymbol {F}_{L} = 0$|⁠, and the Helmholtz Decomposition is unique. For real momentum fields encountered in the interstellar medium (ISM), we must rely on finding clouds that are sufficiently isolated that boundary conditions on |${\rho }{\boldsymbol {v}}$| are not a concern. Our method is therefore best-suited to molecular clouds that are bounded in space (i.e. by a fall-off in density ρ that ensures no non-zero |${\rho }{\boldsymbol {v}}$| values on their boundary).

The more widely distributed atomic medium is less suited to application of our method, since it will be difficult to ensure the absence of large-scale gradients in any finite field. Finally, it should be noted that |${\rho }{\boldsymbol {v}}$| will be continuous across the atomic/molecular transition and restriction to the molecular component is a necessarily limited description of the ISM as a single fluid – not to mention restriction to practically observable regions using trace molecules, and the possible influence of inter-mixed atomic/molecular zones. However, these problems are not in principle insurmountable (with sufficient data) though they may pose considerable challenges if a complete description of the ISM fluid is desired.

2.2 Projection from 3D to 2D

Consider the case where we have access to only one scalar component of the vector field |${\boldsymbol {F}}$|⁠. We assume that the accessible component is oriented along the line of sight (as in, for example, a spectral line observation), which we take as the z-direction. The observable component of |${\boldsymbol {F}}$| is then F_z = F_z⊥ + F_z|| where F_z⊥ and F_z|| are the z components of the transverse and longitudinal parts of |${\boldsymbol {F}}$|⁠. The contribution of F_z to the Fourier-transformed field – |$\tilde{\boldsymbol {F}}(\boldsymbol {k})$| – is |$\tilde{F}_{z}$| and this is oriented along the k_z-direction (i.e. |$F_{z} \hat{\boldsymbol {z}}$| transforms into |$\tilde{F}_{z} \hat{\boldsymbol {k}}_{z}$| where |$\hat{\boldsymbol {z}}$| and |$\hat{\boldsymbol {k}}_{z}$| are unit vectors in the z and k_z directions, respectively).

The condition |$\boldsymbol {k} \cdot \tilde{\boldsymbol {F}}_{\perp } = 0$| requires that |$\tilde{F}_{z\perp } = 0$| along the k_z axis (where k_x = k_y = 0). Clearly we must also find that |$\tilde{F}_{z\perp } = \tilde{F}_{z}$| everywhere in the plane k_z = 0, since there the condition |$\boldsymbol {k} \times \tilde{\boldsymbol {F}}_{||} = 0$| requires that |$\tilde{F}_{z||} = 0$| in this plane.

In this paper, we only consider fields for which the transverse and longitudinal fields, |${\boldsymbol {F}}_{\perp }$| and |${\boldsymbol {F}}_{||}$|⁠, are statistically isotropic, i.e. that their power spectra (⁠|$\tilde{{\boldsymbol {F}}}_{\perp } \cdot \tilde{{\boldsymbol {F}}}_{\perp }^{*}$| and |$\tilde{{\boldsymbol {F}}}_{||} \cdot \tilde{{\boldsymbol {F}}}_{||}^{*}$|⁠) may be written as functions of |$k=|\boldsymbol {k}|$| alone, with no explicit angular dependence in Fourier space. In this case, we may write

\begin{equation} \tilde{{\boldsymbol {F}}}_{\perp } \cdot \tilde{{\boldsymbol {F}}}_{\perp }^{*} = F^{2}_{{\perp }0} f_{\perp }(k), \end{equation}

(12)

\begin{equation} \tilde{{\boldsymbol {F}}}_{||} \cdot \tilde{{\boldsymbol {F}}}_{||}^{*} = F^{2}_{||0} f_{||}(k), \end{equation}

(13)

where |$F^{2}_{{\perp }0}$| and |$F^{2}_{||0}$| are scaling factors and f_⊥(k) and f_||(k) describe the k-dependent power distributions. Note that for such isotropic fields the transverse and longitudinal power distributions for a single scalar component are not isotropic, but have predictable anisotropic structure determined by the following equations:

\begin{equation} \tilde{F}_{z||} \tilde{F}_{z||}^{*} = \tilde{{\boldsymbol {F}}}_{||} \cdot \tilde{{\boldsymbol {F}}}_{||}^{*} \frac{k^{2}_{z}}{k^{2}}, \end{equation}

(14)

\begin{equation} \tilde{F}_{z\perp } \tilde{F}_{z\perp }^{*} = \tilde{{\boldsymbol {F}}}_{\perp } \cdot \tilde{{\boldsymbol {F}}}_{\perp }^{*} \frac{k^{2}_{x}+k^{2}_{y}}{2k^{2}}, \end{equation}

(15)

\begin{equation} \tilde{F}_{x||} \tilde{F}_{x||}^{*} = \tilde{{\boldsymbol {F}}}_{||} \cdot \tilde{{\boldsymbol {F}}}_{||}^{*} \frac{k^{2}_{x}}{k^{2}}, \end{equation}

(16)

\begin{equation} \tilde{F}_{x\perp } \tilde{F}_{x\perp }^{*} = \tilde{{\boldsymbol {F}}}_{\perp } \cdot \tilde{{\boldsymbol {F}}}_{\perp }^{*} \frac{k^{2}_{y}+k^{2}_{z}}{2k^{2}}, \end{equation}

(17)

\begin{equation} \tilde{F}_{y||} \tilde{F}_{y||}^{*} = \tilde{{\boldsymbol {F}}}_{||} \cdot \tilde{{\boldsymbol {F}}}_{||}^{*} \frac{k^{2}_{y}}{k^{2}}, \end{equation}

(18)

\begin{equation} \tilde{F}_{y\perp } \tilde{F}_{y\perp }^{*} = \tilde{{\boldsymbol {F}}}_{\perp } \cdot \tilde{{\boldsymbol {F}}}_{\perp }^{*} \frac{k^{2}_{x}+k^{2}_{z}}{2k^{2}}. \end{equation}

(19)

It is worth briefly mentioning here that, as can be checked in the above equations, the individual scalar components of the vector field are only fully isotropic if, at each |${\boldsymbol {k}}$|⁠, we have |$\tilde{{\boldsymbol {F}}}_{\perp }$| · |$\tilde{{\boldsymbol {F}}}_{\perp }^{*}$| = |$2\tilde{{\boldsymbol {F}}}_{||} \cdot \tilde{{\boldsymbol {F}}}_{||}^{*}$|⁠. This is the case for (e.g.) fractional Brownian motion fields – i.e. there is twice as much power in the transverse component as the longitudinal component. In general, this will not be true.

To illustrate the above power spectra, we show example 3D renderings of |$\tilde{F}_{z\perp } \tilde{F}_{z\perp }^{*}$| and |$\tilde{F}_{z_{||}} \tilde{F}_{z_{||}}^{*}$| in Fig. 1. Isotropic power spectra of |$\tilde{{\boldsymbol {F}}}_{\perp } \cdot \tilde{{\boldsymbol {F}}}_{\perp }^{*} \propto k^{-4}$|⁠, and |$\tilde{{\boldsymbol {F}}}_{||} \cdot \tilde{{\boldsymbol {F}}}_{||}^{*} \propto k^{-4}$| have been assumed. The power spectrum of the longitudinal component of F_z has a characteristic ‘hourglass’ appearance, resulting from the suppression (nulling) of power near (at) k_z = 0. Since |${\tilde{\boldsymbol {F}}}_{||}$| is aligned with |${\boldsymbol {k}}$|⁠, this means that |$\tilde{F}_{z_{||}}$| must be zero in the k_z = 0 plane. Conversely, the transverse power is maximized (at fixed k) in the plane k_z = 0, and diminishes as |k_z/k| approaches unity. An instructive reference point can be obtained by considering the power distributions along the line (k_x = 0, k_y = 0, k_z ≠ 0). Here, the longitudinal power is entirely contained in the z-component, while the transverse power is equally split (assuming isotropy) between the x-component and the y-component (equations 17 and 19 with k_x = k_y = 0).

$3D renderings of constant power surfaces in the power spectra of Fz⊥ (left) and $F_{z_{||}}$ (right), given by equations (15) and (14), respectively. Isotropic power spectra: $\tilde{{\boldsymbol {F}}}_{\perp } \cdot \tilde{{\boldsymbol {F}}}_{\perp }^{*} \propto k^{-4}$, and $\tilde{{\boldsymbol {F}}}_{||} \cdot \tilde{{\boldsymbol {F}}}_{||}^{*} \propto k^{-4}$ have been assumed. In both panels, the origin of coordinates ($\boldsymbol {k} = 0$) lies in the centre of the image.$

Figure 1.

3D renderings of constant power surfaces in the power spectra of F_z⊥ (left) and |$F_{z_{||}}$| (right), given by equations (15) and (14), respectively. Isotropic power spectra: |$\tilde{{\boldsymbol {F}}}_{\perp } \cdot \tilde{{\boldsymbol {F}}}_{\perp }^{*} \propto k^{-4}$|⁠, and |$\tilde{{\boldsymbol {F}}}_{||} \cdot \tilde{{\boldsymbol {F}}}_{||}^{*} \propto k^{-4}$| have been assumed. In both panels, the origin of coordinates (⁠|$\boldsymbol {k} = 0$|⁠) lies in the centre of the image.

Open in new tab Download slide

The power spectra of the other (x, y) components share the same form, but with different orientations, such that (as we assume here) the total transverse and longitudinal power spectra are isotropic (functions of k alone):

\begin{equation} \tilde{{\boldsymbol {F}}}_{\perp } \cdot \tilde{{\boldsymbol {F}}}_{\perp }^{*} = \tilde{F}_{x\perp } \tilde{F}_{x\perp }^{*} + \tilde{F}_{y\perp } \tilde{F}_{y\perp }^{*} + \tilde{F}_{z\perp } \tilde{F}_{z\perp }^{*} = F^{2}_{{\perp }0}f_{\perp }(k), \end{equation}

(20)

\begin{equation} \tilde{{\boldsymbol {F}}}_{||} \cdot \tilde{{\boldsymbol {F}}}_{||}^{*} = \tilde{F}_{x||} \tilde{F}_{x||}^{*} + \tilde{F}_{y||} \tilde{F}_{y||}^{*} + \tilde{F}_{z||} \tilde{F}_{z||}^{*} = F^{2}_{{||}0}f_{||}(k), \end{equation}

(21)

as can be checked by equations (12)–(19). In practice, power spectra, whether observed or simulated, will not be precisely isotropic but we assume that they are statistically isotropic – i.e. that the power values fluctuate randomly around the form assumed above. Note that if the isotropy is only statistical, the geometry still requires that |$\tilde{F}_{z}$| has no longitudinal power in the k_z = 0 plane.

An important property of the above power distributions is that in the plane k_z = 0:

\begin{equation} \tilde{F}_{z} \tilde{F}_{z}^{*} = \tilde{F}_{z\perp } \tilde{F}_{z\perp }^{*} =\frac{1}{2} \tilde{{\boldsymbol {F}}}_{\perp } \cdot \tilde{{\boldsymbol {F}}}_{\perp }^{*}. \end{equation}

(22)

It is straightforward to show that if F_z is spatially averaged over the z-axis via

\begin{equation} F_{z,p}(x,y) = \frac{1}{L} \displaystyle \int _{-L/2}^{L/2} {\rm{d}}z \; F_{z}(x, y, z), \end{equation}

(23)

then the Fourier transform of the line-of-sight-averaged field, F_{z, p}, is given by

\begin{equation} \tilde{F}_{z,p}(k_{x},k_{y}) = \frac{1}{L} \tilde{F}_{z}(k_{x},k_{y},k_{z}=0) \end{equation}

(24)

(see BFP). This means that the projected field, F_{z, p}, contains only contributions from transverse structure in 3D, since its Fourier transform is directly proportional to a k_z = 0 cut through the 3D Fourier transform of F_z (which contains only transverse contributions). Therefore, its power spectrum is a direct measure of the k_z = 0 plane of the 3D power spectrum of F_z⊥ (and, indeed, of |${\boldsymbol {F}}_{\perp }$|⁠). The power spectrum, as a representation of the variance of F_{z, p}, can be used to obtain an estimate of the variance of F_z⊥ in 3D, if the power spectra conform to the structures given in equations (14)–(19). Explicitly, the power spectrum of F_{z, p} is, for k_z = 0:

\begin{eqnarray} \tilde{F}_{z,p} \tilde{F}_{z,p}^{*} (k_{x}, k_{y}) & =& \frac{1}{L^{2}}\tilde{F}_{z\perp } \tilde{F}_{z\perp }^{*} (k_{x}, k_{y}, k_{z}=0) \nonumber \\ & =& \frac{1}{L^{2}} \frac{\tilde{{\boldsymbol {F}}}_{\perp } \cdot \tilde{{\boldsymbol {F}}}_{\perp }^{*}}{2} (k_{x}, k_{y}, k_{z}=0). \end{eqnarray}

(25)

Using Parseval's Theorem, and referring to equation (15), the variance of F_{z, p} is given by

\begin{equation} \sigma _{F_{z,p}}^{2}= {\frac{1}{L^{4}}} \displaystyle \sum _{k_{x}=-\infty }^{\infty } \sum _{k_{y}=-\infty }^{\infty } \; \tilde{F}_{z,p} \tilde{F}_{z,p}^{*}, \end{equation}

(26)

or making use of equation (25):

\begin{equation} \sigma _{F_{z,p}}^{2}= {\frac{1}{L^{6}}} \displaystyle \sum _{k_{x}=-\infty }^{\infty } \sum _{k_{y}=-\infty }^{\infty } \; \frac{\tilde{{\boldsymbol {F}}}_{\perp } \cdot \tilde{{\boldsymbol {F}}}_{\perp }^{*}}{2}. \end{equation}

(27)

The variance of F_z⊥ in 3D is given by

\begin{equation} \sigma _{F_{z\perp }}^{2} = {\frac{1}{L^{6}}} \displaystyle \sum _{k_{x}=-\infty }^{\infty } \sum _{k_{y}=-\infty }^{\infty } \sum _{k_{z}=-\infty }^{\infty } \; \tilde{{\boldsymbol {F}}}_{\perp } \cdot \tilde{{\boldsymbol {F}}}_{\perp }^{*} \frac{k^{2}_{x}+k^{2}_{y}}{2k^{2}}, \end{equation}

(28)

so that, under the assumption of isotropy (i.e. equation 12), we may then use the measured variance of F_{z, p} to estimate the variance of F_z⊥ in 3D via

\begin{equation} \sigma _{F_{z\perp }}^{2} = \sigma _{F_{z,p}}^{2} \times \frac{\displaystyle \sum _{k_{x}=-\infty }^{\infty } \sum _{k_{y}=-\infty }^{\infty } \sum _{k_{z}=-\infty }^{\infty } \; f_{\perp }(k) \frac{k^{2}_{x}+k^{2}_{y}}{k^{2}}}{\displaystyle \sum _{k_{x}=-\infty }^{\infty } \sum _{k_{y}=-\infty }^{\infty } \; f_{\perp }(k)}, \end{equation}

(29)

where we note that f_⊥(k) can be directly measured using the power spectrum of the projected field, |$\tilde{F}_{z,p} \tilde{F}_{z,p}^{*}$|⁠. The overall scaling of the fields (controlled by F_⊥0 in equation 12) is unimportant in determining the ratio|$\sigma _{F_{z\perp }}^{2}/\sigma _{F_{z,p}}^{2}$| but is important if the absolute variance |$\sigma _{F_{z\perp }}^{2}$| is desired. Even more straightforwardly, noting that for an isotropic field, |$\sigma _{F_{z\perp }}^{2} = \sigma _{{\boldsymbol {F}}_{\perp }}^{2}\,\!/ 3$|⁠, we can also write

\begin{equation} \sigma _{F_{z\perp }}^{2} = \frac{2}{3} \sigma _{F_{z,p}}^{2} \times \frac{\displaystyle \sum _{k_{x}=-\infty }^{\infty } \sum _{k_{y}=-\infty }^{\infty } \sum _{k_{z}=-\infty }^{\infty } \; f_{\perp }(k) }{\displaystyle \sum _{k_{x}=-\infty }^{\infty } \sum _{k_{y}=-\infty }^{\infty } \; f_{\perp }(k)}. \end{equation}

(30)

If, in addition, a measurement of the total variance of F_z is available, then the fractional power in transverse modes can be calculated. Even if only the ratio of projected-to-total variance, |$\sigma _{F_{z,p}}^{2}\,\!/\!\sigma _{F_{z}}^{2}$|⁠, is known, the fractional power in transverse modes can still be calculated via

\begin{equation} \frac{\sigma _{F_{z\perp }}^{2}}{\sigma _{F_{z}}^{2}} = \frac{\sigma _{F_{z,p}}^{2}}{{\sigma _{F_{z}}^{2}}} \times \frac{\displaystyle \sum _{k_{x}=-\infty }^{\infty } \sum _{k_{y}=-\infty }^{\infty } \sum _{k_{z}=-\infty }^{\infty } \; f_{\perp }(k) \frac{k^{2}_{x}+k^{2}_{y}}{k^{2}}}{\displaystyle \sum _{k_{x}=-\infty }^{\infty } \sum _{k_{y}=-\infty }^{\infty } \; f_{\perp }(k)}, \end{equation}

(31)

or

\begin{equation} \frac{\sigma _{F_{z\perp }}^{2}}{\sigma _{F_{z}}^{2}} = \frac{2}{3} \frac{\sigma _{F_{z,p}}^{2}}{{\sigma _{F_{z}}^{2}}} \times \frac{\displaystyle \sum _{k_{x}=-\infty }^{\infty } \sum _{k_{y}=-\infty }^{\infty } \sum _{k_{z}=-\infty }^{\infty } \; f_{\perp }(k) }{\displaystyle \sum _{k_{x}=-\infty }^{\infty } \sum _{k_{y}=-\infty }^{\infty } \; f_{\perp }(k)}, \end{equation}

(32)

and again under the assumption of isotropy, this should be equal to the fractional power in transverse modes for the full vector field |${\boldsymbol {F}}$|⁠, i.e.

\begin{equation} \frac{\sigma _{F_{\perp }}^{2}}{\sigma _{F}^{2}} \approx \frac{\sigma _{F_{z\perp }}^{2}}{\sigma _{F_{z}}^{2}}. \end{equation}

(33)

Note that there is no requirement that the longitudinal power spectrum is to be known, nor be in any way dependent on the transverse power spectrum.

A simpler version of equation (29) was used by BFP to estimate the variance of normalized density, |$\sigma ^{2}_{\rho /\rho _{0}}$| from the observationally accessible normalized column density variance |$\sigma ^{2}_{N/N_{0}}$| (where ρ and N are density and column density, and ρ₀ and N₀ are their mean values, respectively). The appropriately modified form of equation (29) for this purpose is

\begin{equation} \sigma ^{2}_{\rho /\rho _{0}} = \sigma ^{2}_{N/N_{0}} \times \frac{\left(\displaystyle \sum _{k_{x}=-\infty }^{\infty } \sum _{k_{y}=-\infty }^{\infty } \sum _{k_{z}=-\infty }^{\infty } \; f(k)\right) - f(0)}{\left(\displaystyle \sum _{k_{x}=-\infty }^{\infty } \sum _{k_{y}=-\infty }^{\infty } \; f(k)\right) - f(0)}, \end{equation}

(34)

where it was assumed that the density power spectrum was isotropic, described by the function f (k), and the mean value of the fields (i.e. the zero-frequency component f(0)) is explicitly subtracted from the summations. Above we had assumed that a frame could be chosen in which |$\langle {\boldsymbol {F}} \rangle = 0$|⁠, whereas BFP used the positive-definite nature of ρ and N to provide a suitable normalization. Obviously, equation (34) does not include the |$(k^{2}_{x}+k^{2}_{y})/k^{2}$| factor either, since ρ is a scalar and not subject to the vector projection factors that determine F_z⊥. For practical application of equations (29)–(34) the sums only extend to the maximum wavevector observable in the field, the consequences of which are discussed in BFP and Brunt (2010).

2.3 Observational considerations

We now review the information provided by spectral line observations of the ISM, with a view to determining which physical field(s) may be analysed by the above system. The principal requirements for such a field are (1) that it is a 3D vector field, and (2) that one of its components can be projected (or averaged) over the line of sight with no weighting by other variable physical fields. It turns out, as explained below, that the momentum density field, |${\boldsymbol {p}} = \rho {\boldsymbol {v}}$|⁠, (hereafter simply ‘momentum’) most closely satisfies these requirements. Application of the method to velocity fields is impossible except under conditions of uniform density, which are essentially never encountered in the ISM. Below we therefore develop a scheme whereby the fraction of momentum power in transverse modes may be estimated observationally – i.e. evaluation of equation (31).

Here, we only consider an optically thin isothermal medium with uniform excitation, so that the infinitesimal contribution to the spectral line intensity generated by density ρ along an infinitesimal path-length dz at position z is given by

\begin{equation} {\rm{d}}I(v) = e \rho \phi (v - v_{z}(z))\, {\rm{d}}z, \end{equation}

(35)

where e is a constant and ϕ(v − v_z) is the normalized profile function, which can usually be represented by a Gaussian:

\begin{equation} \phi (v - v_{z}) = \frac{1}{\sqrt{2 {\pi } \sigma ^{2}_{t,i}}} \exp \left( - \frac{(v - v_{z})^{2}}{2\sigma ^{2}_{t,i}} \right), \end{equation}

(36)

where |$\sigma ^{2}_{t,i}$| represents the dispersion caused by thermal and instrumental broadening. Integrating along the line of sight, the observed spectral line intensity is then

\begin{equation} I(x, y, v) = e \displaystyle \int _{-L/2}^{L/2} {\rm{d}}z \; \rho (x, y, z) \phi (v - v_{z}(x, y, z)), \end{equation}

(37)

where we have assumed the emitting medium lies within a region of spatial size L. For molecular spectral line observations, |$\sigma ^{2}_{t,i}$| is usually very small compared to the overall velocity dispersion, so that a reasonable approximation is ϕ(v − v_z) = δ(v − v_z), where δ(v − v_z) is the Dirac delta function. In this case, equation (37) effectively describes the intensity as a ‘density-weighted histogram’ of line-of-sight velocity, and indeed, this is a standard method to represent optically thin spectral line observations of numerical simulations (e.g. Falgarone et al. 1994; Ostriker, Stone & Gammie 2001). (Note that the delta function may be satisfied multiple times along the line of sight for a turbulent medium.)

With ϕ(v − v_z) = δ(v − v_z), consider the observationally accessible integral (the first velocity-moment of the intensity):

\begin{eqnarray} W_{1}(x, y)&=& \displaystyle \int _{-\infty }^{\infty } {\rm{d}}v \; I(x, y, v) v \nonumber \\ &=& \displaystyle \int _{-\infty }^{\infty } {\rm{d}}v \; e \displaystyle \int _{-L/2}^{L/2} {\rm{d}}z \; \rho (x, y, z) \delta (v - v_{z}) v \nonumber \\ &=& e \displaystyle \int _{-L/2}^{L/2} {\rm{d}}z \; \rho (x, y, z) v_{z}(x, y, z) \nonumber \\ &=& e \displaystyle \int _{-L/2}^{L/2} {\rm{d}}z \; p_{z}(x, y, z) \nonumber \\ &=& e L p_{z,p}, \end{eqnarray}

(38)

where p_z = ρv_z is the z-component of the momentum, and p_{z, p} is (cf. equation 23)

\begin{equation} p_{z,p}(x, y) = \frac{1}{L} \int _{-L/2}^{L/2} {\rm{d}}z \; p_{z}(x, y, z). \end{equation}

(39)

The momentum field therefore satisfies the spatial projection requirement (up to constants of proportionality), while the velocity field does so only under conditions of uniform density. Therefore, in equation (31) we will set F_z = p_z and F_{z, p} = p_{z, p}.

At this point, we identify the following ratio: |$\sigma ^{2}_{p_{z,p}}/\sigma ^{2}_{p_{z}}$| – i.e. the fraction of z-momentum power (variance) projected into 2D – as the most relevant quantity to estimate observationally, since measurements of |$\sigma ^{2}_{p_{z,p}}/\sigma ^{2}_{p_{z}}$| and f_⊥(k) are needed for the evaluation of equation (31). We already have f_⊥(k), as this can be derived from the angular average of the power spectrum of W₁(x, y) (note that the overall normalization is unimportant).

While e and L (in equation 38) can in principle be estimated, a better procedure is to normalize them out. The integrated intensity (the zeroth velocity-moment of the intensity) is

\begin{equation} W_{0}(x,y) = \displaystyle \int _{-\infty }^{\infty } {\rm{d}}v \; I(x, y, v) = e L \rho _{p}(x, y) = eN, \end{equation}

(40)

where

\begin{equation} \rho _{p}(x, y) = \frac{1}{L} \int _{-L/2}^{L/2} {\rm{d}}z \; \rho (x, y, z) \end{equation}

(41)

is the line-of-sight average of ρ (N = Lρ_p is the column density).

By computing spatial averages over x and y (denoted by angle brackets) the ratio:

\begin{equation} \frac{\langle W^{2}_{1} \rangle }{\langle W^{2}_{0} \rangle } = \frac{\sigma ^{2}_{p_{z,p}}}{\langle \rho ^{2}_{p} \rangle } \end{equation}

(42)

can be formed, where we have assumed the calculations are done in the zero-momentum frame (〈p_{z, p}〉 = 0, or equivalently, 〈W₁〉 = 0).

We review now what remains to be calculated in order to form the ratio |$\sigma ^{2}_{p_{z,p}}/\sigma ^{2}_{p_{z}}$| needed for application of equation (31). As is evident from the dimensionality of |$\sigma ^{2}_{p_{z,p}}/\langle \rho ^{2}_{p} \rangle$|⁠, it turns out that a measurement of its 3D analogue, |$\sigma ^{2}_{p_{z}}/\langle \rho ^{2} \rangle$|⁠, is the most useful way to proceed. If this quantity was available, then we could form the ratio:

\begin{equation} \frac{\sigma ^{2}_{p_{z,p}}/\langle \rho _{p}^{2} \rangle }{\sigma ^{2}_{p_{z}}/\langle \rho ^{2} \rangle } = \frac{\sigma ^{2}_{p_{z,p}}}{\sigma ^{2}_{p_{z}}} \frac{\langle \rho ^{2} \rangle }{\langle \rho _{p}^{2} \rangle } = \frac{\sigma ^{2}_{p_{z,p}}}{\sigma ^{2}_{p_{z}}} \frac{\langle (\rho /\rho _{0})^{2} \rangle }{\langle (N/N_{0})^{2} \rangle }, \end{equation}

(43)

where in the last step we have written ρ = ρ₀(ρ/ρ₀) and ρ_p = N/L = ρ₀(N/ρ₀L) = ρ₀(N/N₀) since the mean column density is N₀ = ρ₀L.

The quantity 〈(N/N₀)²〉 can be directly calculated via

\begin{equation} \langle (N/N_{0})^{2} \rangle = \frac{\langle N^{2} \rangle }{\langle N_{0} \rangle ^{2}} = \frac{\langle W^{2}_{0} \rangle }{\langle W_{0} \rangle ^{2}}, \end{equation}

(44)

and we have a means to calculate 〈(ρ/ρ₀)²〉 using equation (34) since |$\langle (\rho /\rho _{0})^{2} \rangle = 1 + \sigma ^{2}_{\rho /\rho _{0}}$| and |$\sigma ^{2}_{N/N_{0}} = \langle (N/N_{0})^{2} \rangle - 1$|⁠. Note that f(k) is the angular average of the column density power spectrum and since the overall normalization is unimportant, the integrated intensity power spectrum, |$\tilde{W}_{0} \tilde{W}^{*}_{0}$|⁠, can be used to measure this.

With 〈(N/N₀)²〉 and 〈(ρ/ρ₀)²〉 known, the desired ratio |$\sigma ^{2}_{p_{z,p}}/\sigma ^{2}_{p_{z}}$| is given by (see equation 43)

\begin{equation} \frac{\sigma ^{2}_{p_{z,p}}}{\sigma ^{2}_{p_{z}}} = \left[\frac{\sigma ^{2}_{p_{z,p}}}{\langle \rho ^{2}_{p} \rangle } \right] \left[\frac{\langle (N/N_{0})^{2} \rangle }{\langle (\rho /\rho _{0})^{2} \rangle } \right] \left[ \frac{\sigma ^{2}_{p_{z}}}{\langle \rho ^{2} \rangle } \right]^{-1}. \end{equation}

(45)

It still remains to find an observational estimate of |$\sigma ^{2}_{p_{z}}/\langle \rho ^{2} \rangle$|⁠. This is a 3D quantity – the ratio of mean squared z-momentum to mean squared density. Explicitly,

\begin{equation} \frac{\sigma ^{2}_{p_{z}}}{\langle \rho ^{2} \rangle } = \frac{\frac{1}{L^{3}} \displaystyle \int _{-L/2}^{L/2} {\rm{d}}x \; \displaystyle \int _{-L/2}^{L/2} {\rm{d}}y \; \displaystyle \int _{-L/2}^{L/2} {\rm{d}}z \; p_{z}^{2}}{\frac{1}{L^{3}} \displaystyle \int _{-L/2}^{L/2} {\rm{d}}x \; \displaystyle \int _{-L/2}^{L/2} {\rm{d}}y \; \displaystyle \int _{-L/2}^{L/2} {\rm{d}}z \; \rho ^{2}}, \end{equation}

(46)

where we have again assumed the calculations are performed in the zero momentum frame. From a slightly different perspective, we see that it may also be viewed as the z-velocity dispersion calculated with a ρ² weight:

\begin{equation} \frac{\sigma ^{2}_{p_{z}}}{\langle \rho ^{2} \rangle } = \frac{\frac{1}{L^{3}} \displaystyle \int _{-L/2}^{L/2} {\rm{d}}x \; \displaystyle \int _{-L/2}^{L/2} {\rm{d}}y \; \displaystyle \int _{-L/2}^{L/2} {\rm{d}}z \; \rho ^{2} v_{z}^{2}}{\frac{1}{L^{3}} \displaystyle \int _{-L/2}^{L/2} {\rm{d}}x \; \displaystyle \int _{-L/2}^{L/2} {\rm{d}}y \; \displaystyle \int _{-L/2}^{L/2} {\rm{d}}z \; \rho ^{2}} = \frac{\langle \rho ^{2} v^{2}_{z} \rangle }{\langle \rho ^{2} \rangle }. \end{equation}

(47)

However, we only have access instead to the z-velocity dispersion calculated with a ρ weight, as follows. Making use of the second velocity-moment of intensity, W₂:

\begin{eqnarray} W_{2}(x, y)&=& \displaystyle \int _{-\infty }^{\infty } {\rm{d}}v \; I(x, y, v) v^{2} \nonumber \\ &=& \displaystyle \int _{-\infty }^{\infty } {\rm{d}}v \; e \displaystyle \int _{-L/2}^{L/2} {\rm{d}}z \; \rho (x, y, z) \delta (v - v_{z}) v^{2} \nonumber \\ &=& e \displaystyle \int _{-L/2}^{L/2} {\rm{d}}z \; \rho (x, y, z) v^{2}_{z}(x, y, z), \end{eqnarray}

(48)

a spatial average of W₂, normalized by the spatial average of W₀, gives

\begin{eqnarray} \frac{\langle W_{2} \rangle }{\langle W_{0} \rangle } &=& \frac{ \frac{e}{L^{2}} \displaystyle \int _{-L/2}^{L/2} {\rm{d}}x \; \displaystyle \int _{-L/2}^{L/2} {\rm{d}}y \; \displaystyle \int _{-L/2}^{L/2} {\rm{d}}z \; \rho (x, y, z) v^{2}_{z}(x, y, z)}{ \frac{e}{L^{2}} \displaystyle \int _{-L/2}^{L/2} {\rm{d}}x \; \displaystyle \int _{-L/2}^{L/2} {\rm{d}}y \; \displaystyle \int _{-L/2}^{L/2} {\rm{d}}z \; \rho (x, y, z) } \nonumber \\ &=& \frac{\langle \rho v^{2}_{z} \rangle }{\langle \rho \rangle }. \end{eqnarray}

(49)

Thus, we have a measure of |$\langle \rho v^{2}_{z} \rangle /\langle \rho \rangle$| but require a measure of |$\langle \rho ^{2} v^{2}_{z} \rangle /\langle \rho ^{2} \rangle$|⁠. While there are circumstances in which |$\langle \rho ^{2} v^{2}_{z} \rangle /\langle \rho ^{2} \rangle = \langle \rho v^{2}_{z} \rangle /\langle \rho \rangle$|⁠, in general they will be different. They can be considered practically equal in two special cases: (1) a uniform or only very weakly varying density field, (2) statistically uncorrelated density and velocity fields. Neither of these conditions is expected to hold in the ISM, and corrective measures are necessary. For now, we write, in advance of notation explained in the appendix:

\begin{equation} \frac{\sigma ^{2}_{p_{z}}}{\langle \rho ^{2} \rangle } = \frac{\langle \rho ^{2} v^{2}_{z} \rangle }{\langle \rho ^{2} \rangle } = g_{21} \frac{\langle \rho v^{2}_{z} \rangle }{\langle \rho \rangle }, \end{equation}

(50)

where g₂₁ is a statistical correction factor of order unity. In the appendix, we show that g₂₁ may be written as

\begin{equation} g_{21} = \langle \xi ^{2} \rangle ^{-\epsilon }, \end{equation}

(51)

where ξ = ρ/ρ₀ and ϵ is a small, positive constant. (We find for the numerical simulations that ϵ is Mach number dependent.) In the appendix, we discuss both numerical and observational estimates of ϵ and look at the effects of g₂₁ on the analysis in Section 4.

Above, we have computed |$\langle \rho v^{2}_{z} \rangle /\langle \rho \rangle$| via the observable ratio W₂/W₀, but note that it may also be computed simply as the dispersion of the summed spectral line profile of the data. It is also important to recognize that a finite-width thermal/instrumental broadening term will cause overestimation of |$\langle \rho v^{2}_{z} \rangle /\langle \rho \rangle$|⁠. The simple fix for this is to subtract the thermal/instrumental dispersion from the raw measurement of |$\langle \rho v^{2}_{z} \rangle /\langle \rho \rangle$| (i.e. the influence of finite thermal/instrumental broadening is simply to convolve the δ-function-mapped data with a smoothing kernel on the v_z-axis which will be typically Gaussian in form). We will return to this point, and other practical observational considerations in Section 4 below.

Putting together all the above, equation (45) becomes

\begin{equation} \frac{\sigma ^{2}_{p_{z,p}}}{\sigma ^{2}_{p_{z}}} = \left[\frac{\sigma ^{2}_{p_{z,p}}}{\langle \rho ^{2}_{p} \rangle } \right] \left[\frac{\langle (N/N_{0})^{2} \rangle }{\langle (\rho /\rho _{0})^{2} \rangle } \right] \left[ g_{21} \frac{\langle \rho v^{2}_{z} \rangle }{\langle \rho \rangle } \right]^{-1}, \end{equation}

(52)

so that the fraction of z-momentum power in transverse modes is given by equation (31):

\begin{equation} \frac{\sigma ^{2}_{p_{z\perp }}}{\sigma ^{2}_{p_{z}}} = \frac{\sigma _{p_{z,p}}^{2}}{{\sigma _{p_{z}}^{2}}} \times \frac{\displaystyle \sum _{k_{x}=-k_{{\rm max}}}^{k_{{\rm max}}} \sum _{k_{y}=-k_{{\rm max}}}^{k_{{\rm max}}} \sum _{k_{z}=-k_{{\rm max}}}^{k_{{\rm max}}} \; f_{\perp }(k) \frac{k^{2}_{x}+k^{2}_{y}}{k^{2}}}{\displaystyle \sum _{k_{x}=-k_{{\rm max}}}^{k_{{\rm max}}} \sum _{k_{y}=-k_{{\rm max}}}^{k_{{\rm max}}} \; f_{\perp }(k)}, \end{equation}

(53)

where k_max is the maximum observed 1D wavenumber, and f_⊥(k) is the angular average of the projected momentum power spectrum. Assuming isotropy, equation (53) provides an estimate of the ratio of transverse to total momentum power in 3D, i.e.

\begin{equation} \frac{\sigma ^{2}_{p_{\perp }}}{\sigma ^{2}_{p}} \approx \frac{\sigma ^{2}_{p_{z\perp }}}{\sigma ^{2}_{p_{z}}}. \end{equation}

(54)

2.4 Summary of observational measurements required

Equations (52)–(54) can be written in terms of known observables as

\begin{eqnarray} \frac{\sigma ^{2}_{p_{\perp }}}{\sigma ^{2}_{p}} \approx \left[\frac{\langle W^{2}_{1} \rangle }{\langle W^{2}_{0} \rangle }\right] \left[ \frac{\langle W^{2}_{0} \rangle / \langle W_{0} \rangle ^{2}}{1 + A(\langle W^{2}_{0} \rangle / \langle W_{0} \rangle ^{2} - 1)} \right] \left[ g_{21} \frac{\langle W_{2} \rangle }{\langle W_{0} \rangle } \right]^{-1} B, \nonumber\\ \end{eqnarray}

(55)

where

\begin{equation} A = \frac{\left(\displaystyle \sum _{k_{x}=-k_{{\rm max}}}^{k_{{\rm max}}} \sum _{k_{y}=-k_{{\rm max}}}^{k_{{\rm max}}} \sum _{k_{z}=-k_{{\rm max}}}^{k_{{\rm max}}} \; f(k)\right) - f(0)}{\left(\displaystyle \sum _{k_{x}=-k_{\rm max}}^{k_{\rm max}} \sum _{k_{y}=-k_{\rm max}}^{k_{\rm max}} \; f(k)\right) - f(0)}, \end{equation}

(56)

\begin{equation} B = \frac{\displaystyle \sum _{k_{x}=-k_{\rm max}}^{k_{\rm max}} \sum _{k_{y}=-k_{\rm max}}^{k_{\rm max}} \sum _{k_{z}=-k_{\rm max}}^{k_{\rm max}} \; f_{\perp } (k) \frac{k^{2}_{x}+k^{2}_{y}}{k^{2}}}{\displaystyle \sum _{k_{x}=-k_{\rm max}}^{k_{\rm max}} \sum _{k_{y}=-k_{\rm max}}^{k_{\rm max}} \; f_{\perp }(k)}, \end{equation}

(57)

and

\begin{equation} f(k) = \frac{1}{2 {\pi } k} \displaystyle \int _{0}^{2 {\pi }} {\rm{d}}\phi \; \tilde{W}_{0}(k,\phi ) \tilde{W}^{*}_{0}(k,\phi ), \end{equation}

(58)

\begin{equation} f_{\perp }(k) = \frac{1}{2 {\pi } k} \displaystyle \int _{0}^{2 {\pi }} {\rm{d}}\phi \; \tilde{W}_{1}(k,\phi ) \tilde{W}^{*}_{1}(k,\phi ), \end{equation}

(59)

are, respectively, the angular averages of the power spectra of the zeroth and first velocity moments of the spectral line intensities (in practice, computed by sums rather than integrals).

Finally, we note that in equation (55), it is evident that some factors cancel, suggesting further simplification is possible. We have opted to leave it in the form presented since we believe it makes more logical sense this way, and the individual terms (in square brackets) in the equation will be analysed below, along with the deprojection factor, B, in Section 4.

3 NUMERICAL SIMULATIONS

The simulations used to test the analytic method were performed with the astrophysical code flash (Fryxell et al. 2000; Dubey et el. 2008), which integrates the ideal, 3D, magnetohydrodynamic (MHD) equations of compressible gas. The MHD equations were closed with a polytropic equation of state, |$P_\rm{th}=c_\rm{s}^2\rho$|⁠, such that the gas remains isothermal with a constant sound speed, c_s. We solve the MHD equations in the hydrodynamic limit (⁠|${\boldsymbol {B}} = 0$|⁠) on 3D, uniform grids with a fixed resolution of 256³ grid points and periodic boundary conditions, using a positive-definite Riemann solver for ideal MHD (Waagan, Federrath & Klingenberg 2011).

To drive turbulence, a stochastic forcing term |${\boldsymbol {F}_\rm{stir}}$| is applied as a source term in the MHD momentum equation. Following common practice, the forcing only acts on large scales 1 < k < 3 (where most of the power is injected at the k = 2 mode in Fourier space, which corresponds to half of the box size L in physical space), i.e. the outer scale of molecular clouds, as favoured by observations (Ossenkopf & Mac Low 2002; Brunt, Heyer & Mac Low 2009), such that turbulence develops self-consistently on smaller scales. We use the Ornstein–Uhlenbeck (OU) process to model |${\boldsymbol {F}_\rm{stir}}$|⁠, which is a well-defined stochastic process with a finite autocorrelation time-scale (Eswaran & Pope 1988; Schmidt, Hillebrandt & Niemeyer 2006), leading to a smoothly varying stochastic force field in space and time. Details about the OU process and the forcing applied in this study can be found in Schmidt et al. (2009), Federrath et al. (2010) and Konstandin et al. (2012a). However, the essential point of our forcing approach is that we can adjust the mixture of solenoidal and compressive modes of |${\boldsymbol {F}_\rm{stir}}$| arbitrarily. This is achieved with the projection tensor |$\mathcal {\underline{P}}^{\,\zeta }({\boldsymbol {k}})$| in Fourier space. In index notation, it reads

\begin{equation} \mathcal {P}_{ij}^\zeta = \zeta \,\mathcal {P}_{ij}^\perp +(1-\zeta )\,\mathcal {P}_{ij}^\parallel = \zeta \,\delta _{ij}+(1-2\zeta )\,\frac{k_i k_j}{|k|^2}, \end{equation}

(60)

where δ_ij is the Kronecker symbol, and |$\mathcal {P}_{ij}^\perp = \delta _{ij} - k_i k_j / k^2$| and |$\mathcal {P}_{ij}^\parallel = k_i k_j / k^2$| are the fully solenoidal and the fully compressive projection operators, respectively. The ratio of compressive power to total power in |${\boldsymbol {F}_\rm{stir}}$| can be derived from equation (60) by evaluating the norm of the compressive component of the projection tensor and dividing it by the total injected power, resulting in

\begin{equation} \frac{F_\rm{comp}}{F_\rm{tot}} = \frac{(1-\zeta )^2}{1-2\zeta +3\zeta ^2}\,, \end{equation}

(61)

for 3D space (Schmidt et al. 2009; Federrath et al. 2010). The projection operator serves to construct a purely solenoidal force field by setting ζ = 1, while for ζ = 0, a purely compressive force field is obtained. Any combination of solenoidal and compressive modes can be constructed by choosing ζ ∈ [0, 1]. Here, we compare simulations with ζ = 1 (sol) and ζ = 0 (comp).

Starting from a uniform density distribution and zero velocities, the forcing excites turbulent motions, which approach a statistically steady state after about two turbulent crossing times, |$2T=L/(\mathcal {M}c_\rm{s})$| (e.g. Klessen, Heitsch & Mac Low 2000; Klessen 2001; Heitsch, Mac Low & Klessen 2001; Federrath, Klessen & Schmidt 2009; Federrath et al. 2010; Price & Federrath 2010; Micic et al. 2012; Federrath 2013), where |$\mathcal {M}$| denotes the 3D, root-mean-squared sonic Mach number. We study simulations in both the subsonic and supersonic regimes of turbulence with |$\mathcal {M}\sim 0.1$|⁠, 0.5, 2, 5 and 15, each with the two limiting cases of purely solenoidal and purely compressive forcing to basically cover the whole range of possible solenoidal momentum ratios between zero and unity, in order to test the analytic method. To ensure that the initial transient phase (t ≲ 2T) is not included in the following analysis, we consider snapshots at t = 3, 4 and 5 T. Each of the snapshots used in the analysis is thus separated by at least one crossing time, effectively representing statistically independent turbulent fields at this temporal separation. We thus improve the independent statistical sampling of our results by including these three snapshots for each simulation, providing an estimate of the typical temporal variations of our results.

We also make use of previously conducted simulations for auxiliary information and testing purposes. These simulations have been previously described in BFP and comprise smoothed particle hydrodynamics (SPH) simulations of solenoidally forced turbulence at a range of (supersonic) Mach numbers and fixed grid simulations of solenoidally forced MHD turbulence at a range of (supersonic) Mach numbers and a range of magnetic field strengths (see BFP for more details).

Finally in this section, we compute example power spectra of longitudinal and transverse momentum components for comparison with Fig. 1. The Fourier transform (⁠|$\tilde{{\boldsymbol {p}}}_{||}$|⁠) of the longitudinal component of the momentum (⁠|${\boldsymbol {p}}_{||}$|⁠) is given by

\begin{equation} {\tilde{\boldsymbol {p}}}_{||} = {\hat{\boldsymbol {k}}}({\hat{\boldsymbol {k}}} \cdot \tilde{{\boldsymbol {p}}}), \end{equation}

(62)

where |${\boldsymbol {p}} = \rho {\boldsymbol {v}}$|⁠, and |$\tilde{{\boldsymbol {p}}}$| is its Fourier transform. The longitudinal momentum is then found by |${\boldsymbol {p}}_{||} $| = FT|$^{-1}(\tilde{{\boldsymbol {p}}}_{||})$|⁠. We then extract the transverse component via

\begin{equation} {\boldsymbol {p}}_{\perp } = {\boldsymbol {p}} - {\boldsymbol {p}}_{||}. \end{equation}

(63)

From the vector momenta, we can extract the z-components (p_z|| and p_z⊥) and measure their power spectra. Fig. 2 shows example 3D renderings for a solenoidally driven p_z-field at rms Mach number 5. The power spectra exhibit the same underlying symmetries as the power spectra shown in Fig. 1, including the suppression (nulling) of longitudinal power as k_z approaches (equals) zero. (This geometrical aspect is enforced even if the power spectra are anisotropic.)

$3D renderings of constant power surfaces in the power spectra of pz⊥ (left) and $p_{z_{||}}$ (right) in a snapshot taken from the numerical simulations (solenoidally forced, rms Mach number = 5). In both panels, the origin of coordinates ($\boldsymbol {k} = 0$) lies in the centre of the image (cf. Fig. 1).$

Figure 2.

3D renderings of constant power surfaces in the power spectra of p_z⊥ (left) and |$p_{z_{||}}$| (right) in a snapshot taken from the numerical simulations (solenoidally forced, rms Mach number = 5). In both panels, the origin of coordinates (⁠|$\boldsymbol {k} = 0$|⁠) lies in the centre of the image (cf. Fig. 1).

Open in new tab Download slide

4 APPLICATION TO HYDRODYNAMIC SIMULATIONS

We now apply the above method to the numerical simulations described in Section 3. Before doing so, we examine the fraction of momentum power contained in transverse modes by direct calculation in 3D. These fractions are calculated via Fourier space, making use of the conditions given in equations (5) and (6). In Fig. 3, we plot the ratio of transverse to total momentum power, |$\sigma ^{2}_{p_{\perp }} / \sigma ^{2}_{p}$|⁠, against the density-weighted Mach number, |${\mathcal {M}_{1}}$|⁠. At low Mach number (⁠|${\mathcal {M}_{1}} <\!{<} 1$|⁠), the different forcing methods result in very different values of |$\sigma ^{2}_{p_{\perp }} / \sigma ^{2}_{p}$| – essentially all the power remains in the respective forcing modes. As |${\mathcal {M}_{1}}$| increases, the transverse fractions appear to converge to |$\sigma ^{2}_{p_{\perp }} / \sigma ^{2}_{p} \approx 0.7\pm 0.1$|⁠, irrespective of the nature of the forcing. It is perhaps unfortunate that this situation arises; interestingly, the transverse velocity power fraction is different for different forcing mechanisms, even at high Mach number (Federrath et al. 2011). We discuss the prospects of measuring the transverse velocity power fraction below.

$Plot of the fraction of momentum power in transverse modes versus 3D density-weighted Mach number, ${\mathcal {M}}_{1}$. Dots: solenoidal forcing; open squares: compressive forcing.$

Figure 3.

Plot of the fraction of momentum power in transverse modes versus 3D density-weighted Mach number, |${\mathcal {M}}_{1}$|⁠. Dots: solenoidal forcing; open squares: compressive forcing.

Open in new tab Download slide

Now we attempt to measure |$\sigma ^{2}_{p_{\perp }} / \sigma ^{2}_{p}$| using only observationally accessible quantities. The analysis below is somewhat idealized, in that we assume certain auxiliary pieces of physical information are available, and we do not consider the effects of instrumental noise and beam smearing – though we will make comments and recommendations at appropriate points.

4.1 Normalized density variance

First, we use the BFP method to estimate the 3D normalized density variance, |$\sigma ^{2}_{\rho /\rho _{0}}$|⁠, using information contained solely in the observationally accessible normalized column density field, namely |$\sigma ^{2}_{N/N_{0}}$| and the angular average of the column density power spectrum, f(k). The second term in equation (52) is then formed as

\begin{equation*} \frac{\langle (N/N_{0})^{2} \rangle }{\langle (\rho /\rho _{0})^{2} \rangle } = \frac{1+\sigma ^{2}_{N/N_{0}}}{1+\sigma ^{2}_{\rho /\rho _{0}}}.\nonumber \end{equation*}

In Fig. 4, we plot the observational value of this ratio versus the true quantity measured by privileged access to the 3D density field. To make the data in this plot, we make use of all three possible orientations of the simulation cubes to generate three column density fields per simulation. Note that large values of 〈(N/N₀)²〉/〈(ρ/ρ₀)²〉 originate from the weakly varying density fields (low Mach number) while the smaller values originate from more variable density fields (higher Mach numbers) in which small-scale structure suffers from a greater degree of line-of-sight averaging. In the observational context, such measurements are subject to instrumental noise (resulting in a noise contribution to the column density variance and column density power spectrum) and the effect of beam smearing of the telescope. These factors can be dealt with using the methods outlined in BFP and Brunt (2010).

Observed versus intrinsic values of the ratio 〈(N/N0)2〉/〈(ρ/ρ0)2〉 measured in the numerical simulations. The ‘observed’ values are calculated using the BFP method, which is accurate to about 10 per cent for statistically isotropic fields.

Figure 4.

Observed versus intrinsic values of the ratio 〈(N/N₀)²〉/〈(ρ/ρ₀)²〉 measured in the numerical simulations. The ‘observed’ values are calculated using the BFP method, which is accurate to about 10 per cent for statistically isotropic fields.

Open in new tab Download slide

4.2 Dispersion in z-axis momentum

The next term in equation (52) we examine is the third term, |$g_{21} \langle \rho v_{z}^{2} \rangle / \langle \rho \rangle$|⁠. This term is designed to measure the total dispersion in z-axis momentum (divided by 〈ρ²〉) as described by equation (50) – or equivalently, the z-axis velocity dispersion calculated with a ρ² weight. As mentioned previously, 〈ρv²〉/〈ρ〉 can either be obtained from the ratio 〈W₂〉/〈W₀〉 or by the dispersion of the summed spectral line profile. The latter option is probably best-suited to practical observational work, though the two are equivalent. Corrections for thermal/instrumental broadening should be made. We point out here that for the low Mach number simulations (⁠|${\mathcal {M}} {<} 1$|⁠), the thermal broadening dominates over the turbulent motions of interest. Here, we assume no thermal broadening (or exact accounting for it) which is rather unrealistic. However, our aim here is to test the principle rather than the practice, since the low Mach number simulations extend the range of |$\sigma ^{2}_{p_{\perp }} / \sigma ^{2}_{p}$| available for analysis (see Fig. 3). We do not recommend that the method is applied practically to subsonic media, but expect instead that it will be applied to large molecular clouds where the thermal broadening has a small influence on the dominant supersonic turbulent motions.

A further idealization that we use here is that we work with Mach numbers rather than velocity dispersions. This is just a convenient system in which to make comparisons between velocity dispersions (in units of the squared sound speed), though we note that our estimation of the parameter ϵ needed to derive g₂₁ via equation (A16) requires the Mach number |${\mathcal {M}}_{1}$| to be known. In the observational context, this would require that an accurate measure of the sound speed is available. We have made a further assumption here, which is: the use of equation (A16) assumes the hydrodynamic limit, since the fit was made to simulations that did not include magnetic fields. There is an apparent dependence of ϵ on the Alfvénic Mach number (see Fig. A4) but this is only important in the limit |${\mathcal {M}}_{\rm A} {<} 1$| and low sonic Mach number (⁠|${\mathcal {M}}_{1} \lesssim 5$|⁠). We assume that auxiliary observations have been made to ensure that these conditions are not met, though the consequences of erroneously assuming they are not met when in fact they (at least within the range of conditions covered by the above analysis) are not severe. In the low |${\mathcal {M}}_{\rm A}$|⁠, low |${\mathcal {M}}_{1}$| regime, ϵ is slightly smaller than that which holds in the hydrodynamic limit (ϵ_∞), and therefore the true correction factor lies between unity and g₂₁(ϵ_∞). As shown below, the practical difference between these two cases is tolerably small.

From the simulation data, we can directly extract the quantity |${\mathcal {M}}_{1z} = \langle \rho v^{2}_{z} \rangle / \langle \rho \rangle c^{2}_{\rm s}$|⁠, which is the ρ-weighted z-velocity dispersion in units of the squared sound speed. The correction factor g₂₁ is designed to convert this into the ρ²-weighted z-axis Mach number: |${\mathcal {M}}_{2z} = \langle \rho ^{2} v^{2}_{z} \rangle / \langle \rho ^{2} \rangle c^{2}_{\rm s}$|⁠. To calculate ϵ we assume isotropy and take |${\mathcal {M}}_{1} = \sqrt{3} {\mathcal {M}}_{1z}$|⁠, allowing us to calculate g₂₁ = 〈(ρ/ρ₀)²〉^−ϵ using equation (A15) with 〈(ρ/ρ₀)²〉 already calculated above. In Fig. 5(a), we plot |$M^{2}_{1z}$| versus |$M^{2}_{2z}$| (measured with privileged access) for the supersonic data, to illustrate the consequences of the naive assumption that |$M^{2}_{1z}$| is a good surrogate for |$M^{2}_{2z}$|⁠. Fig. 5(b) shows the effect of assuming that |$M^{2}_{2z} = g_{21} M^{2}_{1z}$|⁠, which is rather better, though there are a couple of outliers that resist correction (these are notably anisotropic fields). In Fig. 5, we have again used all 3 available choices for the z-axis in each simulation. In the analysis below, we will document the effect of including or ignoring the correction factor g₂₁.

$(a) Estimate of ${\mathcal {M}}_{2z}^{2}$ naively assuming ${\mathcal {M}}_{2z}^{2} = {\mathcal {M}}_{1z}^{2}$; (b) estimate of ${\mathcal {M}}_{2z}^{2}$ assuming ${\mathcal {M}}_{2z}^{2} = g_{21} {\mathcal {M}}_{1z}^{2}$.$

Figure 5.

(a) Estimate of |${\mathcal {M}}_{2z}^{2}$| naively assuming |${\mathcal {M}}_{2z}^{2} = {\mathcal {M}}_{1z}^{2}$|⁠; (b) estimate of |${\mathcal {M}}_{2z}^{2}$| assuming |${\mathcal {M}}_{2z}^{2} = g_{21} {\mathcal {M}}_{1z}^{2}$|⁠.

Open in new tab Download slide

4.3 Projected solenoidal fractions

The first factor in equation (52) is trivially computed using the observationally accessible |$\langle W^{2}_{1} \rangle / \langle W^{2}_{0} \rangle$|⁠. We have already mentioned accounting for the contribution of noise variance to |$\langle W^{2}_{0} \rangle$|⁠, and it remains to deal with the role of noise in calculating |$\langle W^{2}_{1} \rangle$|⁠, which is a bit more complicated. We just mention here that prescriptions for accounting for noise in centroid velocity measurements are available (e.g. Kleiner & Dickman 1985; Miesch & Bally 1994; Brunt & Mac Low 2004) and a suitable modification of these should be made.

Now with all three terms in equation (52) assembled, we can construct an observationally accessible estimate of |$\sigma ^{2}_{p_{z,p}} / \sigma ^{2}_{p_{z}}$| (i.e. the fraction of z-momentum power projected into 2D) for comparison with the same ‘true’ quantity directly extracted from the data with privileged access. Fig. 6 shows this comparison, (a) without and (b) with the correction factor g₂₁ applied. In general, equation (52) provides an accurate observational measure of the ratio |$\sigma ^{2}_{p_{z,p}} / \sigma ^{2}_{p_{z}}$|⁠.

$Observationally estimated values of the fraction of total z-momentum power projected into 2D assuming (a) g21 = 1 (i.e. no statistical correction for density–velocity correlation) and (b) g21 = 〈ξ2〉−ϵ, in both cases plotted versus the true value of this fraction.$

Figure 6.

Observationally estimated values of the fraction of total z-momentum power projected into 2D assuming (a) g₂₁ = 1 (i.e. no statistical correction for density–velocity correlation) and (b) g₂₁ = 〈ξ²〉^−ϵ, in both cases plotted versus the true value of this fraction.

Open in new tab Download slide

Note that there are two main effects that control this ratio. First, if the fraction of compressive (longitudinal) momentum modes is high, then the ratio |$\sigma ^{2}_{p_{z,p}} / \sigma ^{2}_{p_{z}}$| will be low since only transverse modes are projected into 2D. This is (partially) why the compressively forced simulations display small values of |$\sigma ^{2}_{p_{z,p}} / \sigma ^{2}_{p_{z}}$|⁠. Secondly, a higher density variance will typically result in a lower fraction of momentum variance projected into 2D since such fields suffer more averaging/smoothing of small-scale structure (e.g. as mentioned above in the discussion of Fig. 4). This means that even though the transverse momentum fractions are approximately the same at high Mach numbers regardless of the forcing mechanism (see Fig. 3) the compressive forcing drives higher density variance (Federrath, Klessen & Schmidt 2008) which results in a greater degree of line-of-sight averaging of small-scale structure, and therefore lower values of |$\sigma ^{2}_{p_{z,p}} / \sigma ^{2}_{p_{z}}$|⁠. The highest values of |$\sigma ^{2}_{p_{z,p}} / \sigma ^{2}_{p_{z}}$| reached (∼1/3) originate from the subsonic solenoidally forced simulations which have very nearly uniform density. In this case, trivially, ∼1/3 of the momentum variance is projected along one of the three spatial axes.

4.4 Deprojection

The final step in the analysis is the de-projection via equation (53) to estimate the fraction of momentum power in transverse modes in 3D (assuming that |$\sigma ^{2}_{p_{\perp }} / \sigma ^{2}_{p} \approx \sigma ^{2}_{p_{z\perp }} / \sigma ^{2}_{p_{z}}$|⁠). The input to the de-projection factor is just f_⊥(k), obtained from the angular average of the power spectrum of W₁. In an observational context, this should be corrected for noise (using a suitable modification of the methods in Brunt 2010; Brunt & Mac Low 2004) and treatment of the effect of the telescope beam pattern should be included (see BFP for a discussion of this).

In Fig. 7, we compare the observationally derived values of |$\sigma ^{2}_{p_{\perp }} / \sigma ^{2}_{p}$| to those measured exactly with privileged access to the momentum field in 3D. Fig. 7(a) shows the comparison with no statistical correction for density–velocity correlation, while Fig. 7(b) shows the results with the g₂₁ correction applied (we have applied the correction to all data, not just the supersonic fields). In the plots, the plotted points are the mean values of |$\sigma ^{2}_{p_{\perp }} / \sigma ^{2}_{p}$| obtained by averaging the results over all three spatial axes, while the error bars represent the standard deviation about the mean. This is to show the recovery of |$\sigma ^{2}_{p_{\perp }} / \sigma ^{2}_{p}$| from the same field seen from different orientations. Overall, the observational recovery of the intrinsic |$\sigma ^{2}_{p_{\perp }} / \sigma ^{2}_{p}$| values is good, albeit with relatively large scatter for the high Mach number compressively forced simulations.

$Observationally estimated values of the fraction of 3D momentum power in transverse modes assuming (a) g21 = 1 (i.e. no statistical correction for density–velocity correlation) and (b) g21 = 〈ξ2〉−ϵ, in both cases plotted versus the true value of this fraction.$

Figure 7.

Observationally estimated values of the fraction of 3D momentum power in transverse modes assuming (a) g₂₁ = 1 (i.e. no statistical correction for density–velocity correlation) and (b) g₂₁ = 〈ξ²〉^−ϵ, in both cases plotted versus the true value of this fraction.

Open in new tab Download slide

Though the solenoidally forced simulations alone span a limited range in |$\sigma ^{2}_{p_{\perp }} / \sigma ^{2}_{p}$|⁠, the recovery is reliable, albeit with a slight overestimation for intrinsic |$\sigma ^{2}_{p_{\perp }} / \sigma ^{2}_{p} \approx 1$| (applicable to the subsonic, ∼ uniform density fields). If the turbulence was strongly magnetized (i.e. |$\mathcal {M}_{\rm s} \lesssim 5$|⁠, |$\mathcal {M}_{\rm A} \lesssim 1$|⁠), but assumed not to be, then the g₂₁ factors derived from ϵ_∞ would lead to slightly overestimated solenoidal fractions, though for the range of Mach numbers studied here, this is only at the ∼10 per cent level. This should be tested directly in future studies, where the (probably) more important question of anisotropy should be assessed.

5 DISCUSSION

While the simulated fields here contain a number of idealizations (Fourier-space driving, periodic boundary conditions, lack of self-gravity, strict isothermality) the above analysis has demonstrated that, in principle, the fraction of momentum power in transverse (solenoidal) modes, |$\sigma ^{2}_{p_{\perp }} / \sigma ^{2}_{p}$|⁠, may be measurable from observations. (We have outlined practical observational considerations at appropriate points above.) For the subsonic fields, to extend the range of intrinsic |$\sigma ^{2}_{p_{\perp }} / \sigma ^{2}_{p}$| available, we have assumed (unrealistically) that the dominant thermal broadening terms have been accounted for.

We envision that the above method is best-suited to the study of nearby supersonically turbulent molecular clouds for which high-sensitivity, high spatial dynamic range spectral line data are available. So, how relevant are the simulated supersonic fields to typical conditions met in molecular clouds that (currently) could be analysed with the above model? We argue that the key condition that must be satisfied is that of statistical isotropy, rather than any shortcomings in the physical simulations – i.e. the simulations simply create an intrinsic ratio |$\sigma ^{2}_{p_{\perp }} / \sigma ^{2}_{p}$| that we can set out to measure. If strong magnetic fields are present, then this could cause significant anisotropy in the density and velocity fields that we have not accounted for here. The assumption of isothermality may be slightly problematic, since the method relies (weakly) on applying corrections for thermal broadening and deriving the Mach number (to estimate g₂₁). The assumption of a single Mach number is also problematic if the method is to be applied to the atomic ISM. It is also possible of course that the detailed physics may affect the statistical correction factors g₂₁, a point to which we will return below.

5.1 Anisotropies and discontinuities

The next steps in the testing of the method should focus on possible complications and biases induced by anisotropies. In fact, there are suggestions that compressive modes may be suppressed in highly anisotropic media (Hansen, McKee & Klein 2011), which would be very interesting to test observationally once the practicalities of doing so are understood more fully. Testing (or extension) of the method to MHD case with different compressive/solenoidal fractions is also worth pursuing. Filamentary structure, by itself, poses no particular problem as long as the filament orientations are statistically isotropic. In practice, good evidence for isotropy in projected 2D is needed, as well as some confirmation that the cloud's line-of-sight depth is comparable to its projected extent. In the field selection, one must also pay attention to the boundary conditions. In our simulations, these are periodic, so there are no problems arising from edge discontinuities. Fields for analysis should be selected so that no significant edge discontinuities exist, though we are aided in this selection by the ρ- (or N-) weighting of the relevant fields so that (to the extent that any ‘cloud’ is truly isolated) a suitable region may be defined where W₀, W₁ and W₂ are sufficiently close to zero that this condition is satisfied. Edge-tapering (e.g. Brunt & Mac Low 2004) or padding (e.g. Brunt 2010) may be applied to sufficiently large fields, with negligible quantitative consequences.

5.2 Variability

In terms of the reliability of the recovered |$\sigma ^{2}_{p_{\perp }} / \sigma ^{2}_{p}$|⁠, of most concern is the relatively high degree of scatter in the compressively forced fields. The scatter is mainly a consequence of the high degree of variability in these fields, so that variances are contributed to strongly by a small number of extreme field values, perhaps magnifying any anisotropic effects. To quantify this, in Fig. 8 we plot the root-mean-square uncertainty on the measured |$\sigma ^{2}_{p_{\perp }} / \sigma ^{2}_{p}$| ratios against a measure of the column density variability, |$\ln {\left(1 + \sigma ^{2}_{N/N_{0}}\right)}$| (which would be the logarithmic variance, |$\sigma _{ln(N/N_{0})}$|⁠, in the case of a lognormal PDF). Fig. 8 shows that typical errors are at the ∼10 per cent level or better for |$\ln {\left(1 + \sigma ^{2}_{N/N_{0}}\right)} \lesssim 1$|⁠, but then increase sharply for |$\ln {\left(1 + \sigma ^{2}_{N/N_{0}}\right)} \gtrsim 1$|⁠.

$The root-mean-square uncertainty on the measured $\sigma ^{2}_{p_{\perp }} / \sigma ^{2}_{p}$ ratios from Fig. 7, plotted against a measure of the column density variability, $\ln {(1 + \sigma ^{2}_{N/N_{0}})}$.$

Figure 8.

The root-mean-square uncertainty on the measured |$\sigma ^{2}_{p_{\perp }} / \sigma ^{2}_{p}$| ratios from Fig. 7, plotted against a measure of the column density variability, |$\ln {(1 + \sigma ^{2}_{N/N_{0}})}$|⁠.

Open in new tab Download slide

5.3 Variations in solenoidal/compressive fractions

It is evident from Fig. 3 that the |$\sigma ^{2}_{p_{\perp }} / \sigma ^{2}_{p}$| ratio does not allow observational discrimination between solenoidal and compressive forcing, since the |$\sigma ^{2}_{p_{\perp }} / \sigma ^{2}_{p}$| ratios appear to converge to 0.7 ± 0.1 at high Mach numbers independent of the nature of the forcing (it is probably unlikely that these curves cross, though we cannot say more at present). However, given that observational estimates of |$\sigma ^{2}_{p_{\perp }} / \sigma ^{2}_{p}$| have not yet been made, it would be an interesting and important test to see if this ratio (∼0.7 ± 0.1) is realized in nature and whether there may be systematic variations in different environments. Note that, though mathematically well defined, the forcing scheme is somewhat idealized physically – assigning accelerations in Fourier space, which generate non-local accelerations in direct space. It is possible that the |$\sigma ^{2}_{p_{\perp }} / \sigma ^{2}_{p}$| ratios seen here may not necessarily be replicated in real molecular clouds, where large-scale driving sources could include supernovae or spiral shocks.

One could also consider the initially point-like energy injections from outflows in localized regions generating different |$\sigma ^{2}_{p_{\perp }} / \sigma ^{2}_{p}$| ratios. It may also be possible that |$\sigma ^{2}_{p_{\perp }} / \sigma ^{2}_{p}$| could evolve with time in decaying conditions (our simulations here are continually driven) or due to an increasing importance of self-gravity over time. Observational application of this method therefore may be more revealing and interesting than simply confirming the ‘default’ ∼0.7 ± 0.1 ratio seen in the current simulations. A framework for investigating the respective roles of solenoidal and compressive forcing in determining the normalized density variance has been recently presented by Konstandin et al. (2012b), and it would be interesting to compare (or combine) that method with ours.

The inclusion of self-gravity will be particularly interesting. Federrath et al. (2011) have shown that gravitational collapse produces a high fraction of longitudinal modes, which are later converted to solenoidal modes. It may therefore be possible to search for this signature in molecular clouds, after we have validated the method for the self-gravitating case. This is beyond the scope of the current work, but we note that the geometrical constraints, that lead to projection of solenoidal modes only, will not be affected by inclusion of self-gravity. It is true, however, that self-gravity will have a notable effect on the density PDF, and presumably on the degree of density–velocity correlation. A number of theoretical studies have looked at the generation of ∼power-law tails in density PDFs (e.g. Klessen 2000; Kritsuk, Norman & Wagner 2011; Federrath & Klessen 2013; Girichidis et al. 2014) and corresponding features have been seen observationally: early-time PDFs appear ∼lognormal in form, while more evolved clouds have more prevalent power-law tails in their density PDFs (e.g. Kainulainen et al. 2009; Schneider et al. 2013, 2014; Kainulainen, Federrath & Henning 2014).

5.4 Improvements

In the further testing of this model (in addition to the key question of anisotropy), there are two main areas that should be investigated. First, in terms of numerics, we need more information on the role of density–velocity correlations (Section 4) by measuring ϵ over a greater range of physical conditions, including the effects of self-gravity, as mentioned above, and a more detailed investigation of magnetic fields. In the analysis completed above, we have introduced the concept of ϵ and the correction factors g_mn and demonstrated their use in practice. Admittedly, there was a degree of internal ‘tuning’ of ϵ employed, though we did support the measured ϵ values with independent data, and the correction terms g₂₁ here are close enough to unity to be of not overwhelming concern. More important, perhaps, is the possible role that the statistical system presented in Section 4 may play in other contexts for understanding density–velocity correlations in turbulent media.

The second main avenue for improvement will be in the investigation of more mundane observational considerations, such as excitation and opacity. The analysis presented above was completed under the assumption of uniform excitation and in the optically thin limit, and therefore serves as a baseline for quantifying how variable excitation and finite opacity affect the method. Naively, this would make the measurements of |$\sigma ^{2}_{p_{\perp }} / \sigma ^{2}_{p}$| less reliable, and possibly biased, but may also have the advantageous effect of taming some of the field variability.

5.5 Other moments

Finally, we have based the method around the recovery of the fraction of momentum power, 〈ρ²v²〉, in transverse modes simply because the momentum field satisfies the observational requirement that it is projected into 2D unweighted by any other physical fields. One could also consider the fraction of velocity dispersion, 〈v²〉, in transverse modes (which is sensitive to the forcing mechanism – Federrath et al. 2011) or indeed the fraction of energy, 〈(1/2)ρv²〉, in transverse modes – either of which, arguably, more naturally spring to mind as relevant quantities to measure. The prospects for measuring either of the two above-named alternatives to momentum power are rather dim however, since neither field is accessible as a projected quantity observationally (except in the case of uniform density when all three definitions are equivalent).

In the general case where the density is variable, neither the velocity nor energy can be isolated for independent study. For example, note that the z-component of the momentum field, which is projected unweighted into 2D, is the product of the density, ρ, and z-velocity, v_z. Therefore, the crucially important (projected) power spectrum of W₁ involves the 3D convolution of the Fourier transforms of ρ and v_z, of which a single plane (k_z = 0) is available for analysis. To isolate the velocity contribution, one could in principle imagine deconvolving the density contribution for a simply structured density field (e.g. a 3D Gaussian or spatial power law), but for the highly complex and variable density fields encountered in the ISM, this appears futile. One possible way to proceed is to form an understanding of the relation between the three quantities (〈v²〉, 〈(1/2)ρv²〉, and 〈ρ²v²〉) by a means similar to the statistical system presented in the appendix – because the different quantities are just different ρ^q-weighted velocity dispersions. The stumbling block is simply that the Fourier space power distributions may be rather different, so that the deprojection factors become unreliable. In this case, an understanding of any scale-dependence of the statistical relationship between 〈v²〉, 〈(1/2)ρv²〉, and 〈ρ²v²〉 would be required, as well as an understanding of their transverse/longitudinal mode dependence. These considerations are beyond the scope of the current paper, but may be a worthwhile pursuit in future studies.

6 SUMMARY

In this paper, we have introduced an observational method for measuring the fraction of momentum power in solenoidal modes in a turbulent cloud, and confirmed its applicability using hydrodynamic numerical simulations. The method is best-suited to application in nearby molecular clouds for which high-sensitivity, high spatial dynamic range spectral line observations are available. The principal limitation of the method at present is its reliance on the assumption of isotropy. Further work is needed to examine the impact of anisotropy imposed by (e.g.) magnetic fields (BFP) or anisotropic driving of turbulence (Hansen et al. 2011). Isotropy aside, the main limiting factor in the accuracy of the model is variability in the physical fields (density, momentum). We have also introduced a statistical framework for describing density–velocity correlations in turbulent media that should be of relevance beyond its application here.

A big thanks to Daniel Price for allowing us use of the auxiliary numerical simulations, to Maria Cunningham for allowing us access to the Delta Quadrant Survey data and to Dave Acreman for much-needed help with Figs 1 and 2. CB is funded in part by the UK Science and Technology Facilities Council grant ST/J001627/1 (‘From Molecular Clouds to Exoplanets’) and the ERC grant ERC-2011-StG_20101014 (‘LOCALSTAR’), both held at the University of Exeter. CF acknowledges funding provided by the Australian Research Council under the Discovery Projects scheme (grant DP110102191). Supercomputing time at the Leibniz Rechenzentrum (project pr32lo) and at the Forschungszentrum Jülich (project hhd20) are gratefully acknowledged. The software used in this work was in part developed by the DOE-supported ASC/Alliance Center for Astrophysical Thermonuclear Flashes at the University of Chicago.

REFERENCES

Ballesteros-Paredes

J.

Mac Low

M.-M.

,

ApJ

,

2002

, vol.

570

pg.

734

Crossref

Search ADS

Ballesteros-Paredes

J.

D'Alessio

P.

Hartmann

L.

,

MNRAS

,

2012

, vol.

427

pg.

2562

Crossref

Search ADS

Brunt

C. M.

,

A&A

,

2010

, vol.

513

pg.

67

Crossref

Search ADS

Brunt

C. M.

Mac Low

M.-M.

,

ApJ

,

2004

, vol.

604

pg.

196

Crossref

Search ADS

Brunt

C. M.

Heyer

M. H.

Mac Low

M.-M.

,

A&A

,

2009

, vol.

504

pg.

883

Crossref

Search ADS

Brunt

C. M.

Federrath

C.

Price

D. J.

,

MNRAS

,

2010a

, vol.

403

pg.

1507

(BFP)

Crossref

Search ADS

Brunt

C. M.

Federrath

C.

Price

D. J.

,

MNRAS

,

2010b

, vol.

405

pg.

66

Crossref

Search ADS

Chapman

N. L.

Goldsmith

P. F.

Pineda

J. L.

Li

D.

Krčo

M.

,

ApJ

,

2011

, vol.

741

pg.

21

Crossref

Search ADS

Dubey

A.

et al.

Pogorelov

N. V.

Audit

E.

Zank

G. P.

,

ASP Conf. Ser. Vol. 385, Numerical Modelling of Plasma Flows

,

2008

San Francisco

Astron. Soc. Pac.

pg.

145

Google Scholar

Google Preview

OpenURL Placeholder Text

WorldCat

Elmegreen

B. G.

Scalo

J.

,

ARA&A

,

2004

, vol.

42

pg.

211

Crossref

Search ADS

Eswaran

V.

Pope

S. B.

,

Comput. Fluids

,

1988

, vol.

16

pg.

257

Crossref

Search ADS

Falgarone

E.

Lis

D. C.

Phillips

T. G.

Pouquet

A.

Porter

D. H.

Woodward

P. R.

,

ApJ

,

1994

, vol.

436

pg.

728

Crossref

Search ADS

Federrath

C.

,

MNRAS

,

2013

, vol.

436

pg.

1245

Crossref

Search ADS

Federrath

C.

Klessen

R. S.

,

ApJ

,

2013

, vol.

763

pg.

51

Crossref

Search ADS

Federrath

C.

Klessen

R. S.

Schmidt

W.

,

ApJ

,

2008

, vol.

688

pg.

79

Crossref

Search ADS

Federrath

C.

Klessen

R. S.

Schmidt

W.

,

ApJ

,

2009

, vol.

692

pg.

364

Crossref

Search ADS

Federrath

C.

Duval

J.

Klessen

R.

Schmidt

W.

Mac Low

M.-M.

,

A&A

,

2010

, vol.

512

pg.

81

Crossref

Search ADS

Federrath

C.

Chabrier

G.

Schober

J.

Banerjee

R.

Klessen

R. S.

Schleicher

D. R. G.

,

Phys. Rev. Lett.

,

2011

, vol.

107

pg.

114504

Crossref

Search ADS

PubMed

Fischera

J.

Dopita

M. A.

,

ApJ

,

2004

, vol.

611

pg.

919

Crossref

Search ADS

Fryxell

B.

et al. ,

ApJS

,

2000

, vol.

131

pg.

273

Crossref

Search ADS

Girichidis

P.

Konstandin

L.

Whitworth

A. P.

Klessen

R. S.

,

ApJ

,

2014

, vol.

781

pg.

91

Crossref

Search ADS

Hansen

C. E.

McKee

C. F.

Klein

R. I.

,

ApJ

,

2011

, vol.

738

pg.

88

Crossref

Search ADS

Heitsch

F.

Mac Low

M.-M.

Klessen

R. S.

,

ApJ

,

2001

, vol.

547

pg.

280

Crossref

Search ADS

Helmholtz

H.

,

J. Reine Angew. Math.

,

1858

, vol.

55

pg.

25

Crossref

Search ADS

Heyer

M. H.

Brunt

C. M.

,

MNRAS

,

2012

, vol.

420

pg.

1562

Crossref

Search ADS

Heyer

M. H.

Krawczyk

C.

Duval

J.

Jackson

J. M.

,

ApJ

,

2009

, vol.

699

pg.

1092

Crossref

Search ADS

Jackson

J. M.

et al. ,

ApJS

,

2006

, vol.

163

pg.

145

Crossref

Search ADS

Kainulainen

J.

Beuther

H.

Henning

T.

Plume

R.

,

A&A

,

2009

, vol.

508

pg.

35

Crossref

Search ADS

Kainulainen

J.

Federrath

C.

Henning

T.

,

Science

,

2014

, vol.

344

pg.

183

Crossref

Search ADS

PubMed

Kirk

H.

Pineda

J. E.

Johnstone

D.

Goodman

A.

,

ApJ

,

2010

, vol.

723

pg.

457

Crossref

Search ADS

Kleiner

S. C.

Dickman

R. L.

,

ApJ

,

1985

, vol.

295

pg.

466

Crossref

Search ADS

Klessen

R. S.

,

ApJ

,

2000

, vol.

535

pg.

869

Crossref

Search ADS

Klessen

R. S.

,

ApJ

,

2001

, vol.

556

pg.

837

Crossref

Search ADS

Klessen

R. S.

Heitsch

F.

Mac Low

M.-M.

,

ApJ

,

2000

, vol.

535

pg.

887

Crossref

Search ADS

Konstandin

L.

Federrath

C.

Klessen

R. S.

Schmidt

W.

,

J. Fluid Mech.

,

2012a

, vol.

692

pg.

183

Crossref

Search ADS

Konstandin

L.

Girichidis

P.

Federrath

C.

Klessen

R. S.

,

ApJ

,

2012b

, vol.

761

pg.

149

Crossref

Search ADS

Kritsuk

A. G.

Norman

M. L.

Wagner

R.

,

ApJ

,

2011

, vol.

727

pg.

L20

Crossref

Search ADS

Larson

R. B.

,

MNRAS

,

1981

, vol.

194

pg.

809

Crossref

Search ADS

Lo

N.

et al. ,

MNRAS

,

2009

, vol.

395

pg.

1021

Crossref

Search ADS

McKee

C. F.

Ostriker

E. C.

,

ARA&A

,

2007

, vol.

45

pg.

565

Crossref

Search ADS

Mac Low

M.-M.

Klessen

R. S.

,

Rev. Mod. Phys.

,

2004

, vol.

76

pg.

125

Crossref

Search ADS

McQuinn

K. B.

Simon

R.

Law

C. J.

Jackson

J. M.

Bania

T. M.

Clemens

D. P.

Heyer

M. H.

,

ApJ

,

2002

, vol.

576

pg.

274

Crossref

Search ADS

Micic

M.

Glover

S. C. O.

Federrath

C.

Klessen

R. S.

,

MNRAS

,

2012

, vol.

421

pg.

2531

Crossref

Search ADS

Miesch

M. S.

Bally

J.

,

ApJ

,

1994

, vol.

429

pg.

645

Crossref

Search ADS

Molina

F. Z.

Glover

S. C. O.

Federrath

C.

Klessen

R. S.

,

MNRAS

,

2012

, vol.

423

pg.

2680

Crossref

Search ADS

Ossenkopf

V.

Mac Low

M.-M.

,

A&A

,

2002

, vol.

390

pg.

307

Crossref

Search ADS

Ostriker

E. C.

Stone

J. M.

Gammie

C. F.

,

ApJ

,

2001

, vol.

546

pg.

980

Crossref

Search ADS

Padoan

P.

Nordlund

Å.

Jones

B. J. T.

,

MNRAS

,

1997

, vol.

288

pg.

145

Crossref

Search ADS

Price

D. J.

Federrath

C.

,

MNRAS

,

2010

, vol.

406

pg.

1659

Schmidt

W.

Niemeyer

J. C.

Hillebrandt

W.

,

A&A

,

2006

, vol.

450

pg.

265

Crossref

Search ADS

Schmidt

W.

Federrath

C.

Hupp

M.

Kern

S.

Niemeyer

J. C.

,

A&A

,

2009

, vol.

494

pg.

127

Crossref

Search ADS

Schneider

N.

Brooks

K.

,

PASA

,

2004

, vol.

21

pg.

290

Crossref

Search ADS

Schneider

N.

et al. ,

ApJ

,

2013

, vol.

766

pg.

L17

Crossref

Search ADS

Schneider

N.

et al. ,

A&A

,

2014

preprint (arXiv:1403.2966)

Simon

R.

Jackson

J. M.

Clemens

D. P.

Bania

T. M.

Heyer

M. H.

,

ApJ

,

2001

, vol.

551

pg.

747

Crossref

Search ADS

Solomon

P. M.

Rivolo

A. R.

Barrett

J.

Yahil

A.

,

ApJ

,

1987

, vol.

319

pg.

730

Crossref

Search ADS

Stutzki

J.

Bensch

F.

Heithausen

A.

Ossenkopf

V.

Zeilinsky

M.

,

A&A

,

1998

, vol.

336

pg.

697

Vázquez-Semadeni

E.

,

ApJ

,

1994

, vol.

423

pg.

681

Crossref

Search ADS

Vázquez-Semadeni

E.

Ballesteros-Paredes

J.

Rodriguez

L. F.

,

ApJ

,

1997

, vol.

474

pg.

292

Crossref

Search ADS

Waagan

K.

Federrath

C.

Klingenberg

C.

,

J. Comp. Phys.

,

2011

, vol.

230

pg.

3331

Crossref

Search ADS

Walsh

A. J.

Myers

P. C.

Burton

M. G.

,

ApJ

,

2004

, vol.

614

pg.

194

Crossref

Search ADS

Williams

J. P.

Blitz

L.

,

ApJ

,

1998

, vol.

494

pg.

657

Crossref

Search ADS

APPENDIX A: DENSITY–VELOCITY CORRELATIONS

In Section 2.3, we introduced a correction factor g₂₁, needed to convert a ρ-weighted velocity dispersion into a ρ²-weighted velocity dispersion. Here, we discuss the procedural steps necessary to evaluate the relationship between |$\langle \rho v^{2}_{z} \rangle /\langle \rho \rangle$| and |$\langle \rho ^{2} v^{2}_{z} \rangle /\langle \rho ^{2} \rangle$|⁠, in order to derive the statistical correction factor g₂₁ introduced in Section 2.3.

As we demonstrate below, if the density and velocity field are statistically uncorrelated (i.e. they have independent probability distributions), then |$\langle \rho v^{2}_{z} \rangle /\langle \rho \rangle = \langle \rho ^{2} v^{2}_{z} \rangle /\langle \rho ^{2} \rangle$|⁠, and therefore g₂₁ = 1. In general, this will not be the case, so in the following we introduce and test a simple method for quantifying the degree of density–velocity correlation and establish a simple prescription for converting between different ρ^q-weighted velocity dispersions. After theoretical development (Sections A1 and A2) and testing (Section A3), we examine existing observational constraints on these corrections (Section A4).

A1 Theoretical development

The velocity dispersion, |$\sigma ^{2}_{q}$|⁠, calculated using a ρ^q weighting can be written as

\begin{eqnarray} \sigma ^{2}_{q} &=& \frac{\frac{1}{V} \displaystyle \int _{V} {\rm{d}}V \; \rho ^{q} v^{2}}{\frac{1}{V} \displaystyle \int _{V} {\rm{d}}V \; \rho ^{q}} = \frac{\langle \rho ^{q} v^{2}\rangle }{\langle \rho ^{q} \rangle } \nonumber \\ &=& \frac{\frac{1}{V} \displaystyle \int _{V} {\rm{d}}V \; \xi ^{q} v^{2}}{\frac{1}{V} \displaystyle \int _{V} {\rm{d}}V \; \xi ^{q}} = \frac{\langle \xi ^{q} v^{2}\rangle }{\langle \xi ^{q} \rangle }, \end{eqnarray}

(A1)

where ρ₀ = 〈ρ〉, ξ = ρ/ρ₀, V is the volume containing the fields and angle brackets denote spatial averages. In what follows, we use the more convenient variable ξ to perform calculations.

The volumetric integration can be replaced by integrals over ξ and v as follows:

\begin{equation} \sigma ^{2}_{q} = \frac{\displaystyle \int _{0}^{\infty }\displaystyle \int _{-\infty }^{\infty } {\rm{d}}\xi \; {\rm{d}}v\; P_{\xi }(\xi ) P_{v}(v) \; \xi ^{q} v^{2}}{\displaystyle \int _{0}^{\infty } {\rm{d}}\xi \; P_{\xi }(\xi ) \; \xi ^{q}}, \end{equation}

(A2)

where P_ξ(ξ) and P_v(v) are the PDFs of ξ and v, respectively. To account for correlations between density and velocity, we consider the case now where P_v(v) is an implicit function of ξ, so that completing the v integral leads to

\begin{equation} \sigma ^{2}_{q} = \frac{\displaystyle \int _{0}^{\infty } {\rm{d}}\xi \; P_{\xi }(\xi ) \; \xi ^{q} \sigma ^{2}_{v}(\xi ) }{\displaystyle \int _{0}^{\infty } {\rm{d}}\xi \; P_{\xi }(\xi ) \; \xi ^{q}}, \end{equation}

(A3)

where |$\sigma ^{2}_{v}(\xi )$| is a density-dependent velocity dispersion. We propose a simple form for this dispersion as follows:

\begin{equation} \sigma ^{2}_{v}(\xi ) = h(\xi ) \sigma ^{2}_{00}, \end{equation}

(A4)

where |$\sigma ^{2}_{00}$| is a constant, and further propose that

\begin{equation} h(\xi ) = \xi ^{-\epsilon }, \end{equation}

(A5)

with the expectation (but not requirement) that ϵ is a small positive constant, so that higher densities are associated with smaller velocity dispersions. Equation (A3) then becomes

\begin{equation} \sigma ^{2}_{q} = \frac{\displaystyle \int _{0}^{\infty } {\rm{d}}\xi \; P_{\xi }(\xi ) \; \xi ^{q-\epsilon } \sigma ^{2}_{00} }{\displaystyle \int _{0}^{\infty } {\rm{d}}\xi \; P_{\xi }(\xi ) \; \xi ^{q}}, \end{equation}

(A6)

with the immediate result:

\begin{equation} \sigma ^{2}_{q} = \frac{\langle \xi ^{q-\epsilon } \rangle \sigma ^{2}_{0} }{ \langle \xi ^{q} \rangle \langle \xi ^{-\epsilon } \rangle }, \end{equation}

(A7)

where we have written the unweighted (q = 0) velocity dispersion, |$\sigma ^{2}_{0}$|⁠, as

\begin{equation} \sigma ^{2}_{0} = \sigma ^{2}_{00} \langle \xi ^{-\epsilon } \rangle . \end{equation}

(A8)

Note that if density and velocity are statistically uncorrelated (ϵ = 0), then equation (A7) just gives |$\sigma ^{2}_{q} = \sigma ^{2}_{0}$| for all q.

Using equations (A1) and (A7) without the normalizing denominator, we may also write

\begin{equation} \langle \xi ^{q} v^{2} \rangle = \frac{\langle \xi ^{q-\epsilon } \rangle \sigma ^{2}_{0} }{ \langle \xi ^{-\epsilon } \rangle }, \end{equation}

(A9)

i.e. that the spectrum of moments 〈ξ^qv²〉(q) is a scaled, shifted version of the spectrum of moments 〈ξ^q〉(q).

A2 Formulation for lognormal density PDFs

Knowledge of 〈ξ^qv²〉/〈ξ^q〉 would be very useful, but it in turn requires knowledge of 〈ξ^q〉(q) which itself is inaccessible observationally. However, progress can be made by employing an analytic form for the PDF of normalized density, P_ξ(ξ). For isothermal turbulence, a lognormal PDF for ξ is a reasonable approximation (Vázquez-Semadeni 1994; Padoan, Nordlund & Jones 1997; Federrath et al. 2008; Kainulainen et al. 2009), for which

\begin{equation} \langle \xi ^{q} \rangle = \exp \left[ q \langle \ln (\xi ) \rangle + \frac{1}{2} q^{2} \sigma ^{2}_{\ln (\xi )} \right]\! . \end{equation}

(A10)

The normalization of the field to 〈ξ〉 = 1 requires that |$\langle \ln (\xi ) \rangle = -\frac{1}{2} \sigma ^{2}_{\ln (\xi )}$|⁠, so that

\begin{equation} \langle \xi ^{q} \rangle = \exp \left[ \frac{1}{2} \sigma ^{2}_{\ln (\xi )} (q^{2} - q) \right] = \langle \xi ^{2} \rangle ^{\frac{1}{2}(q^{2} - q)}, \end{equation}

(A11)

where we have made use of |$\sigma ^{2}_{\ln (\xi )} = \ln (1 + \sigma ^{2}_{\xi }) = \ln (\langle \xi ^{2} \rangle )$| in the last step. The predicted form for 〈ξ^q〉 (in the case of a lognormal density PDF) agrees well with the moment spectra shown in Fig. A1: i.e. ln (〈ξ^q〉) is parabolic against q, equal to unity at q = 0 and q = 1, and reaches a minimum at q = 0.5.

The spectrum of moments for 〈ξq〉 (dots) and 〈ξqv2〉 (triangles) versus q. The curve (described in the text) is fitted to 〈ξq〉 and shifted to best match 〈ξqv2〉 in order to determine ϵ.

Figure A1.

The spectrum of moments for 〈ξ^q〉 (dots) and 〈ξ^qv²〉 (triangles) versus q. The curve (described in the text) is fitted to 〈ξ^q〉 and shifted to best match 〈ξ^qv²〉 in order to determine ϵ.

Open in new tab Download slide

Inserting this result into equation (A9) then gives

\begin{equation} \langle \xi ^{q} v^{2} \rangle = \langle \xi ^{2} \rangle ^{\frac{1}{2}(q^{2} - q - 2q\epsilon )} \sigma ^{2}_{0}, \end{equation}

(A12)

and equation (A7) becomes

\begin{equation} \sigma ^{2}_{q} = \frac{\langle \xi ^{q} v^{2} \rangle }{\langle \xi ^{q} \rangle } = \langle \xi ^{2} \rangle ^{-q\epsilon } \sigma ^{2}_{0}. \end{equation}

(A13)

We define velocity dispersion ratios, g_mn, as follows:

\begin{equation} g_{mn} = \frac{\sigma ^{2}_{m}}{\sigma ^{2}_{n}} = \frac{\langle \rho ^{m} v^{2} \rangle /\langle \rho ^{m} \rangle }{\langle \rho ^{n} v^{2} \rangle /\langle \rho ^{n} \rangle } = \langle \xi ^{2} \rangle ^{-(m-n)\epsilon }. \end{equation}

(A14)

In the next section, we will make use of the following ratio:

\begin{equation} g_{21} = \frac{\sigma ^{2}_{2}}{\sigma ^{2}_{1}} = \frac{\langle \rho ^{2} v^{2} \rangle /\langle \rho ^{2} \rangle }{\langle \rho v^{2} \rangle /\langle \rho \rangle } = \langle \xi ^{2} \rangle ^{-\epsilon }. \end{equation}

(A15)

Note that if the density field is uniform (〈ξ²〉 = 1), then g₂₁ = 1. In general, the assumption ϵ = 0 requires that there be no statistical correlation between density and velocity. To apply this correction observationally, a measurement of 〈ξ²〉 is required. This is can be done, since the 3D normalized density dispersion, |$\sigma ^{2}_{\xi }$|⁠, and therefore, |$\langle \xi ^{2} \rangle = 1 + \sigma ^{2}_{\xi }$| can be estimated by the BFP method.

A3 Numerical testing

The applicability of equation (A9) to the numerical density and velocity fields can be tested by calculating and comparing 〈ξ^qv²〉(q) and 〈ξ^q〉(q). From the simulated data, moments between q = 0 and q = 2 in steps of 0.05 were calculated in the zero momentum frame, using the full 3D velocity field. To ensure that the resulting values of ϵ are not simply ‘tuned’ to the numerical simulations analysed in this paper, we use a much larger sample of density and velocity fields by including the simulations from BFP. These include both hydrodynamic simulations (using SPH) and MHD simulations (using grid calculations) – see Section 3 and BFP for more details. We only consider supersonic fields in this analysis, since the influence of density fluctuations in subsonic fields is minimal – we will return to this point at the end of this section.

After calculation of the moments, a sixth-order polynomial in q (⁠|$\sum _{n=0}^{6} c_{n} q^{n}$|⁠) was fitted to ln (〈ξ^q〉). Following this, the same polynomial with the same (fixed) c_n was fitted to ln (〈ξ^qv²〉) with additional offset (ϵ) and scaling factor (A) – i.e. |$\ln (A) + \sum _{n=0}^{6} c_{n} (q - \epsilon )^{n}$|⁠.

Fig. A1 shows representative fits to data taken from the latest snapshot of the grid simulations for both solenoidal and compressive forcing (the other simulations yield similar results). It is clear that the proposed form h(ξ) = ξ^−ϵ results in a good representation of the moment spectrum 〈ξ^qv²〉(q). It is also apparent from the plots in Fig. A1 that ϵ (i.e. the horizontal shift of the 〈ξ^qv²〉 moment spectrum relative to the 〈ξ^q〉 moment spectrum) decreases with increasing rms Mach number. To quantify this further, the fitted values of ϵ for all hydrodynamic data are shown in Fig. A2, plotted versus the measured 3D density-weighted Mach number, |${\mathcal {M}}_{1} = \sigma _{1}/c_{\rm s}$|⁠, where c_s is the sound speed. A power-law relation is seen, represented by the fitted line, which is

\begin{equation} \epsilon _{\infty } = (0.38\,{\pm }\, 0.05) {\mathcal {M}}_{1}^{-0.72\, {\pm }\, 0.05}, \end{equation}

(A16)

where we have used the subscript ∞ to denote that the fit was obtained in the hydrodynamic limit where the Alfvénic Mach number |${\mathcal {M}}_{\rm A} \rightarrow \infty$|⁠.

$Plot of ϵ versus 3D density-weighted Mach number, ${\mathcal {M}}_{1}$, derived from the hydrodynamic grid simulations (dots: solenoidal forcing; open squares: compressive forcing; triangles: hydrodynamic (SPH) simulations with solenoidal forcing from BFP.$

Figure A2.

Plot of ϵ versus 3D density-weighted Mach number, |${\mathcal {M}}_{1}$|⁠, derived from the hydrodynamic grid simulations (dots: solenoidal forcing; open squares: compressive forcing; triangles: hydrodynamic (SPH) simulations with solenoidal forcing from BFP.

Open in new tab Download slide

The values of ϵ obtained from the MHD simulations are shown in Fig. A3, where they are plotted versus |${\mathcal {M}}_{1}$|⁠. Relative to the hydrodynamic results, represented by the straight line (ϵ_∞), the MHD fields tend to give lower values of ϵ at the lower Mach numbers. We find that in the low |${\mathcal {M}_{1}}$| regime, ϵ appears to decrease with decreasing Alfvénic Mach Number, |${\mathcal {M}}_{\rm A}$|⁠, as shown in Fig. A4, where we have plotted ϵ/ϵ_∞ versus |${\mathcal {M}}_{\rm A}$|⁠. However, an |${\mathcal {M}}_{\rm A}$|-dependence (if any) is less clear at higher |${\mathcal {M}}_{1}$|⁠. In the analysis below, we will adopt ϵ_∞, obtained in the hydrodynamic limit, as the basis for deriving correction factors (g₂₁).

$Plot of ϵ versus 3D density-weighted Mach number, ${\mathcal {M}}_{1}$, derived from the MHD simulations (BFP). For reference, we have shown the fitted line from Fig. A2.$

Figure A3.

Plot of ϵ versus 3D density-weighted Mach number, |${\mathcal {M}}_{1}$|⁠, derived from the MHD simulations (BFP). For reference, we have shown the fitted line from Fig. A2.

Open in new tab Download slide

$Plot of ϵ (relative to the fitted $\epsilon _{\infty }({\mathcal {M}}_{1})$ relation in the hydrodynamic limit) versus 3D Alfvénic Mach number. Symbols denote the rms sonic Mach number in the simulation (open squares: ${\mathcal {M}}_{1} = 20$; dots: ${\mathcal {M}}_{1} = 10$; triangles: ${\mathcal {M}}_{1} = 4$; asterisks: ${\mathcal {M}}_{1} = 2$).$

Figure A4.

Plot of ϵ (relative to the fitted |$\epsilon _{\infty }({\mathcal {M}}_{1})$| relation in the hydrodynamic limit) versus 3D Alfvénic Mach number. Symbols denote the rms sonic Mach number in the simulation (open squares: |${\mathcal {M}}_{1} = 20$|⁠; dots: |${\mathcal {M}}_{1} = 10$|⁠; triangles: |${\mathcal {M}}_{1} = 4$|⁠; asterisks: |${\mathcal {M}}_{1} = 2$|⁠).

Open in new tab Download slide

For a given 〈ξ²〉, a smaller ϵ leads to a correction factor (g₂₁) closer to unity. In fact, at fixed sonic Mach number, numerical simulations suggest that 〈ξ²〉 should be closer to unity for strongly magnetized turbulence than for hydrodynamic turbulence (Molina et al. 2012). Together with lower ϵ, this would predict that g₂₁ should be closer to unity for strongly magnetized turbulence than for the hydrodynamic case.

Finally, as mentioned previously, we have only considered supersonic fields in the above analysis. Since the velocity dispersion ratios (equation A14) depend on 〈ξ^q〉 to the power −(m − n)ϵ, fields for which 〈ξ^q〉 ≈ 1 (i.e. uniform or weakly varying density fields) require large values of ϵ for small corrections. We find that measurement of ϵ in this regime (by the method given above) is rather unstable, with ϵ increasing strongly, with large scatter, as |$\langle \xi ^{q} \rangle \longrightarrow 1$|⁠. In the analysis in Section 4, we will nevertheless investigate values of the correction factor g₂₁ extrapolated into the subsonic regime.

A4 Constraints from observational data

In the preceding analysis, we made the assumption that the velocity dispersion decreases with the density at which it is measured (i.e. equations A4 and A5), and explored the consequences on velocity dispersions calculated with a ρ^q-weight. To do this, privileged access to the density and velocity fields in 3D is required, so this is not possible to do observationally by exactly the same method. An alternative means of constraining the effect of such density–velocity correlations is to interpret equations (A4) and (A5) literally and examine velocity dispersions measured observationally in different density regimes.

A4.1 Constraints from Larson's relations

The Larson (1981) relations between velocity dispersion (or linewidth) and cloud size (σ_v ∝ L^a) and between density and cloud size (ρ ∝ L^−b) can provide a very crude measure of ϵ in molecular clouds. Combining the two relations gives a measure of ϵ ≈ a/b, though with some major caveats. The original a = 0.38 and b = 1.1 derived by Larson (1981) give ϵ = 0.35, while values of a ≈ 0.5 and b = 1 (assuming Virial equilibrium) give ϵ ≈ 0.5 from Solomon et al. (1987). Large-scale CO surveys use size, linewidth and mean density measured in distinct clouds rather than probing the density-dependence of velocity dispersion internal to individual clouds. Extending such cloud/clump-based analyses to ‘subcloud’ scales are highly questionable (e.g. Ballesteros-Paredes & Mac Low 2002; Schneider & Brooks 2004), but nevertheless tend to yield values of a and b roughly in accord with the large-scale values, though in the presence of significant scatter. Using the analysis of Simon et al. (2001) on four inner Galaxy clouds (in the Solomon et al. 1987 survey region), we find values of ϵ between 0.17 and 0.39 (the lower values being mostly due to shallower linewidth–size relations than that found by Solomon et al. 1987).

Criticisms have been levelled at the density-size relation as being a consequence of limited dynamic range in cloud surface density (e.g. Vázquez-Semadeni, Ballesteros-Paredes & Rodriguez 1997). Using the Solomon et al. (1987) cloud sample, Heyer et al. (2009) found that, over the limited surface density range available, velocity dispersions rose with surface density as σ_v ∝ Σ^1/2 at fixed cloud size. If we make the reasonable assumption that higher surface density indicates higher volume density, then a/b should provide an upper limit to ϵ.

The ϵ values derived here by Larson's relations are quoted for reference, and should be compared to the better-motivated (and notably lower) values derived below, using density-selective tracers.

A4.2 Constraints from density-selective tracers

A better-motivated idea than using Larson's relations is to compare velocity dispersions of trace molecules that are excited in different density regimes in the same cloud. We expect that high-density tracers should have smaller velocity dispersions than low-density tracers, if their spectral lines are averaged over the same (large) volume. A rough estimate of ϵ may be arrived at by assuming that spectral line emission from a given trace molecule is dominated by contributions from material near the molecule's critical density. With this assumption, if tracer A has critical density n_{c, A} and velocity dispersion |$\sigma ^{2}_{v,A}$| and tracer B has critical density n_{c, B} and velocity dispersion |$\sigma ^{2}_{v,B}$|⁠, then equations (A4) and (A5) predict that

\begin{equation} \frac{\sigma ^{2}_{v,A}}{\sigma ^{2}_{v,B}} \approx \left(\frac{n_{c,A}}{n_{c,B}}\right)^{-\epsilon }. \end{equation}

(A17)

From this, ϵ may be estimated via

\begin{equation} \epsilon \approx - \log \left(\frac{\sigma ^{2}_{v,A}}{\sigma ^{2}_{v,B}}\right)\bigg / \log \left(\frac{n_{c,A}}{n_{c,B}}\right). \end{equation}

(A18)

In reality, the tracers will sample a range of densities above their effective critical densities. A better motivation for equations (A17) and (A18) can be arrived at by more carefully considering such a system. If we restrict the analysis to two tracers whose effective critical densities lie on the positive tail of the density PDF and assume that in this regime, self-gravity will push the PDF into a power-law form (P_ξ(ξ) ∝ ξ^−α), then we can calculate velocity dispersions via equation (A6) as

\begin{equation} \sigma ^{2}_{v,A} = \frac{\displaystyle \int _{\xi _{c,A}}^{\infty } {\rm{d}}\xi \; \xi ^{-\alpha } \; \xi ^{1-\epsilon } \sigma ^{2}_{00} }{\displaystyle \int _{\xi _{c,A}}^{\infty } {\rm{d}}\xi \; \xi ^{-\alpha } \; \xi } = \left(\frac{2 - \alpha }{2 - \epsilon - \alpha }\right) \xi _{c,A}^{- \epsilon } \sigma ^{2}_{00}, \end{equation}

(A19)

(and similarly for tracer B) and where we have taken the lower limit of the integral as ξ_{c, A} = n_{c, A}/n₀ = ρ_{c, A}/ρ₀ and assumed that α > q + 1 = 2 so that the integrals converge. Note that the velocity dispersion is calculated with a ρ-weight (q = 1).

Taking the velocity dispersion ratio of the two tracers, we find

\begin{equation} \frac{\sigma ^{2}_{v,A}}{\sigma ^{2}_{v,B}} = \left(\frac{\xi _{c,A}}{\xi _{c,B}}\right)^{- \epsilon }, \end{equation}

(A20)

which is exactly the same as equation (A17). If not a power-law PDF, a sufficiently steeply falling PDF (lognormal, exponential) will give similar results, as the dispersions are heavily influenced by densities at the lower threshold (see e.g. Ballesteros-Paredes, D'Alessio & Hartmann 2012, for examination of a mathematically equivalent system).

A4.3 Observational estimates of ϵ

In this section, we will estimate a few contrasting values of ϵ from observational data.

McQuinn et al. (2002) examined CS (J = 2–1) and ¹³CO (J = 1–0) emission in the inner Galaxy observed as part of the Galactic Ring Survey (Jackson et al. 2006) and found no significant difference between CS (J = 2–1) and ¹³CO (J = 1–0) velocity dispersions, as evidenced by roughly constant brightness temperature ratios for composite spectra averaged over large clouds (many parsec scale). The critical densities of CS (J = 2–1) and ¹³CO (J = 1–0) are ∼10³ cm⁻³ and ∼5 × 10⁵ cm⁻³ respectively. Equation (A18) in this case finds ϵ is ‘very small’ (the g₂₁ correction factor in our lognormal model above would therefore simply be ∼ unity – i.e. no correction). However, McQuinn et al. (2002) ultimately concluded that subthermal excitation was a probable factor in the line excitation (especially for CS) so that the effective critical densities would be lower than nominal.

To examine this more closely, we have used multitracer spectral line data, spanning a larger range in critical density, from the Delta Quadrant Survey (see Lo et al. 2009). Fig. A5 shows three regions within the survey where spectral lines are single-component and allow easy fitting of Gaussian functions to estimate their dispersions. We have averaged the intensities over the boxes shown and determined velocity dispersions for each of the four transitions (C¹⁸O 1–0, CS 2–1, HCO⁺ 1–0 and HNC 1–0). The spectra, with fitted Gaussians, are also shown in Fig. A5, along with the variation of velocity dispersion with critical density of the tracer. From this, we determine ϵ ≈ 0.

Three subregions in the Delta Quadrant Survey, in which the intensities of C18O 1–0, CS 2–1, HCO+ 1–0, and HNC 1–0 spectral lines have been averaged and fitted with Gaussians to determine their velocity dispersions – shown in the centre panel, with spectra scaled to their peak and offset for clarity. In the right-hand panel, we plot log velocity dispersion versus log critical density to determine ϵ as the slope of a fitted line for each region.

Figure A5.

Three subregions in the Delta Quadrant Survey, in which the intensities of C¹⁸O 1–0, CS 2–1, HCO⁺ 1–0, and HNC 1–0 spectral lines have been averaged and fitted with Gaussians to determine their velocity dispersions – shown in the centre panel, with spectra scaled to their peak and offset for clarity. In the right-hand panel, we plot log velocity dispersion versus log critical density to determine ϵ as the slope of a fitted line for each region.

Open in new tab Download slide

Williams & Blitz (1998) examined CS (J = 2–1) and ¹³CO (J = 1–0) and (J = 3–2) emission from a star-forming cloud (the Rosette nebula) and a non-star-forming cloud (G216, ‘Maddalena's Cloud’). They found small differences in clump internal velocity dispersion in the (J = 3–2) and (J = 1–0) ¹³CO lines: slightly broader 3–2 lines in the Rosette, relative to 1–0 (implying ϵ < 0 in this case, since the critical density of the 3–2 line is ∼10 times that of the 1–0 line); slightly narrower lines are found for 3–2 in G216, relative to 1–0. They attribute this behaviour to local feedback effects close to star-forming sources in the Rosette, which have no counterparts in the less active G216. Williams & Blitz (1998) do not report values for the linewidth ratios, but they can be estimated from their fig. 19; we will take |$\sigma ^{2}_{v,3-2}/\sigma ^{2}_{v,1-0}$| ≈ 2 and |$\sigma ^{2}_{v,3-2}/\sigma ^{2}_{v,1-0}$| ≈ 0.5 as representative for the Rosette and G216, respectively. These lead, assuming n_{c, 3−2}/n_{c, 1−0} ≈ 10, to ϵ_Rosette ≈ −0.3 and ϵ_G216 ≈ +0.3. These values of ϵ are notably larger in magnitude than those found for our numerical simulations (and in the case of the Rosette involve a sign change) but may not be representative values for globally determined ϵ as we require. A negative value of ϵ probably cannot be maintained over all densities (over all space), though one could potentially imagine high-density clumps moving through a relatively static low-density substrate as a possible configuration for this. The negative value of ϵ for the Rosette more likely comes instead from (e.g.) outflows injecting energy locally in a characteristic density regime near the J=3–2 critical density. Such behaviour could cause problems for our simple power law characterization in equation (A5).

McQuinn et al. (2002) made large-scale averages of line profiles to produce their CS–¹³CO comparisons, which is more closely matched to our requirements. In contrast, Williams & Blitz (1998) examined linewidths from targeted clumps, where the role of density–velocity correlations are likely to be most emphasized. With this provison in mind, we now examine the CS/¹³CO (J = 1–0) velocity dispersion ratios found by Williams & Blitz (1998). These are |$\sigma ^{2}_{v,{\rm{CS}}}/\sigma ^{2}_{v,{\rm{13CO}}}$| ∼ 0.7 (Rosette) and |$\sigma ^{2}_{v,{\rm{CS}}}/\sigma ^{2}_{v,{\rm{13CO}}}$| ∼ 0.5 (G216). Taking a critical density ratio of 500, we then find ϵ_Rosette ≈ +0.06 and ϵ_G216 ≈ +0.11. With a larger baseline in critical density (and comparable hν/k values for the transitions), these CS/¹³CO-derived values of ϵ are a more reliable measure than the 3–2/1–0-derived values, and are more in line with our numerically derived ϵ-values.

We can estimate values of ϵ using C¹⁸O (J = 1–0) and N₂H⁺ (J = 1–0) data from the Perseus molecular cloud reported by Kirk et al. (2010), who averaged spectra over spatially extended regions. The relative linewidths of C¹⁸O and N₂H⁺ vary amongst the targeted regions (see their figs 7–12). We estimate that |$\sigma ^{2}_{v,{\rm{N2H+}}}/\sigma ^{2}_{v,{\rm{C18O}}}$| varies between ∼0.2 and ∼1 from their graphs, and will take an approximate value of |$\sigma ^{2}_{v,{\rm{N2H+}}}/\sigma ^{2}_{v,{\rm{C18O}}}$| ≈ 0.5 for the ensemble. Assuming the ratio of critical densities is ∼1000, we find ϵ ∼+0.1.

Walsh, Myers & Burton (2004) reported ¹³CO, C¹⁸O and N₂H⁺ (all 1–0) linewidths from a sample of nearby cores. These are targeted, single-point spectra towards separated regions, but taking an average over all spectra, we find that |$\langle\sigma ^{2}_{v,{\rm{N2H+}}}/\sigma ^{2}_{v,{\rm{C18O}}}\rangle$| = 0.66 ± 0.37 and |$\langle\sigma ^{2}_{v,{\rm{N2H+}}}/\sigma ^{2}_{v,{\rm{13CO}}}\rangle$| = 0.29 ± 0.16. From these, we derive ϵ values of +0.06 and +0.18, assuming a critical density ratio of 1000 in both cases. Within the uncertainties the dispersion ratios can be reconciled, though it is clear that ¹³CO linewidths are broader than C¹⁸O in most circumstances. The effective critical density for ¹³CO and C¹⁸O is likely to be different, due to abundance differences (i.e. ¹³CO is more abundant in the lower density regions of the cloud than C¹⁸O). In a medium where the density PDF rises sharply towards lower densities, an abundant molecule's emission is in principle subject to some, potentially significant, contribution from subthermally excited regions where radiative trapping is important, thereby lowering the critical density from its nominal value. However, plausible variations (up to an order of magnitude difference) in their effective critical densities cannot equalize the different ϵ values as ϵ is only logarithmically sensitive to the assumed critical density ratio. We propose that the difference in ϵ can potentially be explained by path-length differences due to abundance effects: a velocity dispersion ratio of |$\langle\sigma ^{2}_{v,{\rm{13CO}}}/\sigma ^{2}_{v,{\rm{C18O}}}\rangle$| ≈ 0.66/.29 ≈ 2.28 can be explained by a path-length ratio of 2.28 in a medium where |$\sigma ^{2}_{v} \propto L$|⁠.

The small, positive nature of ϵ in most of the estimates above is in line with our initial expectation set out in Section 4.1. Two factors can contribute to this. First, the velocity dispersion over some spatial scale L in high density may be physically lower than the velocity dispersion over the same spatial scale L in low-density gas. Secondly, in a turbulent medium, the velocity dispersion increases with spatial scale, and the lower density structures are necessarily more spatially extended than the higher density structures. For molecules with the same (nominal) critical density, the more abundant molecule will have a lower effective critical density, and therefore be more spatially extended and have a higher velocity dispersion. (In principle, the first of these two factors may be reversed (i.e. higher velocity dispersion in denser gas at fixed spatial scale) and still yield a positive ϵ as long as the second factor dominates.)

If a targeted, single-point measurement is made towards an atypical position (e.g. a core) within a larger medium, as was done for some of the observations reported above, this may not provide a reliable measure of the density-dependence of velocity dispersion in the medium as a whole. Instead, if a velocity dispersion measurement is made using a high-density tracer averaged over sufficiently large scales, then it will sample many density enhancements that are spatially distributed, and therefore result in a larger overall velocity dispersion (more comparable to the extended medium seen by a lower density tracer).

Though the ¹³CO (J = 3–2) and (J = 1–0) results above raise some questions, the observationally estimated values of ϵ over large critical density baselines are in reasonable accord with those found in our numerical simulations. We will take 0.05 ≲ ϵ ≲ 0.3 as defining the probable range of ϵ from the above calculations, with values less than ∼0.1 being favoured (i.e. as derived in cases where large-scale spatial averaging is conducted, and when the critical density span is larger). The 3D ρ-weighted Mach numbers in the Rosette, G216 and Perseus are ∼10–20, judging by the ¹³CO (J = 1–0) linewidths, for kinetic temperatures of ∼10–20 K. The observational estimates compare reasonably well with the numerical results (ϵ versus |${\mathcal {M}}_{1}$|⁠) in Fig. A2.

Download all slides

Month:	Total Views:
January 2017	4
March 2017	4
April 2017	2
July 2017	3
August 2017	2
September 2017	1
October 2017	5
November 2017	4
December 2017	6
January 2018	5
February 2018	8
March 2018	14
April 2018	8
May 2018	8
June 2018	7
July 2018	9
August 2018	10
September 2018	4
October 2018	16
November 2018	12
December 2018	10
January 2019	9
February 2019	18
March 2019	9
April 2019	27
May 2019	19
June 2019	7
July 2019	11
August 2019	15
September 2019	7
October 2019	11
November 2019	6
December 2019	7
January 2020	10
February 2020	6
March 2020	7
April 2020	13
May 2020	2
June 2020	5
July 2020	2
August 2020	3
September 2020	6
November 2020	3
December 2020	17
January 2021	10
February 2021	9
March 2021	17
April 2021	9
June 2021	4
July 2021	3
August 2021	3
September 2021	5
October 2021	10
November 2021	3
December 2021	4
January 2022	6
February 2022	3
March 2022	3
April 2022	1
May 2022	9
June 2022	8
July 2022	11
August 2022	6
September 2022	13
October 2022	3
November 2022	3
January 2023	5
February 2023	2
March 2023	3
April 2023	10
May 2023	3
June 2023	3
July 2023	4
August 2023	7
September 2023	4
October 2023	22
November 2023	8
December 2023	6
January 2024	8
February 2024	5
March 2024	13
April 2024	13
May 2024	27
June 2024	6
July 2024	6
August 2024	10
September 2024	6
October 2024	6
November 2024	5
December 2024	8
January 2025	8
February 2025	10
March 2025	12
April 2025	7
May 2025	11

Article Contents

An observational method to measure the relative fractions of solenoidal and compressible modes in interstellar clouds

Abstract

1 INTRODUCTION

2 SOLENOIDAL AND COMPRESSIBLE MODES

2.1 General considerations

2.2 Projection from 3D to 2D

2.3 Observational considerations

2.4 Summary of observational measurements required

3 NUMERICAL SIMULATIONS

4 APPLICATION TO HYDRODYNAMIC SIMULATIONS

4.1 Normalized density variance

4.2 Dispersion in z-axis momentum

4.3 Projected solenoidal fractions

4.4 Deprojection

5 DISCUSSION

5.1 Anisotropies and discontinuities

5.2 Variability

5.3 Variations in solenoidal/compressive fractions

5.4 Improvements

5.5 Other moments

6 SUMMARY

REFERENCES

APPENDIX A: DENSITY–VELOCITY CORRELATIONS

A1 Theoretical development

A2 Formulation for lognormal density PDFs

A3 Numerical testing

A4 Constraints from observational data

A4.1 Constraints from Larson's relations

A4.2 Constraints from density-selective tracers

A4.3 Observational estimates of ϵ

Citations

Views

Altmetric

Email alerts

Astrophysics Data System

Citing articles via

Latest

Most Read

Most Cited

Article Contents

An observational method to measure the relative fractions of solenoidal and compressible modes in interstellar clouds

Abstract

1 INTRODUCTION

2 SOLENOIDAL AND COMPRESSIBLE MODES

2.1 General considerations

2.2 Projection from 3D to 2D

2.3 Observational considerations

2.4 Summary of observational measurements required

3 NUMERICAL SIMULATIONS

4 APPLICATION TO HYDRODYNAMIC SIMULATIONS

4.1 Normalized density variance

4.2 Dispersion in z-axis momentum

4.3 Projected solenoidal fractions

4.4 Deprojection

5 DISCUSSION

5.1 Anisotropies and discontinuities

5.2 Variability

5.3 Variations in solenoidal/compressive fractions

5.4 Improvements

5.5 Other moments

6 SUMMARY

REFERENCES

APPENDIX A: DENSITY–VELOCITY CORRELATIONS

A1 Theoretical development

A2 Formulation for lognormal density PDFs

A3 Numerical testing

A4 Constraints from observational data

A4.1 Constraints from Larson's relations

A4.2 Constraints from density-selective tracers

A4.3 Observational estimates of ϵ

Citations

Views

Altmetric

Email alerts

Astrophysics Data System

Citing articles via

Latest

Most Read

Most Cited

This Feature Is Available To Subscribers Only