EDP Sciences
Free Access
Issue
A&A
Volume 516, June-July 2010
Article Number A25
Number of page(s) 14
Section Interstellar and circumstellar matter
DOI https://doi.org/10.1051/0004-6361/200913904
Published online 22 June 2010
A&A 516, A25 (2010)

Numerical and semi-analytic core mass distributions in supersonic isothermal turbulence

W. Schmidt1,2 - S. A. W. Kern2,3 - C. Federrath3,4 - R. S. Klessen4

1 - Institut für Astrophysik, Universität Göttingen, Friedrich-Hund-Platz 1, 37077 Göttingen, Germany
2 - Lehrstuhl für Astronomie, Institut für Theoretische Physik und Astrophysik, Universität Würzburg, Am Hubland, 97074 Würzburg, Germany
3 - Max-Planck-Institut für Astronomie, Königstuhl 17, 69117 Heidelberg, Germany
4 - Zentrum für Astronomie, Institut für Theoretische Astrophysik, Universität Heidelberg, Albert-Ueberle-Str. 2, 69120 Heidelberg, Germany

Received 18 December 2009 / Accepted 19 February 2010

Abstract
Context. Supersonic turbulence in the interstellar medium plays an important role in the formation of stars. The origin of this observed turbulence and its impact on the stellar initial mass function (IMF) still remain open questions.
Aims. We investigate the influence of the turbulence forcing on the mass distributions of gravitationally unstable cores in simulations of isothermal supersonic turbulence.
Methods. Data from two sets of non-selfgravitating hydrodynamic FLASH3 simulations with external stochastic forcing are analysed, each with static grid resolutions of 2563, 5123 and 10243 grid points. The first set applies solenoidal (divergence-free) forcing, while the second set uses purely compressive (curl-free) forcing to excite turbulent motions. From the resulting density field, we compute the mass distribution of gravitationally unstable cores by means of a clump-finding algorithm. Using the time-averaged probability density functions of the mass density, semi-analytic mass distributions are calculated from analytical theories. We apply stability criteria that are based on the Bonnor-Ebert mass resulting from the thermal pressure and from the sum of thermal and turbulent pressure.
Results. Although there are uncertainties in applying of the clump-finding algorithm, we find systematic differences in the mass distributions obtained from solenoidal and compressive forcing. Compressive forcing produces a shallower slope in the high-mass power-law regime compared to solenoidal forcing. The mass distributions also depend on the Jeans length resulting from the choice of the mass in the computational box, which is freely scalable for non-selfgravitating isothermal turbulence. If the Jeans length corresponding to the density peaks is less than the grid cell size, the distributions obtained by clump-finding show a strong resolution dependence. Provided that all cores are numerically resolved and most cores are small compared to the length scale of the forcing, the normalised core mass distributions are close to the semi-analytic models.
Conclusions. The driving mechanism of turbulence has a potential impact on the shape of the core mass function. Especially for the high-mass tails, the Hennebelle-Chabrier theory implies that the additional support due to turbulent pressure is important.

Key words: hydrodynamics - ISM: clouds - ISM: kinematics and dynamics - methods: numerical - stars: formation - turbulence

1 Introduction

The observed supersonic turbulence in the interstellar medium (ISM) plays an important role in the process of star formation (e.g., Scalo & Elmegreen 2004; Elmegreen & Scalo 2004; McKee & Ostriker 2007; Mac Low & Klessen 2004; Ballesteros-Paredes et al. 2007). The origin of this turbulence and the characteristics of different turbulence-driving mechanisms are still a matter of debate. Owing to the ability to counterbalance gravity globally and to provoke gravitational collapse locally, turbulence is expected to have a strong impact on the shape of the stellar initial mass function (IMF). The high-mass end ( $m\ga 0.5~ {M_{\hbox{$\odot$ }}}$) is typically observed to exhibit a power law of the form $\mathcal{N}(m){\rm d}m \propto m^{-\alpha}{\rm d}m$ (e.g., Lada & Lada 2003; Kroupa 2001; Salpeter 1955; Chabrier 2003), where $\mathcal{N}(m)$ is the number of stars in the linear mass interval ${\rm d}m$ with Salpeter (1955) power-law index $\alpha=2.35$. Observations of dense cores embedded in star-forming regions of turbulent giant molecular clouds show a similar but slightly shallower power-law distribution of mass with an exponent $\alpha$ in the range of 1.5-2.5 (Motte et al. 1998; Elmegreen & Falgarone 1996). The origin of the IMF and the Salpeter power law for the high-mass tail are still under debate, but there is a possible connection between the core mass function (CMF) and the IMF Motte et al. 1998; Johnstone & Bally 2006; Nutter & Ward-Thompson 2007; Johnstone et al. 2001; Williams et al. 2000; Testi & Sargent 1998; Johnstone et al. 2000; Alves et al. 2007; see however, the critical discussion by Clark et al. 2007. Padoan & Nordlund (2002) propose a theoretical explanation for the CMF/IMF based on scaling properties of turbulence and the Jeans criterion for gravitational instability. Combining the formalism of Press & Schechter (1974) for cosmological structure formation with the notion of a turbulent Jeans mass, a new analytic theory of the IMF was formulated by Hennebelle & Chabrier (2008).

Numerical simulations of turbulent molecular clouds reported in the literature (for instance, Federrath et al. 2010b; Banerjee et al. 2009; Hennebelle et al. 2008; Schmidt et al. 2009; Kritsuk et al. 2007) inject the turbulent energy via different methods, but apart from an early study by Klessen (2001), there has been no systematic analysis of how different driving schemes affect the CMF. Possible physical mechanisms for maintaining ISM turbulence range from stellar feedback (supernovae, outflows, ionizing radiation) to large-scale dynamical instabilities in the Galaxy (Elmegreen & Scalo 2004; Norman & Ferrara 1996; Mac Low & Klessen 2004), in this work we investigate the influence of the compressibility of the turbulence forcing on the mass spectra of gravitationally unstable cores.

We apply a clump-finding algorithm to the density fields from recent simulations of supersonic isothermal gas, using the FLASH3 code to solve the equations of hydrodynamics on static grids with resolutions of 2563, 5123 and 10243 grid points and periodic boundary conditions (Federrath et al. 2010b,2009,2008b). The turbulent energy is continuously injected by a specific force f on length scales corresponding to wavenumbers 1<k<3, where k is normalised by $2\pi/X$ for a box of size X. Each component of this force is modelled by an Ornstein-Uhlenbeck process (Federrath et al. 2010b; Eswaran & Pope 1988; Schmidt et al. 2006,2009). We define the integral scale L by the mean wavelength of the forcing, i.e., L=X/2. The relative importance of solenoidal and compressive modes of the stochastic forcing is set by the parameter $\zeta \in [0,1]$. The simulations represent two extreme cases: Purely solenoidal forcing ( $\nabla\cdot\textbf{\textit{f}}=0$, $\zeta=1$), on the one hand, and purely compressive forcing ( $\nabla~\times~\textbf{\textit{f}}=0$, $\zeta=0$), on the other, driven to an RMS Mach number of ${\fam=2 M}_{{\rm rms}}\approx 5.5$. Turbulence becomes statistically stationary for $t \geq 2\ \mathcal{T}$, where $\mathcal{T}$ is the dynamical timescale. In order to take the intermittent nature of the density field into account, we analyse each of the 81 density snapshots available in the time interval $2~\mathcal{T} \leq t \leq 10~\mathcal{T}$ to compute time-averaged core statistics. For a detailed description of the simulation setup and numerical techniques, see Federrath et al. (2010b). Recent results of isothermal turbulence simulations indicate a strong dependence of the density field statistics on the compressibility of the stochastic forcing (Federrath et al. 2010b; Schmidt et al. 2009; Federrath et al. 2009,2008b). This effect makes it is a reasonable assumption that the resulting core mass distributions should show different properties for different turbulence driving mechanisms. Furthermore, the datasets enable us to study in which way the effective resolution of the simulation affects the resulting mass distributions.

This article is structured as follows. In the next section, we review the analytical theories for the CMF by Padoan & Nordlund (2002) and Hennebelle & Chabrier (2008). For the sake of clarity, we express the main results of these theories in a consistent notation. In Sect. 3, we discuss the parameters determining the global mass scale, which plays a key role in computations of the core mass functions. The clump-finding algorithm is briefly described in Sect. 4, and the resulting mass distributions are presented. Particular emphasis is put on the question of numerical convergence. The distributions obtained by clump-finding are compared to semi-analytic computations of the CMF in Sect. 5, and the results are summarised and discussed in Sect. 6.

2 Analytical theories for the core mass function

In this work, we follow Hennebelle & Chabrier (2008) by writing $\mathcal{N}(M)={\rm d}N/{\rm d}M$, where N(M) is the cumulative number of cores with masses below M. The function $\mathcal{N}(M)$ is also called the mass spectrum. Because of the wide range of different masses, it is useful to define the CMF by the number of cores per logarithmic mass interval. To that end, we define two dimensionless mass variables:

  • $m:=M/M_{\hbox{$\odot$ }}$ is the mass relative to the solar mass.
  • $\tilde{M}:=M/M_{{\rm J}}^{0}$ is the mass relative to the thermal Jeans mass at the mean density $\bar{\rho}$.
While it is useful to express the CMF as a function of m in an observational context, $\tilde{M}$ is the most natural choice for theoretical considerations. The linear mass distribution is related to the distributions with respect to $\log m$ or $\log\tilde{M}$ as follows[*]:

\begin{displaymath}
\mathcal{N}(M) =
(M_{\hbox{$\odot$ }}m)^{-1}\frac{{\rm d}N...
...tilde{M})^{-1}\frac{{\rm d}N}{{\rm d}{\rm log}~\tilde{M}}\cdot
\end{displaymath} (1)

Integration of over M yields the total number of cores,

\begin{displaymath}
\begin{array}{rl}
N_{{\rm tot}} =& \int_{0}^{\infty}\mathca...
... d}{\rm log}~\tilde{M}}~{\rm d}{\rm log}~\tilde{M}.
\end{array}\end{displaymath} (2)

2.1 Padoan-Nordlund theory

Padoan & Nordlund (2002) based their analytical approach on the following assumptions. The core size is comparable to the thickness of isothermal postshock gas, the turbulent velocity fluctuations follow a power law, and the minimal mass that is unstable against gravitational collapse is given by the thermal Jeans mass. They show that the number of cores per mass interval is given by an expression of the form

\begin{displaymath}
\mathcal{N}(M) \propto
M^{-(x+1)}\int_{0}^{M}{\rm pdf}_{{\rm J}}(M')~{\rm d}M',
\end{displaymath} (3)

where the parameter x depends on the scaling properties of turbulence. Integration of the probability density function ${\rm pdf}_{{\rm J}}(M')$ of the Jeans mass for a given distribution of the mass density yields the fraction of cores with masses $M'\le M$ that are unstable according to the Jeans criterion, $M'\propto\rho^{-1/2}$. For cores of sufficiently high mass, this integral asymptotically approaches a constant. For this reason, the high-mass tail of $\mathcal{N}(M)$ obeys a power law with exponent -(x+1).

Equation (3) can be written in terms of the probability density function ${\rm pdf}(\delta)$ of the logarithmic density fluctuation $\delta=\log(\rho/\bar{\rho})$. Since the Jeans stability criterion implies $\log\tilde{M}=-\frac{1}{2}\delta$, we have ${\rm pdf}_{{\rm J}}(M')=-2{\rm pdf}(\delta)/M'$, and it follows that

\begin{displaymath}
\frac{{\rm d}N}{{\rm d}{\rm log}~\tilde{M}} \propto
\tild...
...int_{-2\log\tilde{M}}^{\infty}{\rm pdf}(\delta)~{\rm d}\delta.
\end{displaymath} (4)

The upper bound on the core mass corresponds to a lower bound on the density fluctuation. The distribution (4) is invariant under a change in the global scales. Obviously, this cannot be the physical CMF if $M_{{\rm J}}^{0}\ga M_{{\rm tot}}$, where $M_{{\rm tot}}:=(2L)^{3}\bar{\rho}$ is the total mass in the computational domain (also see the discussion of the mass scale in Sect. 3). For this reason, the Padoan-Nordlund theory has to be considered as an asymptotic description for $M_{{\rm J}}^{0}\ll M_{{\rm tot}}$. The analytic calculation of the mass spectrum is straightforward if the probability density function of the mass-density is lognormal, i.e.,

\begin{displaymath}{\rm pdf}(\delta)=
\frac{1}{\sqrt{2\pi\sigma^{2}}}\exp\left[-\frac{(\delta-\bar{\delta})^{2}}{2\sigma^{2} }\right]\cdot
\end{displaymath} (5)

where $\bar{\delta}=-\sigma_{0}^{2}/2$.

Originally, the power-law exponent x was obtained from the jump conditions for shocks in magnetohydrodynamic turbulence, and it was shown that the resulting slope of the high-mass tail is close to the Salpeter slope $x=\alpha-1=1.35$ (see Padoan & Nordlund 2002). Padoan et al. (2007), hereafter referred to as PN07, also determined x for isothermal hydrodynamic (HD) turbulence:

\begin{displaymath}
x_{{\rm HD}}=\frac{3}{5-2\beta},
\end{displaymath} (6)

where $\beta$ is the slope of the turbulence energy spectrum. Using data from a driven isothermal HD turbulence simulation, PN07 calculated numerical mass distributions by means of clump-finding for purely thermal support. The spectral index $\beta=1.9$ (Kritsuk et al. 2007) implies $x_{{\rm HD}}=2.5$ for the slope of the high-mass tail. This value recovers the general trend of their numerical mass distributions well (see Figs. 2 and 3 in PN07). However, they did not perform time-averaging for their mass distributions, so no estimate of the statistical significance of their results was provided. It was concluded that the high-mass tail is generally too stiff for HD turbulence to be consistent with the observed CMF.

2.2 Hennebelle-Chabrier theory

The theory of Hennebelle & Chabrier (2008) provides a general framework for mass distributions of cores on the basis of the Press-Schechter statistical formalism. In contrast to the theory by Padoan & Nordlund (2002), there is an inherent notion of a core scale R that is connected to its mass. Assuming a log-normal distribution of the gas density, analytic formulae for the mass distribution $\mathcal{N}(M)$ can be derived. In the simplest case, M is given by the thermal Jeans mass for a given density. While the Padoan-Nordlund theory only uses thermal support, the Hennebelle-Chabrier theory considers a combination of thermal and turbulent support. On the one hand, turbulence is promoting star formation through the broadening of the density PDF, while on the other, it also quenches star formation through the scale-dependent support of the cores.

By applying an approximation that excludes the largest cores, it is possible to generalise their approach to arbitrary PDFs of the mass density. For cores on a length scale R comparable to the integral scale L, the mass distribution depends on the derivative of ${\rm pdf}(\delta)$ with respect to R, which in turn depends on the slope of the spectrum of the density fluctuations. It would be worthwhile investigating the influence of the spectral properties of the density field on the CMF. Because of the difficulties in the numerical determination of the derivative of the $\delta$-PDF, however, we restrict our investigation to length scales $R\ll L$, for which the general form of the mass spectrum (Eq. (33) in Hennebelle & Chabrier 2008) simplifies to

\begin{displaymath}
\mathcal{N}(M) \simeq
-\bar{\rho}\left(M\frac{{\rm d}M}{{\...
...
\frac{{\rm d}\delta}{{\rm d}R}\exp(\delta){\rm pdf}(\delta).
\end{displaymath} (7)

For a comparison with the Padoan-Nordlund theory, the case of purely thermal support is treated as follows. The core mass $M=M_{{\rm J}}^{0}\tilde{M}$ is related to the logarithmic density fluctuation $\delta$ by the Jeans criterion for stability against gravitational collapse, and R is given by the thermal Jeans length, hence,

\begin{displaymath}
\tilde{M}=\tilde{R}=\exp(-\delta/2).
\end{displaymath} (8)

The dimensionless scale parameter is defined by $\tilde{R}:=R/\lambda_{{\rm J}}^{0}$, where $\lambda_{{\rm J}}^{0}:=a_{{\rm J}}^{1/3}c_{{\rm s}}/(G\bar{\rho})^{1/2}$ is the thermal Jeans length for the mean density ($c_{\rm s}$ is the isothermal speed of sound, G the gravitational constant, and $a_{{\rm J}}$ a geometry factor). Substituting the inverse mass-density relation (8), it follows from Eq. (7) that

\begin{displaymath}
\frac{{\rm d}N}{{\rm d}{\rm log}~\tilde{M}}=
-\frac{\bar{\rho}}{M_{{\rm J}}^{0}}\tilde{M}^{-3}{\rm pdf}(-2\log\tilde{M}).
\end{displaymath} (9)

This expression for the CMF is fully determined by the lower bound of the integral in Eq. (4), and, generally, the high-mass tail does not obey a power law.

If turbulence pressure contributes to the support against gravitational collapse, an implicit relation between the gravitationally unstable mass and the density fluctuation results:

$\displaystyle \tilde{M} = \tilde{R}\left(1+\mathcal{M}_{\ast}^{2}\tilde{R}^{2\eta}\right),$         $\textstyle \quad \delta = \log\left(\frac{1+\mathcal{M}_{\ast}^{2}\tilde{R}^{2\eta}}{\tilde{R}^{2}}\right)\cdot$ (10)

The intensity of turbulence is specified by the Mach number $\mathcal{M}_{\ast}$ of turbulent velocity fluctuations on the length scale $\lambda_{{\rm J}}^{0}$,

\begin{displaymath}
\mathcal{M}_{\ast} := \frac{\mathcal{M}_{{\rm rms}}}{\sqrt{3}}
\left(\frac{\lambda_{{\rm J}}^{0}}{L}\right)^{\eta},
\end{displaymath} (11)

where $\eta=(\beta-1)/2$, and $\beta$is the slope of the turbulence energy spectrum function in the inertial subrange. Since the turbulent pressure on the length scale $\lambda_{{\rm J}}^{0}$ equals $\rho(\mathcal{M}_{{\rm\ast}}c_{{\rm s}})^{2}$, the parameter $\mathcal{M}_{\ast}$ measures the relative significance of turbulent vs. thermal support for cores of size $\sim \lambda_{{\rm J}}^{0}$. The turbulent Mach number on the length scale R is given by $\mathcal{M}_{\ast}\tilde{R}^{\eta}$. Note that $\mathcal{M}_{\ast}\tilde{R}^{\eta}$ is about unity if R is close to the sonic length scale (Federrath et al. 2010b; Schmidt et al. 2009). The CMF including turbulent support is given by

\begin{displaymath}
\begin{array}{rl}
\displaystyle \frac{{\rm d}N}{{\rm d}{\rm ...
...\rm pdf}\left[\log(\tilde{M}/\tilde{R}^{3})\right],
\end{array}\end{displaymath} (12)

where $\tilde{R}$ is numerically determined by the inversion of Eq. (10) for a given value of $\tilde{M}$.

3 Clump-finding and scaling

For the calculation of the mass distributions we applied the implementation of the clump-finding algorithm by PN07. As a first step the algorithm divides the three-dimensional density field into k discrete density levels $\rho_{i}$ based upon the settings of the input parameters $\rho_{{\rm min}}$ and $f=\rho_{i+1}/\rho_{i}$ with $i=1\ldots k$. The former parameter defines the minimum density level in units of the average density $\bar{\rho}$, while the latter defines the spacing between two adjacent levels. In the next step each density level is scanned for regions of connected cells with density values higher than the current density level. A connected region is counted as an object (a core) if its mass exceeds the local Bonnor-Ebert mass (Ebert 1955; Bonnor 1956). The ratio of the Bonnor-Ebert mass to the solar mass is given by

\begin{displaymath}
\begin{array}{rl}
m_{{\rm BE}}&= \displaystyle \frac{1.18 c...
...ght)^{-1/2}\left(\frac{T}{10~{\rm K}}\right)^{3/2},
\end{array}\end{displaymath} (13)

where the gravitational constant $G=6.67\times 10^{-8}~ {\rm cm^3~ g^{-1}~ s^{-2}}$, $c_{\rm s}$ is the isothermal sound speed, T the temperature, and $\mu$ the mean molecular weight of the gas. The number density $n=\rho/(\mu m_{{\rm H}})$, where $m_{{\rm H}}$ is the mass of the hydrogen atom, is averaged over the connected region. In the last step, each object that can be split into two or more objects is rejected. For a more detailed description of the clump-finding method, see PN07.

For this analysis, it is crucial to consider the degrees of freedom in the scaling of isothermal, non-self-gravitating turbulence. For a given temperature T and box size L, the hydrodynamical scales of the system are fixed, and the forcing magnitude determines the RMS Mach number $\mathcal{M}_{{\rm rms}}$of turbulence in the statistically stationary state. If self-gravity is not included, the flow properties are invariant with respect to the choice of the mean mass density $\bar{\rho}$. This is palpable if the compressible Euler equations are written in terms of the logarithmic density $\delta=\ln(\rho/\bar{\rho})$ (see, for instance, Federrath et al. 2010b). Therefore, the mass scale is arbitrary.

On the other hand, if gravitationally unstable cores are identified via postprocessing, then the value of $m_{{\rm BE}}$ (13) resulting from the chosen values of n and T fixes the mass scale. In this paper, we characterise the mass scale by the total number $N_{{\rm BE}}$ of Bonnor-Ebert masses with respect to the mean density contained in the simulation box:

\begin{displaymath}
\begin{array}{rl}
N_{{\rm BE}} =&\displaystyle \frac{m_{{\rm...
...^{-3/2}\left(\frac{2L}{1.782\ {\rm pc}}\right)^{3},
\end{array}\end{displaymath} (14)

where $m_{{\rm BE}}(\bar{\rho})$ is the Bonnor-Ebert mass for the mean density $\bar{\rho}$, and $\bar{n}=\bar{\rho}/(\mu m_{{\rm H}})$. The ratio of the thermal Jeans length with the geometry factor $a_{{\rm J}}=1.18$ corresponding to the definition of the Bonnor-Ebert mass to the integral length is given by

\begin{displaymath}\frac{\lambda_{{\rm J}}^{0}}{L}=\frac{2}{N_{{\rm BE}}^{1/3}}\cdot
\end{displaymath} (15)

The Hennebelle-Chabrier theory accounts for the dependence of the core statistics on the parameter $N_{{\rm BE}}$ via the factor $\bar{\rho}/M_{{\rm J}}^{0}\sim N_{{\rm BE}}^3/(2L)^3$ in the distributions (9) and (12). Apart from that, variations of the mass scale correspond to changes of the scale-dependent Mach number,

\begin{displaymath}
\mathcal{M}_{\ast} =
\frac{\mathcal{M}_{{\rm rms}}}{\sqrt{3}}
\left(\frac{N_{{\rm BE}}}{8}\right)^{-\eta/3},
\end{displaymath} (16)

where $\mathcal{M}_{\ast}$ characterises the ratio of turbulent pressure to thermal pressure on the length scale $\lambda_{{\rm J}}^{0}$ (see Eq. (11)). To include the influence of the turbulent pressure on the stability against gravitational collapse in the clump-finding tool, we modified the stability criterion by replacing the thermal speed of sound in the definition of the Bonnor-Ebert mass with an effective speed of sound (see, for instance, Bonazzola et al. 1987):

\begin{displaymath}
c_{{\rm eff}}^{2} = c_{{\rm s}}^{2} + \frac{1}{3}\sigma_{{\rm core}}^{2}(v).
\end{displaymath} (17)

The turbulent velocity dispersion $\sigma_{{\rm core}}(v)$ is computed from the mass-weighted RMS velocity fluctuation with respect to the centre-of-mass velocity for a particular core.

The physical scales that were chosen by PN07 are $n = 10^4\;{\rm cm^{-3}}$, $X = 2L = 6\;{\rm pc}$ and $T = 10\;{\rm K}$. In this case, $N_{{\rm BE}}\approx 1.2\times 10^{5}$ and $\lambda_{{\rm J}}^{0}/L\approx 0.04$. For the simulations with the highest resolution, we have $\lambda_{{\rm J}}^{0}\approx 20.8\Delta$, where $\Delta$ is the size of the grid cells. As we see in the next section, this introduces a serious resolution issue, because many cores only extend over a few cells for this parameter set. Another difficulty is that the above parameters are atypical for observed molecular cloud properties. According to the Larson (1981) relation, $\bar{n}\approx 3000\;{\rm cm^{-3}}(L/1~{\rm pc})^{-1}$ (see also Heyer et al. 2009; Falgarone et al. 1992). This suggests that the assumed forcing scale $L=X/2=3\;{\rm pc}$is too large for molecular clouds of mean number density $\bar{n}=10^4\;{\rm cm^{-3}}$. To analyse the impact of the physical scales, we computed the core mass distributions for $L = 0.6\;{\rm pc}$ again without changing the temperature and the mean number density. For these parameters, the number density is consistent with the Larson relation within a factor of two, which is reasonable given the scattering of observed molecular cloud properties (Falgarone et al. 1992). In this case, $N_{{\rm BE}}=10^{3}$, and $\lambda_{{\rm J}}^{0}/L\approx 0.2$. We also note, however, that the validity of Larson-type relations for the mass density was questioned by Ballesteros-Paredes & Mac Low (2002).

4 Core statistics

We begin with a detailed analysis of the mass distributions obtained for purely thermal support and the PN07 mass scale corresponding to $N_{{\rm BE}}\approx 1.2\times 10^{5}$. In particular, the tuning of the clump-finding algorithm and influence of numerical resolution are considered in Sects. 4.1 and 4.2. Finally, the mass distributions for fewer Bonner-Ebert masses in the box, $N_{{\rm BE}}\approx 1.0\times 10^{3}$, as well as the influence of turbulent pressure are discussed in Sect.  4.3.

\begin{figure}
\par\subfigure[Solenoidal forcing.]{\includegraphics[width=9cm,cl...
...{\includegraphics[width=9cm,clip]{13904fg1b.eps} }
\vspace*{2.8mm}
\end{figure} Figure 1:

The core mass distributions for the 5123 solenoidal a) and compressive b) simulations as a function of the clump-finding algorithm parameter f (increasing from top curve to bottom curve) which sets the relative spacing between two adjacent density levels. Error bars contain the 1$\sigma $ temporal fluctuations and are only indicated for f=1.04 for the sake of clarity.

Open with DEXTER

Table 1:   Time-averaged properties of the mass distributions for solenoidal and compressive forcing calculated with the physical scalings applied in PN07.

4.1 Tuning of the clump-finding algorithm

As described in Sect. 3 the clump-finding routine needs two technical input parameters: the relative spacing f between two adjacent density levels and the minimum density level $\rho_{{\rm min}}$ above which the datacube is divided into discrete levels. The influence of these parameters has already been described in PN07. For the minimum density level, we choose the mean density $\rho_{0}$ of the simulation box, since we observe no significant differences in the resulting core mass distributions for values of $\rho_{{\rm min}} < \rho_{0}$, which is also found in PN07. The influence of the level spacing f on the resulting mass distributions is shown in Fig. 1. The slope of the high-mass range and the total number of cores is very sensitive to the chosen spacing of the density levels. For higher f-values the algorithm detects less but more massive cores. For decreasing spacing between density levels and therefore finer ``scanning'' of the datacube, the massive cores are split into many more lower mass cores, and this leads to a steepening of the high-mass tail of the core mass distributions. For f=1.04, 1.08, and 1.16, this effect lies within the temporal fluctuations of the mass distributions as indicated by the error bars in Fig. 1. PN07 also conclude that the mass distributions are basically converged for $f \leq 1.16$. For the calculations in the present study, we chose a value of f=1.04. Nevertheless, our numerical mass distributions should be compared to the observed IMF/CMF with caution, especially, in the high-mass range. The parameter dependency of core mass distributions is a common problem for all clump-finding methods in the literature. Smith et al. (2008) apply different clump-finding methods to data of SPH simulations of molecular clouds with the result that the definition of an individual core and its mass depend strongly on the method and the parameter settings (see also Klessen & Burkert 2000; Klessen 2001). In the observers' community, the clump-finding routines CLUMPFIND (Williams et al. 1994) and GAUSSCLUMPS (Stutzki & Guesten 1990) are commonly used to decompose position-position-velocity information of molecular clouds into individual cores. These techniques suffer from the same parameter dependency, as shown by Schneider & Brooks (2004) and more recently by Pineda et al. (2009).

\begin{figure}
\par\subfigure[Solenoidal forcing, $N_{{\rm BE}}=1.2\times10^{5}$...
...]{\includegraphics[width=9cm,clip]{13904fg2d.eps} }
\vspace*{2.5mm}
\end{figure} Figure 2:

Core mass distributions for solenoidal (left) and compressive forcing (right) at numerical resolutions of 2563 (dashed), 5123 (dot-dashed), and 10243 (solid), normalised to the total number of cores, $N_{{\rm tot}}$. Error bars indicate 1$\sigma $ temporal fluctuations.

Open with DEXTER

4.2 Numerical resolution study for purely thermal support

Figure 2 indicates that the mass distributions are very sensitively to the grid resolution of the simulation. For $N_{{\rm BE}}=1.2\times 10^{5}$, the high-mass tails steepen systematically with increasing grid resolution for both solenoidal and compressive forcing, showing no sign of convergence in the sense of an asymptotic limit, even for the highest resolution of 10243 grid points. The mass distributions for $N_{{\rm BE}}=10^{3}$appear to be better converged, at least in the case of solenoidal forcing. The statistical properties of the core mass distributions are summarised in Table 1. The mean mass $\overline{m}$ shifts to lower values with increasing resolution, approaching a time-averaged value of $\left< \overline{m}\right> \approx 1~{M_{\hbox{$\odot$ }}}$ for the 10243 runs, while the width $\sigma_{m}$ of the distributions decreases. This is caused by the fragmentation of the cores found in the low resolution simulations into smaller and denser cores in the high resolution simulations which also satisfy the condition for gravitational instability. The sum of all core masses, $m_{{\rm cores}}$, is a few percent of the total box mass $m_{{\rm tot}}$.

An important property of the mass distributions is the slope x of the high-mass wings of ${\rm d}N/{\rm d}{\rm log}~m$, which is related to the exponent $\alpha$ of the linear distribution via $x=\alpha-1$. Assuming that there is a power-law range, we applied error-weighted least squares fits to the high-mass power-law regime (2563: $m \geq 15~{M_{\hbox{$\odot$ }}}$, 5123: $m \geq 8~{M_{\hbox{$\odot$ }}}$, 10243: $m \geq 3~{M_{\hbox{$\odot$ }}}$) in Fig. 2. The power-law exponents obtained by this method are referred to as $x_{{\rm LS}}$. For the highest resolution, we find $x_{{\rm LS}}=3.1\pm 0.2$ for solenoidal forcing and a shallower slope of $x_{{\rm LS}}=2.1\pm 0.2$ for compressive forcing. The mass distributions of the solenoidal runs tend to show a higher degree of convergence. The differences in the slopes between the three solenoidal simulations decrease with increasing resolution. The slope $x_{{\rm LS}}$ obtained from the 5123 simulation is steeper by a value of 1.0 than the 2563 while the values of $x_{{\rm LS}}$ for the resolutions of 5123 and 10243 are nearly identical in the solenoidal case. In contrast, the slopes of the mass distributions obtained from the pure compressive forcing increase by a value of 0.6 ( $256^3\rightarrow 512^3$) and by a value of 0.5 ( $512^3\rightarrow 1024^3$), indicating slower convergence than in the solenoidal forcing case.

We also used the method of maximum likelihood estimation (MLE) (for details see Appendix A) inspired by the work of Clauset et al. (2007) to determine the slope $x_{{\rm MLE}}$ of the core mass distributions and to check that the distribution really follows a power law. This method avoids biases that can result from least-square fits (e.g., Bauke 2007) and provides an estimate of the mass value $m_{{\rm min}}$ above which the power-law assumption is made. Furthermore, the probability p that the power-law assumption holds can be calculated. Following Clauset et al. (2007), we rule out power laws if $p \leq 0.1$. The results of the MLE method applied to the data of our core masses are summarised in Table 1. For the highest resolution, we find a power-law slope of $x_{{\rm MLE}}=2.2\pm0.3$ for $m\geq m_{{\rm min}}=2.7\pm 0.8~{M_{\hbox{$\odot$ }}}$ in the case of compressive forcing and $x_{{\rm MLE}}=3.2\pm0.7$ for $m\geq m_{{\rm min}}=3.5\pm 0.9~{M_{\hbox{$\odot$ }}}$ in the solenoidal forcing case. The errors are the 1$\sigma $temporal fluctuations of these values. The slopes are consistent with the values obtained via least squares fits. At numerical resolutions of 2563 and 5123, the time-averaged p-values are less than 0.1 for both solenoidal and compressive forcing. Thus, they are not consistent with power laws. However, at a numerical resolution of 10243, we find $p=0.21\pm0.26$ for solenoidal and $p=0.23\pm0.21$ for compressive forcing. The standard deviations of these p-values are very large, such that instantaneous CMFs can either be consistent with power laws or do not exhibit power-law behaviour at the high-mass end of the distribution. We emphasise that it is absolutely necessary to compute time-averaged quantities in order to obtain statistically meaningful results. The p-value can only rule out the power law hypothesis, while other distributions could be a better fit even for $p\geq0.1$. The possible inconsistency of the high-mass CMFs with a power law-model was also discussed by Ballesteros-Paredes et al. (2006).

Table 2:   Least-square estimates of the power-law exponent $x_{{\rm LS}}$obtained from the 10243 simulations for different values of $N_{{\rm BE}}$ with/without turbulent pressure included in the core stability criterion.

\begin{figure}
\subfigure[Solenoidal forcing, $N_{{\rm BE}}=1.2\times10^{5}$ , t...
...t.]{\includegraphics[width=9cm,clip]{13904fg3b.eps} }
\vspace*{2mm}
\end{figure} Figure 3:

PDF of core size r defined in Eq. (18) for solenoidal a) and compressive b) forcing for a numerical resolution of 2563 (dashed), 5123 (dot-dashed), and 10243 (solid) and PN07 scaling. Error bars indicate 1$\sigma $ temporal fluctuations of the PDFs.

Open with DEXTER

For $N_{{\rm BE}}=10^{3}$, the scatter around the tails of the mass distributions is too large for the MLE method to be applicable. Therefore, we only list the power-law exponents $x_{{\rm LS}}$ of the high-mass tails of the CMFs in Table 2. Although it is quite inaccurate, especially because of the low number statistics, this method can be used to describe the general trend of the CMFs. For the solenoidal simulation with turbulent support and $N_{{\rm BE}}=10^{3}$, we do not give a value for the power-law exponent because the error bars are simply too large. Decreasing the total number of Bonner-Ebert masses in the box from $N_{{\rm BE}}$ from $1.2\times 10^{5}$ to 103 leaves the fitted value of $x_{{\rm LS}}$ relatively unchanged (taking the large error bars into account). The effect of adding turbulent support on the slope of the CMF is discussed in Sect. 4.3.

To further investigate the effects of numerical resolution on the core properties, we looked at the size r of the cores. We defined the approximate size of a core as

\begin{displaymath}
r=n^{1/3}_{{\rm cells}}\Delta,
\end{displaymath} (18)

where $n_{{\rm cells}}$ is the number of grid cells contained in a core and $\Delta=X/N$ with N=256, 512, 1024 depending on the numerical grid resolution. We first consider the mass scale chosen by PN07, i.e., $N_{{\rm BE}}\approx 1.2\times 10^{5}$. Then we comment on the changes arising from a lower number of Jeans masses in the box (corresponding to $N_{{\rm BE}}\approx \times10^{3}$). The typical size of large cores should be $\sim \lambda_{{\rm J}}^{0}$, which is about $5.2\Delta$ for N=256 and $20.8\Delta$ for N=1024. The fragmentation properties of the gas on length scales less than about $10\Delta$, which encompasses most of the range of core sizes, are affected by numerical smoothing. Consequently, the probability density functions (PDFs) of r, which are plotted in Fig. 3, change substantially with resolution. In particular, one can also see that the fractions of cores with $r\ga \lambda_{{\rm J}}^{0}$ are relatively small, especially for the high-resolution data. At lower resolution, however, these fractions increase; i.e., more big cores are found relative to the number of smaller cores. Thus, the maxima of the mass distributions shown in Fig. 2 are shifted towards higher masses with decreasing resolution. In addition, the discrepancies in the high-mass wings might be caused by a bias of the clump-finding algorithm to select cores that contain several gravitationally unstable cores.

Table 3:   Time-averaged properties of the core size PDFs for solenoidal and compressive forcing.

For a certain variance, $\sigma^{2}$, of the density fluctuations, we expect that cores of size $\sim\sigma^{-1/2}\lambda_{{\rm J}}^{0}$ are most frequent, because the Jeans length changes with the inverse square root of the mass density. From the values of $\sigma $ for solenoidal and compressive forcing (see Federrath et al. 2010b, Table 1), it follows that $(\sigma_{\rm comp}/\sigma_{\rm soln})^{1/2}=\sqrt{5.9/1.9}\approx1.76$. This agrees with the relative peak positions of the distributions for N=1024 in Fig. 3. Moreover, the mean size $\overline{r}$ of the cores for solenoidal forcing is about the same factor more than $\overline{r}$ for compressive forcing, as one can see from the values listed in Table 3.

The minimal core size is roughly given by the peak densities $\rho_{{\rm max}}$ in the turbulent gas. The definition of the Jeans length implies $r_{{\rm min}}\sim(\bar{\rho}/\rho_{{\rm max}})^{1/2}\lambda_{{\rm J}}^{0}$. Since $\rho_{{\rm max}}$ is roughly $500\bar{\rho}$ for solenoidal forcing and $10^{4}\bar{\rho}$ for compressive forcing (Federrath et al. 2009), it follows that $r_{{\rm min}}\sim0.9\Delta$ and $r_{{\rm min}}\sim 0.2\Delta$, respectively, for the highest resolution. Consequently, the smallest cores are marginally resolved in the 10243 simulation with solenoidal forcing, while they are definitely below the resolution limit for compressively driven turbulence. As an indicator, we calculated the fractions $f_{r\leq 2}$ of cores with $r \leq 2$ (see Table 3). Indeed, significant fractions $f_{r\leq 2}$ were obtained for all resolutions in the case of compressive forcing. In the solenoidal simulations, on the other hand, $f_{r\leq 2}$ drops from $15.8\%$ (2563) to $0.037\%$ in the 10243 simulation. These trends are also visible in Fig. 3.

Setting $N_{{\rm BE}}=10^{3}$, the estimate of the minimal core size from the peak density yields $r_{{\rm min}}\approx 4.6\Delta$ for solenoidal forcing and $1.0\Delta$ for compressive forcing. In the solenodial case, $f_{r\leq 2}$ drops to zero for the highest resolution. Compared to the PN07 setting, the cores are in general larger (see Table 3 and the top panels of Fig. 4). Accordingly, the resolution dependence of the mass distributions is less pronounced, particularly in the case of solenoidal forcing. However, comparing the mass distributions without normalisation in Figs. 5a and 5c, one can see that the total number of cores decreases by two orders of magnitude if the mass scale is chosen such that $N_{{\rm BE}}=10^{3}$. As in the case $N_{{\rm BE}}\approx 1.2\times 10^{5}$, we find that the mean core size for solenoidal forcing is roughly twice as large as for compressive forcing.

\begin{figure}
\par\subfigure[Solenoidal forcing, $N_{{\rm BE}}=10^{3}$ , therma...
...{\includegraphics[width=8.9cm,clip]{13904fg4d.eps} }
\vspace*{-1mm}
\end{figure} Figure 4:

PDFs of core size distribution as in Fig. 3 but for $N_{{\rm BE}}=10^{3}$ with/without turbulent pressure support (see Sect. 2.2 and 3).

Open with DEXTER

4.3 Influence of the turbulent pressure

As explained in Sect. 3, we used the effective speed of sound according to Eq. (17) in the definition of the Bonnor-Ebert mass, Eq. (13), to compute the mass distributions of cores with thermal and turbulent support. For brevity, we use the term turbulent support, for which it is understood that the instability criterion is based on the sum of thermal and turbulent pressure.

\begin{figure}
\par\subfigure[$N_{{\rm BE}}=1.2\times10^{5}$ , thermal support.]...
...bulent support.]{\includegraphics[width=9cm,clip]{13904fg5d.eps} }\end{figure} Figure 5:

CMFs of compressive (solid) and solenoidal (dot-dashed) forcing for different values of $N_{{\rm BE}}$ with/without turbulent support (see Sects. 2.2 and 3) and a grid resolution of 10243. The least-square fits to the high-mass tails are shown as dashed lines.

Open with DEXTER

Figures 5b and 5d compare the CMFs for $N_{{\rm BE}}=1.2\times 10^{5}$ and $N_{{\rm BE}}=10^{3}$ with turbulent support for the highest grid resolution of 10243. While turbulent support has no significant impact on the CMF for $N_{{\rm BE}}=1.2\times 10^{5}$, it clearly changes the CMF for $N_{{\rm BE}}=10^{3}$. For solenoidal forcing, turbulent support leads to a large shift in the peak in the CMF, corresponding to the pronounced change of the core size distribution. As one can see in Fig. 5d, however, the temporal fluctuations are larger than the time average. This is a consequence of how few cores there are. In the compressive case with turbulent support, the peak of the CMF is relatively robust, but we observe a flattening of the high-mass tail from $x_{{\rm LS}}\approx 2.0$ for the CMF without turbulent support to $x_{{\rm LS}}\approx 1.2$ (see Table 2), which is close to the Salpeter value of 1.35.

Table 3 gives a summary of the basic properties of the core size distribution for all possible combinations of grid resolution, turbulent/thermal support, and $N_{{\rm BE}}$. If turbulent support is included, $f_{r\leq 2}$ becomes negligible both for solenoidal and for compressive forcing if N=1024. The core size PDFs are plotted in Fig. 4. The large error bars are the result of to find large cores increases substantially for solenoidal forcing, while the PDF for compressive forcing does not change appreciably. In the case of solenoidal forcing, the greatest effect on the core size is observed for the highest resolution of 10243 grid cells. Thus, the support of cores against gravitational collpase is greatly enhanced by the turbulent pressure on length scales that are sufficiently large compared to the grid scale. In Sect. 5.2, we show that the turbulent pressure exceeds the thermal pressure on length scales above $30\Delta$, i.e., outside the range of length scales that are subject to significant numerical dissipation.

5 Comparison to theoretical models

In this section, the core mass distributions computed with the clump-finding algorithm are compared to theoretical predictions of Padoan & Nordlund (2002) as well as Hennebelle & Chabrier (2008). The influence of the width of the density PDF resulting from different forcing was investigated by Hennebelle & Chabrier (2009). Because Federrath et al. (2010b,2008b) showed that the density PDFs are not exactly log-normal, we evaluated Eqs. (4), (9), and (12) for the numerical PDFs calculated from the simulation data. The PDFs for the three numerical resolutions are plotted in Federrath et al. (2010b, Figs. 4 and 6). The resulting CMFs are thus semi-analytic. This is particularly important for the case of compressively driven turbulence, for which the core mass distributions calculated from log-normal fits to the PDFs of the mass density deviate substantially from the semi-analytic distributions based on the numerical PDFs (see also Schmidt et al. 2009). For the comparison with the data from our clump-finding analysis, we normalise the distributions with the total number of cores $N_{{\rm tot}}$, which is calculated from Eq. (2). Allowing for an arbitrary geometry factor $a_{{\rm J}}$ in the definition of the thermal Jeans mass (see Hennebelle & Chabrier 2008), we set $a_{{\rm J}}=1.18$ and identify $M_{\rm J}^{0}$ with the Bonnor-Ebert mass $m_{{\rm BE}}(\rho)M_{\hbox{$\odot$ }}$ (see Eq. (13)). With this, we have a criterion for gravitational instability on the basis of the thermal pressure that is consistent with the clump-finding algorithm.

5.1 Purely thermal support

Substitution of the time-averaged PDFs obtained from 10243 simulations (Federrath et al. 2010b,2008b) into Eqs. (4) and (9), yields the semi-analytic CMFs plotted in the top panels (a,b) of Fig. 6. These distributions are independent of the parameter $N_{{\rm BE}}$ because of the normalisation by the total number of cores, $N_{{\rm tot}}$. As one can see, both the Padoan-Nordlund and the Hennebelle-Chabrier theory imply that the peak of the CMF is shifted towards lower masses for compressively driven turbulence. This is a direct consequence of the shape of the density PDFs. Since compressive forcing produces higher density peaks, there is a higher fraction of low Jeans masses compared to solenoidal forcing. This results in the formation of smaller cores, as noted in Sect. 4.2). As a consequence, the high-mass tail of the CMF (9) is significantly less stiff in the compressive case (the asymptote for $\tilde{M}\gg 1$ is mainly given by the factor ${\rm pdf}(-2\log\tilde{M})$. The slope of the CMF (4), on the other hand, is determined by the factor $\tilde{M}^{-x}$, where x is defined by Eq. (6), in the limit of high masses. The turbulence energy spectra computed from the simulation data yield the power-law exponents $\beta=1.86\pm0.05$ for solenoidal forcing and $\beta=1.94\pm0.05$ for compressive forcing (Federrath et al. 2010b). Thus, the slopes of the high-mass tails predicted by the Padoan-Nordlund theory are $x\approx 2.3$ for solenoidal forcing and $x\approx 2.7$ for compressive forcing.

\begin{figure}
\par\subfigure[Thermal support, $\lambda_{{\rm J}}^{0}/L=0.04$ , ...
...BE}}=10^{3}$ .]{\includegraphics[width=9cm,clip]{13904fg6d.eps} }
\end{figure} Figure 6:

Comparison of the semi-analytic mass distribution following from the Hennebelle-Chabrier (HC) theory and the Padoan-Nordlund (PN) theory with the corresponding distributions obtained via clump-finding (large dots) for two different choices of global mass scale. The thin dotted lines are the tangents to the mass distributions with the Salpeter slope x=1.35. For an explanation of turbulent support, see Sects. 2.2 and 3.

Open with DEXTER
\begin{figure}
\par\subfigure[$\lambda_{{\rm J}}^{0}/L=0.04$ , $N_{{\rm BE}}=2\t...
...]{\includegraphics[width=9cm,clip]{13904fg7b.eps} }
\vspace*{-2mm}
\end{figure} Figure 7:

Normalised gravitationally unstable mass as a function of the logarithmic density fluctuation for purely thermal support (thick dashed lines) and with additional turbulent support (solid lines), as defined by Eqs. (8) and (10). For solenoidal forcing (see Table 4), the mass corresponding to the scale for which thermal pressure equals turbulent pressure ( $\mathcal{M}_{\ast}\tilde{R}^{\eta}=1$) is indicated by the dot-dashed horizontal lines, and the dashed horizontal lines specify the mass for which the core scale equals the box size ( $\tilde{R}=N_{{\rm BE}}^{1/3}$).

Open with DEXTER

Figure 6 shows also the results from the clump-finding analysis. Let us first consider the case $\lambda _{{\rm J}}^{0}/L=0.04$ ( $N_{{\rm BE}}=2\times 10^{5}$). For compressively driven turbulence, the clump-finding data match the theoretical predictions very poorly. As discussed in Sect. 4.2, the size of the smallest cores is less than the grid resolution in this case. Since the numerically unresolved cores are missing in the core statistics, the distribution is biased towards higher masses in comparison to the semi-analytic distributions. For solenoidal forcing, on the other hand, the smallest cores are at least marginally resolved (see Table 4), and the low-mass wing, as well as the peak position following from the clump-finding analysis, is reasonably close to the theoretical predictions. While the distribution obtained by clump-finding roughly agrees with the CMF following from the Padoan-Nordlund theory also for $\tilde{M}\gg 1$, the slope of $x\approx 3.1$ (see Table 2) is significantly steeper than the theoretical value 2.3. In Sect. 5.2, it is shown that this discrepancy can be resolved by including turbulent pressure. Remarkably, the mass distributions from the clump-finding analysis are matched quite well by distributions of the form (9) for $\lambda _{{\rm J}}^{0}/L=0.2$ ( $N_{{\rm BE}}=10^{3}$). There are small deviations in the high-mass tails, which are, however, well within the error bars (see Fig. 5). In comparison to the distribution (4) predicted by the Padoan-Nordlund theory, large discrepancies become apparent. Since the cores are resolved enough both for solenoidal and for compressive forcing for $\lambda _{{\rm J}}^{0}/L=0.2$ (see Table 4), the most likely explanation is that $N_{{\rm BE}}$ is too small and, consequently, the high-mass cores are not within the asymptotic regime to which the Padoan-Nordlund theory applies.

Table 4:   Dependence of various parameters on the mass scale and the forcing.

5.2 Turbulent support

Since Padoan & Nordlund (2002) consider only the thermal Jeans mass, we concentrate on the Hennebelle-Chabrier theory for the case that includes both thermal and turbulent support. Eliminating $\tilde{R}$ from Eqs. (10) by means of numerical root finding yields the mass-density relations plotted in Fig. 7. For comparison, also the relation (8) for purely thermal support is shown. One can see that low-density cores of size R greater than $\lambda_{{\rm J}}^{0}$ can maintain a much higher mass against gravitational collapse if turbulent pressure is included in the Jeans criterion. The high-density asymptote, on the other hand, coincides with the thermally supported branch, because cores of high enough density are associated with small length scales R, for which the turbulent velocity dispersion becomes negligible compared to the speed of sound, i.e., $\mathcal{M}_{\ast}\tilde{R}^{\eta}\ll 1$. This implies $\tilde{M}\simeq\tilde{R}\simeq \exp(-\delta/2)$. For the length scale $R_{{\rm s}}:=L(\mathcal{M}_{{\rm rms}}/\sqrt{3})^{-1/\eta}$, we have $\mathcal{M}_{\ast}\tilde{R}_{{\rm s}}^{\eta}=1$, i.e., the turbulent pressure of the gas equals its thermal pressure. Substituting the values of the RMS Mach numbers and the scaling exponents from Federrath et al. (2010b), the values of $\mathcal{M}_{{\rm\ast}}$ listed in Table 4 are obtained. It follows that $R_{{\rm s}}\approx 0.074L$ and 0.082L for solenoidal and compressive forcing, respectively. These scales are close to the sonic length scales $\lambda_{{\rm s}}\approx 0.077L$ and 0.074L, which Federrath et al. (2010b) determined from the turbulence energy spectra. The corresponding dimensionless masses, $\tilde{M}_{{\rm s}}$, are also listed in Table 4[*].

An important limitation of calculations based on the approximation (7) is that R has to be small compared to the box size 2L, i.e., $\tilde{R}\ll N_{{\rm BE}}^{1/3}$. The mass parameters $\tilde{M}_{R~=~2L}$ corresponding to $\tilde{R}=N_{{\rm BE}}^{1/3}$ are indicated in Fig. 7. The values of $\tilde{M}_{R~=~2L}$ are listed in Table 4 for all parameter sets. As one can see, we have the constraints $\tilde{M}\ll 50$ and $\tilde{M}\ll 250$ for $\lambda _{{\rm J}}^{0}/L=0.2$ and 0.04, respectively. Since these limits are much higher than the peak positions of the mass distributions, there is a sufficient margin to study the high-mass wings.

The PDF data for the density fluctuation $\delta$ yield the mass distributions plotted in the bottom panels (c,d) of Fig. 6. Comparing the mass distributions with and without turbulent support, the high-mass tails are flatter in the former. As expected, the difference is more pronounced for $\lambda _{{\rm J}}^{0}/L=0.2$, because of the higher contribution of turbulent pressure for large core masses (see Fig. 7). Remarkably, the tail of the mass distribution for compressively driven turbulence that is plotted in Fig. 6d is much less stiff than what PN07 report for hydrodynamic turbulence and is in very good agreement with the Salpeter power law. The peak positions, on the other hand, are nearly unaffected, because the turbulent pressure on the length scales corresponding to the peaks of the mass distributions is low compared to the thermal pressure.

Regarding the mass distributions obtained by clump-finding with turbulent support (bottom panels of Fig. 6), similar trends to the distributions with purely thermal support can be seen. For $\lambda _{{\rm J}}^{0}/L=0.04$ ( $N_{{\rm BE}}\approx 1.2\times 10^{5}$), the clump-finding distribution is shifted towards higher masses for compressively driven turbulence because the smallest cores cannot be resolved (see Sect. 5.1). The discrepancy is much less in the case of solenoidal forcing, for which the overall shape of the semi-analytic and clump-finding distributions agree quite closely. However, the slopes of the tails following from the clump-finding data are steeper in both cases. Consequently, it appears that the clump-finding algorithm underestimates the turbulent support of high-mass cores, provided that the theoretical description of the CMF is correct.

In the case $\lambda _{{\rm J}}^{0}/L=0.2$ ( $N_{{\rm BE}}=10^{3}$), the mass distributions agree for compressive forcing (except for a small shift), whereas the mass distribution obtained by clump-finding is markedly different from the semi-analytic model for solenoidal forcing. The analysis in Sect 4.2 showed that the cores tend to be smaller for compressively driven turbulence. In accordance with the clump-finding results, the Hennebelle-Chabrier theory implies a significant flattening of the high-mass wing in this case, and the slope is found to be close to the Salpeter value. For solenoidal forcing, however, a significant number of cores have sizes that are of the same order of magnitude as the integral scale L (see Table 3 and Fig. 4). This entails a violation of the basic assumption that was made in the derivation of Eq. (12). Apart from that, it is important to realize that both the clump-finding algorithm and the semi-analytic approach rely on the notions of the Bonnor-Ebert or Jeans mass and turbulence pressure. Since it is based on the collapse of a spherical cloud, the former is highly idealised and does not apply to the fully non-linear regime, while turbulence pressure is an ensemble-average property of isotropic turbulence. One mechanism that might account for the massive cores is merging; indeed, since the gas is less compressed, the cores are more extended and should merge more easily than in the compressible forcing case. Also they should live longer, as it takes more time for them to be assembled, and it is likely that these cores are destroyed by violent passing shocks.

6 Conclusions

After omputing the distributions of gravitationally unstable cores by means of a clump-finding algorithm as well as semi-analytic methods for hydrodynamic simulations of forced supersonic turbulence, we have found significant differences resulting from the turbulence forcing, the choice of the global mass scale, and the effects of turbulent support. For a comparison with theoretical predictions, the numerical probability density functions of the mass density fluctuations were used to evaluate the formulae for the mass distribution. Our analysis completes a comprehensive study of the influence of forcing on the density and velocity statistics of isothermal supersonic turbulence (Schmidt et al. 2008; Federrath et al. 2010b,2009,2008b). In the following, we summarise the main results:

1.
For hydrodynamic simulations without explicit treatment of self-gravity, the CMF depends on the choice of the global mass scale, which determines the number of Bonnor-Ebert masses, $N_{{\rm BE}}$, with respect to the mean density in the computational domain. For clump-finding, one has to observe two constraints: $N_{{\rm BE}}$ must be large enough to allow for many cores, but if $N_{{\rm BE}}$ is chosen too large, a significant fraction of cores will be numerically unresolved and the resulting CMF will be shifted towards higher masses. In our analysis, these constraints are satisfied for two cases: solenoidal forcing with $N_{{\rm BE}}\approx 1.2\times 10^{5}$ and compressive forcing with $N_{{\rm BE}}=10^{3}$. The lower value of $N_{{\rm BE}}$ is consistent with Larson-type relations.

2.
Comparing the results for solenoidal and compressive forcing, the most noticeable trends are that the CMF for compressively driven turbulence displays more low-mass cores, while the high-mass tails are flatter, particularly if turbulent support is significant. These trends are implied by Eqs. (8) and (10) for the stronger density contrast produced by compressive forcing (also see Fig. 7). We also found that the CMF peaks at lower mass in the case of compressive forcing. These results follow both from the clump-finding analysis and the semi-analytic computation of the mass distributions.

3.
Because of the considerable scatter of core masses and possible biases of the clump-finding algorithm, it is difficult to ascertain power laws for the mass distributions. Nevertheless, we attempted to determine power-law exponents for the high-mass tails from least-square fits and by means of MLE. The results agree within the statistical uncertainties for both methods. However, the results for purely thermal support are at odds with the theory of Padoan & Nordlund (2002), which asymptotically applies to $N_{{\rm BE}}\rightarrow\infty$ and predicts power-law tails for the CMF. In particular, the difference between the exponents following from the clump-finding analysis is more pronounced than the theoretical prediction.

4.
The mass distributions for cores with purely thermal support agree well with the theory of Hennebelle & Chabrier (2008) for $N_{{\rm BE}}=10^{3}$, although the tails are less stiff towards high masses than what is expected on the basis of the semi-analytic models. There are large discrepancies for $N_{{\rm BE}}=1.2\times 10^{5}$, but the disagreement might be spurious because of the strong dependence of the high-mass wings on the numerical resolution. Apparently, the mass distributions for high values of $N_{{\rm BE}}$ are in closer agreement with Padoan & Nordlund (2002) if purely thermal support is assumed. However, it should be noted that the mass-density pdfs will be strongly affected by self-gravity in this case.

5.
For $N_{{\rm BE}}=10^{3}$, turbulent pressure is important for a wider range of core masses than for $N_{{\rm BE}}=1.2\times 10^{5}$. For this reason, the CMFs change substantially if turbulent support is included. The effect is particularly strong for turbulence that is driven by solenoidal forcing, because the mass and the size of virtually all cores is increased by turbulent support. However, the mass distribution obtained by clump-finding cannot be reproduced theoretically, because the basic assumption that the size of the cores is much smaller than the integral scale is clearly not fulfilled in this case. Moreover, a stability criterion that is based on the notions of Jeans mass and turbulent pressure might fail if the cores become too large. On the other hand, we find very good agreement between the clump-finding analysis and the Hennebelle-Chabrier for compressively driven turbulence with $N_{{\rm BE}}=10^{3}$. Compared to purely thermal support, the high-mass tail is considerably flatter in this case.
One of the difficulties in assessing the applicability of the semi-analytic models is that we cannot disentangle limitations of the underlying theoretical concepts from shortcomings of the clump-finding algorithm. When studying the influence of the clump-finding parameters, it becomes clear that the shape of the resulting mass distribution varies significantly with these parameters. While the peak position appears to be quite robust, the tails are more affected by the tuning of the clump-finding algorithm. Apart from these uncertainties and the large scatter of the core masses, there might be systematic biases. In particular, it is not entirely clear whether clump-finding traces down the smallest gravitationally unstable cores.

Even so, we have been able to shed more light on the gravitational fragmentation of turbulent gas. The influence of the forcing is irrefutable for a range of length scales that certainly extends beyond the energy-containing range. An open question is which regime in the interstellar medium is suitably modelled by a particular mode of forcing. If compressive excitation of turbulence occurs on length scales that tend to be larger than the size of molecular clouds, then the core statistics on much smaller scales, i.e., inside the cores of the clouds, still might be universal. Federrath et al. (2010b), however, suggests that different forcing mechanisms can produce genuine differences in molecular clouds, which affect the properties of turbulence even on the length scales of molecular cloud cores. This calls for further advances from both observational and the theoretical directions.

Another important result is the influence of turbulent pressure on the gas fragmentation. This is not new at all from the theoretical point of view. But we have been able to demonstrate numerically that turbulent support entails a flattening of the CMF towards high masses. In one extreme case, the slope of the high-mass tail was found to be close to the Salpeter slope. In this respect, the effect of turbulent pressure on gravitational fragmentation is analogous to the effect of magnetic pressure in magnetohydrodynamical turbulence. This can be understood on the basis of the dispersion relation resulting from a linear stability analysis, which shows that the magnetic pressure, as well as the turbulent pressure, contributes to the effective pressure of the gas. We do not suggest that the observed CMF might be explicable in terms of hydrodynamical turbulence alone, because it is known that magnetic fields play an important role in the interstellar medium. However, our results show that turbulent pressure might significantly contribute to the slope of the high-mass tail of the CMF. Apart from that, this effect is important in numerical simulations, where unresolved turbulent velocity fluctuations are close to the speed of sound or higher. In this case, the corresponding turbulent pressure can be computed with a subgrid-scale model (Schmidt 2009).

To achieve further progress with the theoretical explanation of the CMF, it is paramount to account for processes that modify the PDFs of the mass density. One such process is the thermal instability induced by radiative cooling in the interstellar medium. It has already been noticed by Passot & Vázquez-Semadeni (1998) that polytropic equations of state, which are simple models for cooling, lead to power-law tails in the density PDFs. This is confirmed by Li et al. (2003) for three-dimensional turbulence-in-a-box simulations and by Audit & Hennebelle (2009) for three-dimensional simulations of colliding flows with a polytropic exponent $\gamma=0.7$, although a power law was not obtained for thermally bistable turbulence. Seifried et al. (2010), on the other hand, find power-law tails for thermally bistable turbulence produced by compressive forcing, using the cooling function of Audit & Hennebelle (2005). The effects of a polytropic equation of state on the CMF have been analysed by Hennebelle & Chabrier (2009).

Another process that gives rise to density PDFs with power-law tails is the accretion and contraction of gas in self-gravitating turbulence (e.g. Klessen et al. 2000; Federrath et al. 2008a). Moreover, dense structures can merge because of their mutual gravitational attraction. Consequently, the statistical characteristics will evolve with time (Klessen 2000), and the corresponding core mass distribution can differ significantly from the non-selfgravitating case (Klessen 2001), depending on the relative strength of selfgravity and the evolutionary stage of the star-forming cloud. Indeed, power-law, high-density tails are observed in high-resolution extinction maps of nearby molecular clouds (Kainulainen et al. 2009).

In numerical simulations of self-gravitating turbulence, the mass distribution of gravitationally collapsing objects can be determined directly. However, the sink particles that are usually applied to numerically capture the collapsing gas (Padoan & Nordlund 2009; Federrath et al. 2010a; Bate et al. 1995; Jappsen et al. 2005; Krumholz et al. 2004) do not directly correspond to the cores we have considered here. While sink particles are dynamic objects that accrete gas and thus can acquire different masses, they are obviously not associated with a variable length scale, as is the case for cores[*]. Theoretically, the relation between length scale and mass becomes apparent in the Hennebelle-Chabrier theory, but also the clump-finding algorithm identifies regions of variable size as cores. For this reason, relating the mass distributions of sink particles and cores is not trivial. Nevertheless, the approach via sink particles has the considerable advantage that concepts such as the Bonner-Ebert mass, which is rather arbitrary in the highly non-linear regime, are not required. In either case, including self-gravity in turbulence simulations is indispensable to improve our understanding of gravoturbulent fragmentation.

Acknowledgements
We thank Paolo Padoan (UCSD) for making his clump-finding algorithm available for this work and for critical comments. We are grateful to Patrick Hennebelle (ENS Paris) for numerous discussions and plenty of advice. Thanks also to Jens C. Niemeyer (IAG) for his support. The simulations and data analysis used resources from HLRBII project h0972 at the Leibniz Supercomputer Centre in Garching, Germany. CF is grateful for financial support from the International Max Planck Research School for Astronomy and Cosmic Physics (IMPRS-A) and the Heidelberg Graduate School of Fundamental Physics (HGSFP), which is funded by the Excellence Initiative of the Deutsche Forschungsgemeinschaft (DFG) GSC 129/1. R.S.K. acknowledges financial resources provided by the Deutsche Forschungsgemeinschaft under grants no. KL 1358/1, KL 1358/4, KL 1359/5. R.S.K. furthermore is grateful for subsidies from a Frontier grant of Heidelberg University sponsored by the German Excellence Initiative and for support from the Landesstiftung Baden-Württemberg via their programme International Collaboration II, grant no. P-LS-SPII/18. R.S.K. also acknowledges financial support from the German Bundesministerium für Bildung und Forschung via the ASTRONET project STAR FORMAT (grant 05A09VHA).

Appendix A: Maximum likelihood estimation for power-law scaling parameters

Estimating the parameters $\alpha$ and $x_{{\rm min}}$ of a power-law distribution of the form

\begin{displaymath}
d(x)=C x^{-\alpha},\;\;\;x\geq x_{{\rm min}}
\end{displaymath} (A.1)

with a constant C and exponent $\alpha$ by simply fitting straight lines on a log-log plot introduces strong biases (e.g. Clauset et al. 2007; Newman 2005). Furthermore, this technique gives no information whether assuming a power law is reasonable or not. To accurately estimate the scaling parameters of the high-mass range in our mass distributions we followed the MLE method of Clauset et al. (2007). In this appendix we give a brief overview of the method. For a detailed description please refer to the original paper.

Clauset et al. (2007) use a maximum likelihood estimator

\begin{displaymath}
\hat{\alpha}=1+n\left[\sum_{i=1}^{n}{\ln \frac{x_{i}}{x_{{\rm min}}}}\right]^{-1}
\end{displaymath} (A.2)

for the power-law exponent $\alpha$ assuming that the dataset, consisting of n values xi with $i=1 \ldots n$ is drawn from a continuous power-law distribution for $x_{{\rm i}}\geq x_{{\rm min}}$, and $\hat{\alpha}$ is the estimated value for $\alpha$ of the underlying distribution (Eq. (A.1)). The value of the lower bound $x_{{\rm min}}$ of the distribution is the crucial factor in this calculation. If the estimated $\hat{x}_{{\rm min}}$ is too small ( $\hat{x}_{{\rm min}} < x_{{\rm min}}$), we fit a power-law to non-power-law data and the parameter estimation will be biased. Too high values ( $\hat{x}_{{\rm min}} > x_{{\rm min}}$) will throw out a certain amount of legitimate, power-law-distributed data points. A reasonable value for $\hat{x}_{{\rm min}}$ minimises the Kolmogorov-Smirnov (KS) statistics defined by

$\displaystyle D= \max\limits_{x \geq x_{{\rm min}}}\ \vert S(x)-P(x)\vert,$   (A.3)

where S(x) is the cumulative distribution function (CDF) of the observed data points and P(x) the CDF of the fitted power-law model. The optimal value for $\hat{x}_{{\rm min}}$ gives the minimum difference between the data and the fitted model, which is just the minimised KS-statistics D of Eq. (A.3).

To check that the observed dataset is likely to be drawn from a power-law distribution, we calculate a p-value using a semi-parametric bootstrap approach. Therefore, we create 2500 synthetic datasets per snapshot, which follow a real power law with the estimated parameters $\hat{x}_{{\rm min}}$ and $\hat{\alpha}$ of the observed core mass values. For $x_{{\rm i}} < \hat{x}_{{\rm min}}$, the synthetic distributions follow the same non-power law behaviour as the observed dataset. For each of the synthetic distributions, we again estimate the parameters $\alpha$ and $x_{{\rm min}}$ and compute the KS-statistics. The p-value is the fraction of the KS-statistics of the synthetic datasets whose value is higher than the KS-statistics of the observed dataset. This means that p is the probability of getting a goodness of fit, e.g. the KS-statistics, for a real power-law distributed synthetic dataset that is at least as bad as the goodness of fit of our observed dataset. For $p \leq 0.1$ the power-law assumption has to be ruled out. A p-value higher than 0.1 does not necessarily mean that the underlying distribution follows a power law. The p-value can only rule out the hypothesis of a power-law distribution.

References

Footnotes

... follows[*]
In this article, $\log$ denotes the natural logarithm.
...[*]
These values implicitly depend on the integral scale L, which is set to half the box size (see Sect. 1). More accurately, L can be determined from the velocity structure functions (Kritsuk 2010, private communication). Indeed, this would entail a somewhat closer match between the semi-analytic core mass distributions and the clump-finding data. Since the trends remain unaltered, however, we keep the simple definition of L in this article.
... cores[*]
There is a fixed accretion radius extending over a few grid cells, which can be interpreted as the length scale associated with sink particles.

All Tables

Table 1:   Time-averaged properties of the mass distributions for solenoidal and compressive forcing calculated with the physical scalings applied in PN07.

Table 2:   Least-square estimates of the power-law exponent $x_{{\rm LS}}$obtained from the 10243 simulations for different values of $N_{{\rm BE}}$ with/without turbulent pressure included in the core stability criterion.

Table 3:   Time-averaged properties of the core size PDFs for solenoidal and compressive forcing.

Table 4:   Dependence of various parameters on the mass scale and the forcing.

All Figures

  \begin{figure}
\par\subfigure[Solenoidal forcing.]{\includegraphics[width=9cm,cl...
...{\includegraphics[width=9cm,clip]{13904fg1b.eps} }
\vspace*{2.8mm}
\end{figure} Figure 1:

The core mass distributions for the 5123 solenoidal a) and compressive b) simulations as a function of the clump-finding algorithm parameter f (increasing from top curve to bottom curve) which sets the relative spacing between two adjacent density levels. Error bars contain the 1$\sigma $ temporal fluctuations and are only indicated for f=1.04 for the sake of clarity.

Open with DEXTER
In the text

    \begin{figure}
\par\subfigure[Solenoidal forcing, $N_{{\rm BE}}=1.2\times10^{5}$...
...]{\includegraphics[width=9cm,clip]{13904fg2d.eps} }
\vspace*{2.5mm}
\end{figure} Figure 2:

Core mass distributions for solenoidal (left) and compressive forcing (right) at numerical resolutions of 2563 (dashed), 5123 (dot-dashed), and 10243 (solid), normalised to the total number of cores, $N_{{\rm tot}}$. Error bars indicate 1$\sigma $ temporal fluctuations.

Open with DEXTER

In the text

    \begin{figure}
\subfigure[Solenoidal forcing, $N_{{\rm BE}}=1.2\times10^{5}$ , t...
...t.]{\includegraphics[width=9cm,clip]{13904fg3b.eps} }
\vspace*{2mm}
\end{figure} Figure 3:

PDF of core size r defined in Eq. (18) for solenoidal a) and compressive b) forcing for a numerical resolution of 2563 (dashed), 5123 (dot-dashed), and 10243 (solid) and PN07 scaling. Error bars indicate 1$\sigma $ temporal fluctuations of the PDFs.

Open with DEXTER

In the text

      \begin{figure}
\par\subfigure[Solenoidal forcing, $N_{{\rm BE}}=10^{3}$ , therma...
...{\includegraphics[width=8.9cm,clip]{13904fg4d.eps} }
\vspace*{-1mm}
\end{figure} Figure 4:

PDFs of core size distribution as in Fig. 3 but for $N_{{\rm BE}}=10^{3}$ with/without turbulent pressure support (see Sect. 2.2 and 3).

Open with DEXTER

In the text

      \begin{figure}
\par\subfigure[$N_{{\rm BE}}=1.2\times10^{5}$ , thermal support.]...
...bulent support.]{\includegraphics[width=9cm,clip]{13904fg5d.eps} }\end{figure} Figure 5:

CMFs of compressive (solid) and solenoidal (dot-dashed) forcing for different values of $N_{{\rm BE}}$ with/without turbulent support (see Sects. 2.2 and 3) and a grid resolution of 10243. The least-square fits to the high-mass tails are shown as dashed lines.

Open with DEXTER

In the text

  \begin{figure}
\par\subfigure[Thermal support, $\lambda_{{\rm J}}^{0}/L=0.04$ , ...
...BE}}=10^{3}$ .]{\includegraphics[width=9cm,clip]{13904fg6d.eps} }
\end{figure} Figure 6:

Comparison of the semi-analytic mass distribution following from the Hennebelle-Chabrier (HC) theory and the Padoan-Nordlund (PN) theory with the corresponding distributions obtained via clump-finding (large dots) for two different choices of global mass scale. The thin dotted lines are the tangents to the mass distributions with the Salpeter slope x=1.35. For an explanation of turbulent support, see Sects. 2.2 and 3.

Open with DEXTER
In the text

  \begin{figure}
\par\subfigure[$\lambda_{{\rm J}}^{0}/L=0.04$ , $N_{{\rm BE}}=2\t...
...]{\includegraphics[width=9cm,clip]{13904fg7b.eps} }
\vspace*{-2mm}
\end{figure} Figure 7:

Normalised gravitationally unstable mass as a function of the logarithmic density fluctuation for purely thermal support (thick dashed lines) and with additional turbulent support (solid lines), as defined by Eqs. (8) and (10). For solenoidal forcing (see Table 4), the mass corresponding to the scale for which thermal pressure equals turbulent pressure ( $\mathcal{M}_{\ast}\tilde{R}^{\eta}=1$) is indicated by the dot-dashed horizontal lines, and the dashed horizontal lines specify the mass for which the core scale equals the box size ( $\tilde{R}=N_{{\rm BE}}^{1/3}$).

Open with DEXTER
In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.