EDP Sciences
Free Access
Issue
A&A
Volume 510, February 2010
Article Number A110
Number of page(s) 12
Section Interstellar and circumstellar matter
DOI https://doi.org/10.1051/0004-6361/200912236
Published online 18 February 2010
A&A 510, A110 (2010)

The thermodynamics of molecular cloud fragmentation

Star formation under non-Milky Way conditions

S. Hocuk - M. Spaans

Kapteyn Astronomical Institute, University of Groningen, PO Box 800, 9700 AV Groningen, The Netherlands

Received 31 March 2009 / Accepted 23 November 2009

Abstract
Context. Properties of candidate stars, forming out of molecular clouds, depend on the ambient conditions of the parent cloud. We present a series of 2D and 3D simulations of fragmentation of molecular clouds in starburst regions, as well as of clouds under conditions in dwarf galaxies, leading to the formation of protostellar cores.
Aims. We explore in particular the metallicity dependence of molecular cloud fragmentation and the possible variations in the dense core mass function, as the expression of a multi-phase ISM, due to dynamic and thermodynamic effects in starburst and metal-poor environments.
Methods. The adaptive mesh refinement code FLASH is used to study the level of fragmentation during the collapse. With this code, including self-gravity, thermal balance, turbulence, and shocks, collapse is simulated with four different metallicity-dependent cooling functions. Turbulent and rotational energies are considered as well. During the simulations, number densities of 10 $^8{-}10^9\rm ~cm^{-3}$ are reached. The influences of dust and cosmic ray heating are investigated and compared to isothermal cases.
Results. The results indicate that fragmentation increases with metallicity, while cosmic ray and gas-grain collisional heating counteract this. We also find that modest rotation and turbulence can affect the cloud evolution as far as fragmentation is concerned. In this light, we conclude that radiative feedback in starburst regions will inhibit fragmentation, while low-metallicity dwarfs should also enjoy a star formation mode in which fragmentation is suppressed.

Key words: stars: formation - hydrodynamics - dust, extinction - equation of state - turbulence - ISM: clouds

1 Introduction

The initial mass function (IMF) in our local neighborhood is observed to be a power-law function, nicely following a Salpeter slope. Recent studies (Elmegreen et al. 2008; Sabbi et al. 2007; Figer 2005a,b) confirm this universality at low and high masses. When we turn to observations of more radical environments like the Galaxy center, a ``universal'' shape is again confirmed and only hints of different IMFs appear. Measurements of abundance patterns (Cunha et al. 2008; Ballero et al. 2007) indicate a flatter IMF. Observations of the Arches cluster initially showed a shallow IMF (Stolte et al. 2002; Figer et al. 1999) and a turnover mass, hence a typical mass, at a few solar masses (6-7 $M_{\odot}$, Stolte et al. 2005). However, Kim et al. (2006) do detect low-mass stars with their deep photometry and merely find a local bump at $\sim$6.3 $~M_{\odot}$ and a somewhat shallower IMF ( $\Gamma = -1.0$ to -1.1). Dib et al. (2007) confirm these results and claim a top-heavy IMF. Again, there appear to be perfectly reasonable explanations for the higher mass stars detected in the Arches cluster, as Portegies Zwart et al. (2007) argue with their idea of a still collapsing cluster core. Of course, the idea that the conditions and the environment play an important role in shaping the IMF is worth studying. Looking at the mass-to-light ratios of ultra-compact dwarf galaxies, Dabringhausen et al. (2009) conclude from their models that it is most likely caused by a top-heavy stellar IMF. If so, then they attribute this to feedback effects induced by massive stars. Still, these studies do not benefit from resolved stellar population data. One can associate the conditions and the environment of dwarf galaxies with those in the early Universe (Salvadori & Ferrara 2009; Hunter 2008). It is also believed that the primordial, high-redshift, and zero-metallicity IMF must have been top-heavy, $\sim$100 $~M_{\odot}$ (Omukai & Palla 2001; Jappsen et al. 2009; Abel et al. 2002; Omukai & Palla 2003; Bromm et al. 2002; Abel et al. 2000). It thus appears that the IMF is universal, but that it is worthwhile to explore strong environmental (feedback) influences, such as starburst and dwarf galaxies, to understand the IMF better (Klessen et al. 2007). Since star formation, fragmentation, and the IMF all depend on initial conditions (Greif et al. 2008; Skillman et al. 2003), the question is how strong the impact of these initial conditions is. Also, how different should the environment be or how strong the feedback to cause a significant deviation from a Salpeter IMF?

Molecular clouds are the birth places of stars. These clouds are formed (rapidly) when interstellar gas and dust gather and cool down under the influence of hydrodynamics, thermodynamics, and gravity (Ballesteros-Paredes et al. 1999; Heitsch & Hartmann 2008; Hartmann et al. 2001; Klessen et al. 2005a; Dobbs 2008; Dobbs & Bonnell 2008). Molecules form in the gas phase and on dust grains and allow low temperatures to be reached (Omukai et al. 2005; Dulieu et al. 2010; Cazaux et al. 2005; Cazaux & Spaans 2009). Dense regions, like the center of the cloud, are able to form many stars, with a modest efficiency ($\sim$5%; Clark et al. 2005; Elmegreen 2000). The distribution of the dense regions is thought to be the precursor of a stellar initial mass function (Goodwin et al. 2008; Motte et al. 1998; Simpson et al. 2008; Motte & André 2001), and perhaps the stellar IMF is directly linked to this core mass function, as it seems for the Pipe nebula (Lada et al. 2008). Following a typical molecular cloud from its infancy to a mature state, for different starting conditions, gives insight into the processes that determine stellar masses and is the focus of this work. Specifically, we consider low-metallicity environments (dwarfs) and warm dusty systems (starbursts).

In the next section we present the initial conditions for our 38 simulations and present the results in the following section. Twenty of these simulations form the basis of this research. They comprise 4 different metallicities with 5 combinations of turbulent and rotational energies each. We add to these basis simulations comparison studies that include extra heating sources and cases where the primary heating and cooling terms are removed.

2 Simulations

2.1 Numerical model and simulation setup

To model a collapsing molecular cloud, we solved the hydrodynamical equations using the adaptive mesh refinement (AMR) code FLASH (Fryxell et al. 2000). FLASH uses the PARAMESH library (MacNeice et al. 2000; Olson et al. 1999) for grid refinement and parallelization. Hydrodynamic equations were solved using the piecewise parabolic method (PPM, Colella & Woodward 1984).

For our 2D simulations we used a grid with a maximum resolution that amounts to 40962 cells when completely refined, which is a refinement level of 10 in FLASH terms. To minimize computational demands, we utilized the power of adaptive meshes by initializing the simulation at refinement level 5, i.e., 1282 cells, and refine with a density threshold of $n = 600 \rm ~cm^{-3}$ in order not to violate the Truelove criteria (Truelove et al. 1997). This density was found to be the optimal balance between following the physics of cloud collapse and computational demand. The most dynamic ``central'' area, however, always started with the highest refinement to maximize precision and minimize resolution dependent effects in this region. The border resolution between the refinement levels was set-up in such a way that there is no sharp transition between the minimum and the maximum refinement.

We are aware of the problems that might arise if resolution is not sufficient. Truelove et al. (1997) state that the Jeans length should at least be resolved by 4 cells, i.e., Nj>4 in order to ensure that the collapse is of a physical rather than numerical nature. Another pitfall is that the transfer of gravitational energy to rotational energy can lead to unphysical values if the initial resolution is poor, i.e., the number of initial cells Ni<2145 (Commerçon et al. 2007). We are well beyond these limits. Our minimum number of cells, extending beyond the central region, is at least $N_i\geqslant16384$ cells. The Jeans length is resolved by more than Nj>77 cells. For even greater precision, we keep the tree structure of the top level hierarchy intact in FLASH. In other words, we initiate all simulations with one top-level block. Every simulation is continued up to twelve dynamical timescales, which is computed as

\begin{displaymath}t_{\rm dyn}=\sqrt{3\pi/32G\rho},
\end{displaymath} (1)

where G is the gravitational constant and $\rho$ is the mass density.

We do not include magnetohydrodynamic (MHD) effects in our models. Tilley & Pudritz (2007,2005) argue that these should not matter for the core formation as much as the turbulence does. However, MHD certainly is an important aspect, see Banerjee et al. (2007,2009), which we will investigate further in the future.

2.2 Initial conditions

We start the basis simulations with four metallicities $Z/Z_{\odot}=1,10^{-1},10^{-2},10^{-3}$ with respect to solar metallicity and take $Z=Z_{\odot }$ as a benchmark value. The dependence of the model outcomes on metallicity is the focus of this work.

For all simulations, we impose a Gaussian turbulence with a characteristic dispersion of $\sigma=2$ km s-1, comparable to observed values (Caselli et al. 2002; Falgarone et al. 2001) and on the order of 0.3 $\sqrt{GM/R}$. We did not assume the clouds to be rotationally supported but rather that turbulent motions, following the Larson laws, dominate. We adopted a Larson-like scale relation for the turbulence of the form $\sigma = 2~(R/ \rm 1~pc)^{0.5}~\rm km~s^{-1}$ (Larson 1981). Bonnell et al. (2006) also show that molecular clouds exhibit a $\sigma \propto R^{0.5}$ velocity dispersion law. The turbulence was introduced as part of the initial conditions of the cloud and not maintained by external forcing. The virialization of the gas that falls into local potential wells did maintain a level of super-thermal random motions for several dynamical timescales.

We ran each metallicity condition with five rotational energies, i.e., $\beta _0=10^{-1},10^{-2},10^{-3},10^{-4},0$, where $\beta _0$ is the initial ratio of rotational to gravitational energy. A Keplerian rotation was initiated around the z-axis multiplied by the ratio $\beta _0$. We define

\begin{displaymath}\beta_0 \equiv \frac{I\omega^2} {GM^2/r} = \frac{v^2 r} {GM}\cdot
\end{displaymath} (2)

Here, I is the moment of inertia, $\omega$ the angular velocity, v the radial velocity, and M and r are the total cloud mass and its corresponding radius. Studies show that a typical ratio for molecular cloud cores is $\beta_0=0.02$ (Caselli et al. 2002). These authors have observed 57 cloud cores and find energies within the range of $\beta_0 \approx 10^{-1}{-}10^{-4}$. Cosmological simulations indicate that star-forming host clouds have rotation energies on the order of $\beta_0 \simeq 0.1$ (Yoshida et al. 2006; Bromm et al. 2002).

We simulated an initially spherical cloud with a radius of 10 pc within a box that is periodic in both gravity and spatial boundaries. Given the larger box size of 26 pc and because of symmetric initial conditions, we neither expect nor see that the periodicity is an issue. We initialized the cloud with a mean constant density of $n\simeq10^3 \rm ~cm^{-3}$ (Frieswijk 2008), since it is not likely that a well-developed density profile exists before virialization. These initial conditions mean that we start with an initial cloud mass of $28.3\times 10^{3}$ solar masses. The initial gas temperature is put at 100 K throughout the cloud, simulating starburst regions (Ott et al. 2005; Spaans & Silk 2000).

2.3 The simulations

We present twenty simulations, with flat initial density profiles, where we focus on the effect of metallicity on cloud fragmentation. In a similar fashion to the work by Machida (2008) and Machida et al. (2009), we also investigate the importance of rotation and turbulence in these. We present 15 more simulations with the same setup to test and compare specific physical cases.

Metallicity:

four sets of five simulations are done to follow the collapse and fragmentation of metallicity-dependent cooling functions for turbulent, rotating, and nonrotating cases. These simulations have an initial cloud temperature of 100 K, which represents the high temperatures within the active environments where the molecular clouds are formed.

Dust and CRs in starbursts:

five solar metallicity runs with additional heating sources; cosmic rays and dust. In massive star-forming regions prominent in starbursts, there can be a lot of warm dust affecting the thermal balance of the cloud. Dust can act as a cooling mechanism when the gas temperature is higher than the dust temperature. Once the temperature of the gas is lower than the dust, this is the case for photon-dominated regions (PDRs) with AV > 5 mag, it will become a source of heating.

Cosmic rays can ionize and heat the gas, but their contribution is generally thought to be low. In starbursts, the impact of CR-heating is stronger because of the increased production by massive stars, i.e., supernovae.

Isothermal at 25 K:

five runs kept the temperature nearly constant at a starting cloud temperature of 25 K, representing an environment where the temperature is higher than what we observe in the Milky Way. For these simulations, the temperature is not allowed to go below 25 K and the compressional heating term is removed. The results are compared to cases where the thermal balance is incorporated.

Isothermal at 10 K:

a second series of near-isothermal runs are initiated at T=10 K in order to compare them to Milky Way conditions and to check the impact of the isothermal condition.

2.4 Chemistry, cooling and heating

In our simulations we used a detailed cooling function created with the Meijerink & Spaans (2005) code. The metallicity-dependent cooling includes fine-structure emissions from carbon and oxygen, as well as molecular lines from species like CO and H2O. The level populations were corrected for LTE effects above 103 cm-3. The chemistry includes gas phase and grain surface formation of H2 and HD (Cazaux & Spaans 2004,2009) and line trapping by the initial column ( $N_{\rm H}\sim \rm 10^{22.5}~cm^{-2}$ and $\sigma=2$ km s-1) through the cloud, as described in Meijerink & Spaans (2005) and using the multi-zone escape probability method of Poelman & Spaans (2005). Line trapping during cloud collapse was ignored (see the discussion). A background radiation field was included that conforms to 0.1 G0 (in units of the Habing field, Habing 1968) in the dwarfs and 30 G0 for the starbursts (Spaans & Norman 1997; Klessen et al. 2007), see Fig. 1.

\begin{figure}
\par\includegraphics[width=8cm,clip]{12236fg1.eps}
\end{figure} Figure 1:

Four cooling curves by Meijerink & Spaans corresponding to four metallicities $Z = Z_{\odot}, 10^{-1}~Z_{\odot}, 10^{-2}~Z_{\odot}$, and $10^{-3}~Z_{\odot}$. The cooling rates ( $\Lambda _{\rm cr}$) correspond to a number density of 104 cm-3.

Open with DEXTER

At high column densities, cosmic ray heating can play an important role. For our heating by cosmic rays, we used the definition of Meijerink & Spaans (2005) as presented in Eq. (3), which scales linearly with the density.

\begin{displaymath}\Gamma_{\rm CR} = 1.5 \times 10^{-11} ~\zeta n({\rm H}_2) \rm ~erg~cm^{-3}~s^{-1}.
\end{displaymath} (3)

We took the cosmic ray ionization rate per H2 molecule $\zeta$ as  $9\times 10^{-16}$ s-1. This corresponds to a starburst environment with a background star formation rate of $\sim$30 $M_{\odot}$ per year, with the appropriate scaling for dwarf galaxies (Spaans & Silk 2005).

\begin{figure}
\par\mbox{\includegraphics[width=7.5cm,clip]{12236fg2.eps}\hspace*{4mm}
\includegraphics[width=7.5cm,clip]{12236fg3.eps} }
\end{figure} Figure 2:

A 2D slice of a 3D cloud collapse simulation at $t \simeq t_{\rm dyn}$. Plotted variable is density. Left: a face-on view of the disk. Right: an edge-on view of the disk.

Open with DEXTER

Besides the heating by cosmic rays, heat can be transferred through collisions between dust and gas. The gas-dust temperature difference is important here. We add gas-grain collisional heating to our simulations as given by Hollenbach & McKee (1989,1979) and Meijerink & Spaans (2005);

                                $\displaystyle \Gamma_{\rm coll.}$ = $\displaystyle 1.2 \times 10^{-31} n^2 \left( \frac{T_{\rm k}} {1000} \right)^{\frac{1}{2}} \left( \frac{100\ {\rm\AA}} {a_{\rm min}} \right)^{\frac{1}{2}}$ (4)
    $\displaystyle \times \left[ 1 - 0.8\exp \left( \frac{-75} {T_{\rm k}} \right) \right] (T_{\rm d} - T_{\rm k})
\rm ~erg ~cm^{-3} ~s^{-1}.$  

In this equation, $T_{\rm k}$ is the gas kinetic temperature, n is the average number density and $a_{\rm min}$ is the minimum grain size that is taken as 10 Å. We keep the grain size distribution fixed (MRN, Mathis et al. 1977). The dust temperature we use, $T_{\rm d} = 39$ K (Wiklind & Henkel 2001), is for a starburst region (Meijerink & Spaans 2005) and is kept constant. The heat transfer between gas and grains will act as a cooling mechanism for the gas when it has a higher temperature than the dust.

Because line trapping beyond a column density 10 $^{22.5} \rm ~cm^{-2}$ is not taken into account, the temperatures presented in this paper are lower limits. Above $\sim$10 $^{4.5}~(Z_{\odot}/Z)$ particles per $\rm cm^{-3}$, heating by dust dominates. In this case, line trapping becomes less important because the n2 dependence of the gas-dust coupling prevails over the n dependence, or weaker, of line emission in LTE. The fine intricacies of cooling and heating together with gravity determine the main results of this paper.

3 Results

First, we performed three 3D simulations at refinement levels 9 and 10 (20483 and 40963 cells). We did these simulations to see if the extra dimension produces significant differences in the results as compared to 2D. Due to one less degree of freedom, 2D simulations can result in a higher merger rate and thus fewer fragments with higher masses in the end. The three 3D turbulent simulations were carried out for modest rotational energies, $\beta _0 = 10^{-1}$, $\beta_0 = 10^{-3}$, and one nonrotating initial condition, to avoid excessive influence by rotational motions on the outcomes. These runs are slightly lower in resolution than the 2D simulations, and we follow them up to 4  $t_{\rm dyn}$. In the two rotating cases, we see a disk forming within a time close to the free fall timescale $t \approx t_{\rm dyn}$, see Fig. 2. After $t_{\rm dyn}$, one effectively has a 2D system. The nonrotating simulation, however, does not exhibit significant flattening. Nevertheless, we do observe that the results, even for this case, is very similar to its 2D counterpart for solar metallicity. It has been verified that fragmentation starts on the same timescales as for the 2D simulations, $t \gtrapprox t_{\rm dyn}$. The number of fragments is somewhat higher, while their masses are comparable to their 2D counterparts. We emphasize that metallicity effects have our interest and feel that 2D simulations suffice to capture any salient features.

The final stage of every simulation is analyzed for the number of fragments (cores) and their masses. The results of all basis simulations in metallicity Z and initial turbulent and rotational energy are shown as density images after 12 theoretical dynamical timescales in Fig. 3.

We find the fragments by using an algorithm called isovolume within the visualization tool VisIt (Childs et al. 2005). A given density threshold value lets the algorithm select the connected components above this threshold. Putting a density cap at 105 cm-3 to select the dense cores (Alves et al. 2001; Caselli et al. 2008), we were able to count the number of fragments and determine their masses. The choice of threshold density matters for the number of fragments one derives and for their masses. However, we find that there is no significant change in the number of fragments for a spread of one order of magnitude in our threshold density, while the fragment masses do decrease with increasing threshold density. Clumps and fragments can also be identified from their deep potential wells like the bound p-cores method by Smith et al. (2009). The fragment sizes range between r=6000 AU and 84 000 AU. A complete overview of the fragments and their masses can be found in Table 1.

\begin{figure}
\par\includegraphics[width=16.5cm,clip]{12236fg4.ps}
\end{figure} Figure 3:

Density plots of the final stage ( $t=12t_{\rm dyn}$) of every simulation of metallicity against rotational energy. Orange to red depicts high density, typically higher than $\rm 10^{4}~cm^{-3}$ and up to $\rm 10^{9}~cm^{-3}$. Light blue portrays those densities lower than $\rm 1 ~cm^{-3}$. Top to bottom: decreasing in metallicity, $Z=Z_{\odot}, 10^{-1}~Z_{\odot}, 10^{-2}~Z_{\odot}, 10^{-3}~Z_{\odot}$. Left to right: decreasing in rotational energy $\beta _0=10^{-1},10^{-2},10^{-3},10^{-4},0$. All images have boxsizes of $4\times 4$ parsec.

Open with DEXTER

Table 1:   The number of fragments and their masses for each case of metallicity versus rotational energy.

\begin{figure}
\par\begin{tabular}{c c}
\thicklines
\includegraphics[scale=0.38...
...g7.ps} &
\par\includegraphics[scale=0.38]{12236fg8.ps}\end{tabular}
\end{figure} Figure 4:

Temperature plots for four different metallicities at $t=12t_{\rm dyn}$, $Z=Z_{\odot }$ ( top left), $10^{-1}~Z_{\odot}$ ( top right), $10^{-2}~Z_{\odot}$ ( bottom left), $10^{-3}~Z_{\odot}$ ( bottom right), for $\beta _0 = 10^{-1}$. The color red represents temperatures above 50 K. Contour plots of number density are overplotted in these figures in black. Contour levels are 10 6, 105, 104, 103, 102, 101 cm-3. The axes are in units of 10 $^{18} \rm ~cm$ in both directions. In the bottom left corner of each image, the Jeans length is illustrated by a white circle. To calculate the Jeans length, the average value of the temperatures and densities of the clumped regions of each image is used. These are: T=(10.0, 13.1 13.4, 19.0)[K] and $\rho =(7.0$, 3.3, 0.41, 0.22) $[\rm \times 10^{-19}~g~cm^{-3}]$.

Open with DEXTER

3.1 Metallicity and multi-phase structure

We evaluated the temperature inside the cores and that of the surrounding gas of every simulation. In most cases, the cooling function is able to cool the gas efficiently. Even though there is compressional heating, the gas can cool down to temperatures as low as 10 K. The difference lies in the rate with which each model cools down or heats up. For the various metallicities, we present temperature plots in Fig. 4. The images are representative for every value of $\beta _0$. The density is overplotted with contours in these images to highlight the high density regions, and specifically, to show the location of the dense cores. The Jeans length is illustrated on the lower left side of the images[*].

We find that the temperature increases inside the fragments as the metallicity decreases. The regions where the density is higher than its surroundings have temperatures of about 10 K on average for $Z=Z_{\odot }$ and 20 K for $Z=10^{-3}~Z_{\odot}$. We see the same trend for the lower density gas, in that the amount of cool gas, i.e. gas with temperatures lower than 50 K, decreases with Z. For a gas density above 10 $^{4} \rm ~cm^{-3}$ and a gas temperature of more than 100 K, fine-structure cooling of [OI] dominates with mid-J CO emission and water at the 20% level. Below 100 K, low-J CO line emission is the main coolant with contributions from water on the 10% level. Diffuse gas below 10 $^{4} \rm ~cm^{-3}$ is mainly cooled by [CII].

In Figs. 5 and 6 we present phase diagrams, gas temperature versus density, which highlight the response of the gas to metallicity driven cooling and additional heating processes. It is immediately obvious that a cold dense phase is supported, but also that more warm gas exists at lower metallicity or higher dust temperatures. In Fig. 5, one immediately notices that for high metallicities, the slope ${\rm d}\log(T)/{\rm d}\log(n)$ is steeper between $\sim$1 and $\sim$10 $^{3}{-}10^{4}~\rm cm^{-3}$ and that the cooling time is shorter. We also note that the maximum density achieved is lower at lower metallicities. The diagrams show that a two-phase medium exists with preferred densities of around 1 and 10 $^{4}~\rm cm^{-3}$. This is independent of the specific grid structure.

Above 10 $^{4}{-}10^{5}~\rm cm^{-3}$, compressional heating starts to dominate the thermal balance and raises the temperature. The increase due to adiabatic heating is weakened for higher metallicities, while the temperature can rise up to 1000 K for lower metallicities, if collapse is not halted completely. For low turbulent/rotational energies and in the case of zero rotation, there is a distinct curvature at a single density present. This is where the temperature rises quickly due to compression so that the inefficient cooling cannot compensate for this fast increase. At this point, the object is thermally supported until it reaches higher temperatures where the cooling gets increasingly more effective.

\begin{figure}
\begin{tabular}{c c}
\par\includegraphics{12236fg9-NEW.eps} &
\in...
...12236f15-NEW.eps} &
\includegraphics{12236f16-NEW.eps}\end{tabular}
\end{figure} Figure 5:

The phase diagrams of each metallicity after twelve dynamical timescales for two $\beta _0$ ratios, 10-1 ( left) and 10-3 ( right). Top to bottom, highest ( $Z=Z_{\odot }$) to lowest ( $Z=10^{-3}~Z_{\odot}$) metallicity.

Open with DEXTER

\begin{figure}
\par\includegraphics[width=8.5cm,clip]{12236f17-NEW.eps}\hspace*{...
...}\hspace*{4mm}
\includegraphics[width=8.5cm,clip]{12236f20-NEW.eps}
\end{figure} Figure 6:

The phase diagrams of several nonrotating simulations with additional physics for different environments after twelve dynamical timescales. The two figures on the right represent isothermal simulations at 25 K ( top) and 10 K ( bottom). The top left one is a phase diagram of the solar metallicity simulation with dust and CR heating. The bottom left figure is for a nonrotating 3D simulation also with $Z=Z_{\odot }$.

Open with DEXTER

There are always accretion shocks with high temperatures in their wake. This is important for the clouds' evolution and fragmentation. We find shocks in low-density gas and occasionally near the fragments. They consist of thin layers of slightly overdense gas with high temperature and a steep velocity gradient. Shock waves are more prominent when the collapse is rapid. Shock heating can dominate under such conditions.

3.2 The effect of metallicity, turbulence, and rotation on fragmentation

\begin{figure}
\par\includegraphics[width=7.8cm,clip]{12236f21.eps}\hspace*{4mm}
\includegraphics[width=8cm,clip]{12236f22.eps}
\end{figure} Figure 7:

Average fragment mass ( $\langle M_{\rm f} \rangle$) with initial conditions as determined after twelve dynamical times. Left: $\langle M_{\rm f} \rangle$ versus metallicity for each initial rotational energy. Right: $\langle M_{\rm f} \rangle$ versus rotational energy for each metallicity.

Open with DEXTER

We notice the immediate result that higher metallicity and higher rotational energy leads to more fragmentation, also found by Clark et al. (2008). A similar in-depth work on this has also been performed by Machida (2008) and Machida et al. (2009), at higher densities.

From a metallicity perspective,

a clear trend is visible, which at higher metallicities show more fragments and the average masses of the fragments are generally lower. The last statement is illustrated in Fig. 7. If we sum all the fragment masses, the total fragmented mass increases with the amount of fragmentation, hence with metallicity. This is caused by the fragments accreting mass over time and more fragments accreting more mass. The compressibility of the gas also depends on the metallicity. One may understand this by using a polytropic equation of state of the form $P\propto\rho^{\gamma}$, for an ideal gas $P\propto\rho T$. Increased metallicity then causes $\gamma$ (Eq. (5)) to be softened to unity or less since cooling is more efficient:

\begin{displaymath}\gamma = 1 + \frac{{\rm d}\log(T)} {{\rm d}\log(\rho)}\cdot
\end{displaymath} (5)

In a few low metallicity cases, we did not find fragments with densities above the threshold criterion. In the same way, the inability to fragment holds for all low metallicity runs. While turbulence initially causes density perturbations, and structures with marginally high density, the cloud does not develop into more compact objects. We see that these undeveloped structures are rather diffuse. The internal pressure is too high for them to collapse on their own and the cooling is not efficient enough to bring down their temperatures. The latter effect causes the thermal Jeans mass, given in Eq. (6), to become very high for these cases and prevents further collapse:

\begin{displaymath}M_J = \left( \frac{3} {4\pi} \right)^{\frac{1}{2}} \left( \fr...
...{\frac{3}{2}} \frac{T^{\frac{3}{2}}} {\rho^{\frac{1}{2}}}\cdot
\end{displaymath} (6)

In solar mass units, the Jeans equation can be formulated as (Frieswijk 2008)

\begin{displaymath}M_J \simeq 7.5 \left(\frac {T} {10~{\rm K}}\right) \left(\fra...
..._2}} {10^4~{\rm cm}^{-3}}\right)^{-\frac {1} {2}}~(M_{\odot}).
\end{displaymath} (7)

Typical Jeans masses are $\sim$10 solar masses for a 10-20 K object with a number density of $\sim$10 $^{5}~\rm cm^{-3}$. For low metallicities, where the temperature can rise to a few 100 K, the Jeans mass lies in the range $10^{2}{-}10^{3}~M_{\odot}$. Nonetheless, the Jeans mass can become as low as a few solar masses at higher densities despite the increase in gas temperature. A typical Jeans mass diagram of a solar metallicity simulation is shown in Fig. 8.

From a rotational perspective,

we see that faster rotation results in more fragments in general. We find fragmentation happening for all given initial rotational energies for the highest metallicity case, because of the imposed turbulence. Even at the lowest rotational energy, collapse to a binary still occurs. Since the turbulence imparts significant angular momentum, the $\beta=0$ case looks similar to the low $\beta$ cases, except that there is more diffuse material surrounding the cores. All simulations are supersonically turbulent over the first 1-2 dynamical timescales and level off to 0.4-2 km s-1 after 3 $t_{\rm dyn}$. Typically, the effective line widths ($\Delta v$) at 12 $t_{\rm dyn}$ are about 1 km s-1 on the scale of the fragment ($\sim$0.1 pc). Higher metallicity tends to lead to more turbulent motions. For the lower metallicity cases, we do not see fragmentation below $\omega_0 t_{\rm dyn} < 0.03$. Here, $\omega_0$ is the angular velocity, which is obtained just by putting $v = \omega_0 r$ into Eq. (2). Combining it with Eq. (1), we get Eqs. (8) (Matsumoto & Hanawa 2003) and (9). In fact, fragmentation only occurs in our studies when the criterion of Eq. (10) is met. This is found by fitting a function through the points where we still see fragmentation:

\begin{displaymath}\omega_0 = \sqrt{\frac{\beta_0 G M} {r^3}}
\end{displaymath} (8)

\begin{displaymath}\omega_0 t_{\rm dyn} = \sqrt{\frac{\pi^2 \beta_0} {8}}
\end{displaymath} (9)

\begin{displaymath}\omega_0 t_{\rm dyn} ~\rm exp \left( \frac{-(\log(Z/Z_{\odot}...
...} {1.74} \right) \sqrt{Z/Z_{\odot}} \approx 1.28\times10^{-2}.
\end{displaymath} (10)

Figure 7, right, shows that the average fragment mass $\langle M_{\rm f} \rangle$ is higher at lower initial rotational energies. We see here that the cloud is fragmenting and the average fragment mass is greatly reduced for lower rotational energies when the cloud has higher metallicity and is turbulent.

The specific angular momenta, j, of the fragments with size scales of 0.03-0.1 pc are typically a few times 10 $^{-21} \rm ~cm^{2}/s$. For larger fragments ($\sim$0.3 pc), j is somewhat higher and $\sim$10 $^{22} \rm ~cm^{2}/s$, in agreement with Bodenheimer (1995).

3.3 The effect of the environment on fragmentation

We present the density plots of the final stage of every run in Fig. 9. The results of these follow-up simulations can be found in Table 2, and the representative average fragment mass is plotted in Fig. 10.

Dust and CRs in starbursts:

when we add an extra heating source due to dust, i.e., gas-grain collisional heating, and one due to elevated cosmic ray ionization, we see noticeable changes in the end results. Initially, the cooling becomes more efficient, because the gas has a higher temperature than the dust so that the dust acts as a coolant for the gas. When the gas temperature drops and falls below the dust temperature, the gas is heated by the dust (and by cosmic rays). In the end, the fragment temperatures are close to the dust temperature of 39 K. The phase diagram for a nonrotating case in Fig. 6 (upper left) shows this steep decrease in temperature, slightly cooling below the dust temperature before it stabilizes around the dust temperature. The strong collisional coupling between the dust and the gas imposes this tight T-n relation.

Heating by cosmic rays and dust is strong enough to prevent further fragmentation on smaller scales and causes the largest fragments to accrete nearby material. Typically, fewer fragments are formed at this point, for the models with $\beta_0 \geq 10^{-2}$. There are more fragments for the case $\beta_0 = 10^{-4}$ where overall gas heating due to adiabatic compression and CRs renders $T_{\rm g} > T_{\rm d}$. This scenario is comparable to the situation for population II.5 stars (Schneider et al. 2006). Thus, dust can enhance fragmentation when the conditions are right.

A remarkable result that is seen when there is zero initial rotation is that the warm dusty cloud stops fragmenting in contrast to the expected trend. The would-be fragments are blurred out and settle in a disk around the sole compact object instead. We find that the cloud is very sensitive to initial rotation in this regime. A rapid turnover from a fragmenting cloud to a nonfragmenting case is found to be around $\beta_0 \sim 10^{-5}$. Interestingly, this is in close agreement with the lower limit of fragmentation ( $\beta_0 = 10^{-4}$) that Machida et al. (2009) find. Starburst systems may have super-solar metallicities (Israel 2009). We have verified that increasing the metallicity of the cloud to twice solar has little impact on our results.

Isothermal:

when we keep the temperature constant, we effectively assume that the pressure scales linearly with the density. Two near-isothermal runs are compared to the thermal balance runs at 10 K and 25 K. In environments that are less active, molecular clouds are known to have these low temperatures. Although we initialize with these temperatures, accretion shocks, which can increase the temperature up to a few hundred degrees in their wake, still affect the temperature of the cloud and the cores. Typically, the temperature remains constant around the initialized value for gas with density higher than $n=1 \rm ~cm^{-3}$ (see the two diagrams on the right in Fig. 6), while the temperature for the diffuse gas can go up one order of magnitude. For the 25 K condition, this makes the Jeans mass higher than simulations with incorporation of the thermal balance and where the cooling is dominant. Compared to the highest metallicity case $Z=Z_{\odot }$ of the metallicity study, for $\beta_{0}\ga 10^{-3}$, the fragments are less well developed and have densities similar to the lower metallicity results. Even for the 10 K isothermal case there are fewer fragments, though it has the same final temperature as the cooling included runs.

At lower rotational energies the difference is less apparent, since compressional heating for non-isothermal runs can overcome the cooling. In fact, for nonrotating clouds, fragmentation for the 10 K isothermal case tends to be higher. This change in fragmentation is caused by the temperature gradient with increasing density. Although the thermal Jeans mass of the metallicity runs is initially higher than in all the near-isothermal models, fragmentation is stimulated due to softening of $\gamma$. For the lowest rotational energy ( $\beta_0 = 10^{-4}$) and the nonrotating case in the solar metallicity model, $\gamma$ is stiff ( $\gamma=1.4$), and the temperature stays relatively high ($\ga$100 K) compared to the near-isothermal case, because compressional heating suppresses fragmentation. The evident difference between the two near-isothermal cases (10 and 25 K) can be explained by the difference in Jeans mass of a factor $\sim$4.

\begin{figure}
\par\includegraphics[width=8.5cm,clip]{12236f23-NEW.eps}
\end{figure} Figure 8:

A diagram that shows the Jeans mass as a function of number density for a nonrotating solar metallicity run at $t=12t_{\rm dyn}$.

Open with DEXTER

Table 2:   Analysis of the fragments for the environmental simulations.

4 Conclusions and discussion

We performed 38 hydrodynamical simulations where we rigorously tested the dependence of molecular cloud fragmentation on the ambient conditions as they pertain to starburst and dwarf galaxy regions. We compared these with the well known conditions of the Milky Way. Our results on the metallicity simulations show that the amount of fragmentation and the compressibility of the gas scales with increasing metallicity. In this, the cooling by gas phase and grain surface chemical species plays a subtle but fundamental role. We find that turbulence and rotationalso play an important role. Turbulence of 0.3 $E_{\rm grav}$, following the Larson laws, is a crucial ingredient in the development of fragments. Although rotation is observed for sub-pc scale molecular clouds, it is unlikely that clouds of these sizes (parsec) have rotational support. All simulations were therefore performed for $\beta _0$ values of much less than unity. These results agree with expectations and are comparable to the work of Machida (2008) and Machida et al. (2009), albeit at lower density. The metallicity has a greater impact on cloud fragmentation than $\beta _0$. For the highest metallicity model, fragmentation (into a binary) is still possible for even the lowest (or zero) rotational energy.

Regarding the initial temperature of the simulations, we note that starting with a high temperature of 100 K, as for starbursts, does not seem to be too important. When we do not consider the extra heating sources, radiative cooling can quickly overcome the high initial temperature and cool the cloud down to $\sim$10 K. This is true for most cases except when the initial rotational energy is very low. In this situation, the adiabatic heating and the heating due to shocks, which are more prominent now, can overcome the cooling. Fragments have higher temperatures >100 K under these circumstances, which suppresses fragmentation. The surrounding gas, though, stays cold ($\sim$10-20 K).

When we use low metallicity, similar to dwarf galaxies or the early universe, the average fragment mass increases in our simulations. This leads to massive cloud cores and can eventually induce massive star formation (Dabringhausen et al. 2009; Lada et al. 2008), relevant to the origin of the IMF.

\begin{figure}
\par\includegraphics[width=16.8cm,clip]{12236f24.ps}
\end{figure} Figure 9:

Density plots of the environmental simulations at $t=12t_{\rm dyn}$. Molecular clouds in environments having dust and cosmic rays are simulated with solar metallicity. Orange to red depicts densities typically higher than $\rm 10^{4}~cm^{-3}$ and up to $\rm 10^{9}~cm^{-3}$. Light blue portrays those densities lower than $\rm 1 ~cm^{-3}$. All images have boxsizes of $4\times 4$ parsec.

Open with DEXTER

We summarize our main conclusions as follows

  • Increased metallicity promotes fragmentation. There is a clear and strong dependence of fragmentation on metallicity, where the average fragment mass decreases with increasing metallicity.

  • Heating by cosmic rays and dust, for initial conditions of high-metallicity, are strong enough to slow down cooling and reduce fragmentation. When the gas and dust are thermally coupled, we find that the gas temperature stabilizes close to the dust temperature, $T_{\rm g} = T_{\rm d}$, at densities >10 $^5 \rm ~cm^{-3}$.

  • Fragmentation is suppressed for isothermal clouds. In almost all of the near-isothermal simulations that we did, we see a smaller number of fragments compared to their non-isothermal counterparts, except when adiabatic heating is very strong. Although the temperatures are exactly the same in several cases, even lower compared to the non-isothermal runs for the lowest $\beta _0$ ratios, and the Jeans masses are comparable, the evolution and fragmentation of the cloud differs. We attribute this to the stiffer effective equation of state, i.e. $\gamma=1$, which evidently has a strong impact (Klessen et al. 2005b).
For the main part in this research, we performed 2D simulations. One can argue that 2D simulations are not representative for a 3D cloud. Certain physics, like MHD, is indeed impossible to do in 2D, but we have not used any kind of physics or initial condition that requires a third dimension. We find that a higher merger rate occurs for 2D simulations thanks to one less degree of freedom and that our fragments have a somewhat higher mass. Still, we do not see many mergers in general. The 3D simulations that we ran show that a disk forms in a time that is comparable to or is less than the fragmentation time scale if even a modest rotation exists, effectively making the problem a 2D one. The results are slightly different, perhaps due to the lower resolution, or possibly owing to the shorter simulation time, but show comparable structures. A disk does not form when there is complete lack of initial rotation. This, however, does not change our results, since the nonrotating 3D simulation also fragments into two compact objects with comparable masses and enjoys a similar temperature evolution (Fig. 6) as its 2D counterpart.

Changes in the density profile do not affect the evolution of the cloud as far as the impact of metallicity is concerned. Testing a power-law density profile (Eq. (11)), as observed in active galactic regions (Mapelli et al. 2008; Genzel et al. 2003), against a flat profile and keeping the mean density constant, we find that the results are similar. The results (not shown) are the same for the number of fragments, but also alike in their temperature evolution. Only the fragment masses are marginally lower when using a flat profile, due to competitive accretion (Bonnell 2005):

\begin{displaymath}\rho = \rho_0 \left( \frac{1} {1 + R/R_{\rm c}} \right)^{\alpha} ~~\rm cm^{-3},
\end{displaymath} (11)

where R is the specific cloud radius, $\alpha=1.4$ the power of the profile, $\rho_0$ the central density of $n=10^4 \rm ~cm^{-3}$, and $R_{\rm c}$ the characteristic length of 1 pc. This is the radius within which the central density remains approximately constant.

\begin{figure}
\par\includegraphics[width=8.5cm,clip]{12236f25.eps}\end{figure} Figure 10:

Average fragment mass versus initial rotational energy for simulations within three distinct physical environments at $t=12t_{\rm dyn}$.

Open with DEXTER

We do not reach star densities, but the fragmentation and formation of dense cores is determined at the much earlier phases that we probe. We have not made use of sink particles, as this is usually done to make a jump in resolution and form protostars from dense cores, but our fragments should be able to form many stars. Opacity effects beyond 10 $^{22.5} \rm ~cm^{-2}$are not considered in our simulations simply because we do not reach very high densities. The maximum densities that we reach are $n\sim10^8{-}10^9$ cm-3. Molecular opacity at high densities can actually have quite a different effect on fragmentation compared to densities below 10 $^{3} \rm ~cm^{-3}$ (Poelman & Spaans 2005,2006). A cloud can become optically thick sooner and trap the lines when there are more metals to form molecules with. This in turn will heat the system and prevent further fragmentation. Incorporating radiative transfer will be the next step in the continuation of this research.

Acknowledgements
We would like to thank Floris van der Tak for useful discussions. MS thanks Henrik Beuther for a stimulating discussion. Special thanks go to Latif for all the long and fruitful conversations. We are very grateful to the anonymous referee and the editor M. Walmsley for insightful and constructive reports that greatly helped this work. The FLASH code was in part developed by the DOE-supported Alliance Center for Astrophysical Thermonuclear Flashes (ACS) at the University of Chicago. Our simulations were carried out on the Donald Smits Center for Information Technology using the HPC Cluster, University of Groningen.

References

Footnotes

... images[*]
The diameter of this circle corresponds to twice the Jeans length, which we refer to here as the Jeans circle. For an object to be able to fragment, it must have a size of at least twice the Jeans length and thus, must be larger than the Jeans circle.

All Tables

Table 1:   The number of fragments and their masses for each case of metallicity versus rotational energy.

Table 2:   Analysis of the fragments for the environmental simulations.

All Figures

  \begin{figure}
\par\includegraphics[width=8cm,clip]{12236fg1.eps}
\end{figure} Figure 1:

Four cooling curves by Meijerink & Spaans corresponding to four metallicities $Z = Z_{\odot}, 10^{-1}~Z_{\odot}, 10^{-2}~Z_{\odot}$, and $10^{-3}~Z_{\odot}$. The cooling rates ( $\Lambda _{\rm cr}$) correspond to a number density of 104 cm-3.

Open with DEXTER
In the text

  \begin{figure}
\par\mbox{\includegraphics[width=7.5cm,clip]{12236fg2.eps}\hspace*{4mm}
\includegraphics[width=7.5cm,clip]{12236fg3.eps} }
\end{figure} Figure 2:

A 2D slice of a 3D cloud collapse simulation at $t \simeq t_{\rm dyn}$. Plotted variable is density. Left: a face-on view of the disk. Right: an edge-on view of the disk.

Open with DEXTER
In the text

  \begin{figure}
\par\includegraphics[width=16.5cm,clip]{12236fg4.ps}
\end{figure} Figure 3:

Density plots of the final stage ( $t=12t_{\rm dyn}$) of every simulation of metallicity against rotational energy. Orange to red depicts high density, typically higher than $\rm 10^{4}~cm^{-3}$ and up to $\rm 10^{9}~cm^{-3}$. Light blue portrays those densities lower than $\rm 1 ~cm^{-3}$. Top to bottom: decreasing in metallicity, $Z=Z_{\odot}, 10^{-1}~Z_{\odot}, 10^{-2}~Z_{\odot}, 10^{-3}~Z_{\odot}$. Left to right: decreasing in rotational energy $\beta _0=10^{-1},10^{-2},10^{-3},10^{-4},0$. All images have boxsizes of $4\times 4$ parsec.

Open with DEXTER
In the text

  \begin{figure}
\par\begin{tabular}{c c}
\thicklines
\includegraphics[scale=0.38...
...g7.ps} &
\par\includegraphics[scale=0.38]{12236fg8.ps}\end{tabular}
\end{figure} Figure 4:

Temperature plots for four different metallicities at $t=12t_{\rm dyn}$, $Z=Z_{\odot }$ ( top left), $10^{-1}~Z_{\odot}$ ( top right), $10^{-2}~Z_{\odot}$ ( bottom left), $10^{-3}~Z_{\odot}$ ( bottom right), for $\beta _0 = 10^{-1}$. The color red represents temperatures above 50 K. Contour plots of number density are overplotted in these figures in black. Contour levels are 10 6, 105, 104, 103, 102, 101 cm-3. The axes are in units of 10 $^{18} \rm ~cm$ in both directions. In the bottom left corner of each image, the Jeans length is illustrated by a white circle. To calculate the Jeans length, the average value of the temperatures and densities of the clumped regions of each image is used. These are: T=(10.0, 13.1 13.4, 19.0)[K] and $\rho =(7.0$, 3.3, 0.41, 0.22) $[\rm \times 10^{-19}~g~cm^{-3}]$.

Open with DEXTER
In the text

  \begin{figure}
\begin{tabular}{c c}
\par\includegraphics{12236fg9-NEW.eps} &
\in...
...12236f15-NEW.eps} &
\includegraphics{12236f16-NEW.eps}\end{tabular}
\end{figure} Figure 5:

The phase diagrams of each metallicity after twelve dynamical timescales for two $\beta _0$ ratios, 10-1 ( left) and 10-3 ( right). Top to bottom, highest ( $Z=Z_{\odot }$) to lowest ( $Z=10^{-3}~Z_{\odot}$) metallicity.

Open with DEXTER
In the text

  \begin{figure}
\par\includegraphics[width=8.5cm,clip]{12236f17-NEW.eps}\hspace*{...
...}\hspace*{4mm}
\includegraphics[width=8.5cm,clip]{12236f20-NEW.eps}
\end{figure} Figure 6:

The phase diagrams of several nonrotating simulations with additional physics for different environments after twelve dynamical timescales. The two figures on the right represent isothermal simulations at 25 K ( top) and 10 K ( bottom). The top left one is a phase diagram of the solar metallicity simulation with dust and CR heating. The bottom left figure is for a nonrotating 3D simulation also with $Z=Z_{\odot }$.

Open with DEXTER
In the text

  \begin{figure}
\par\includegraphics[width=7.8cm,clip]{12236f21.eps}\hspace*{4mm}
\includegraphics[width=8cm,clip]{12236f22.eps}
\end{figure} Figure 7:

Average fragment mass ( $\langle M_{\rm f} \rangle$) with initial conditions as determined after twelve dynamical times. Left: $\langle M_{\rm f} \rangle$ versus metallicity for each initial rotational energy. Right: $\langle M_{\rm f} \rangle$ versus rotational energy for each metallicity.

Open with DEXTER
In the text

  \begin{figure}
\par\includegraphics[width=8.5cm,clip]{12236f23-NEW.eps}
\end{figure} Figure 8:

A diagram that shows the Jeans mass as a function of number density for a nonrotating solar metallicity run at $t=12t_{\rm dyn}$.

Open with DEXTER
In the text

  \begin{figure}
\par\includegraphics[width=16.8cm,clip]{12236f24.ps}
\end{figure} Figure 9:

Density plots of the environmental simulations at $t=12t_{\rm dyn}$. Molecular clouds in environments having dust and cosmic rays are simulated with solar metallicity. Orange to red depicts densities typically higher than $\rm 10^{4}~cm^{-3}$ and up to $\rm 10^{9}~cm^{-3}$. Light blue portrays those densities lower than $\rm 1 ~cm^{-3}$. All images have boxsizes of $4\times 4$ parsec.

Open with DEXTER
In the text

  \begin{figure}
\par\includegraphics[width=8.5cm,clip]{12236f25.eps}\end{figure} Figure 10:

Average fragment mass versus initial rotational energy for simulations within three distinct physical environments at $t=12t_{\rm dyn}$.

Open with DEXTER
In the text


Copyright ESO 2010

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.