EDP Sciences
Free Access
Issue
A&A
Volume 522, November 2010
Article Number A24
Number of page(s) 5
Section Extragalactic astronomy
DOI https://doi.org/10.1051/0004-6361/201015055
Published online 28 October 2010

© ESO, 2010

1. Introduction

At a distance of 100 pc, a massive 107 M black hole is able to emit X-rays with a flux of 100 erg s-1 cm-2 (Meijerink et al. 2007). The presence of X-rays in star forming regions alters fundamental processes in star forming molecular clouds, that eventually determine stellar masses, like ionization and gas and dust heating, and do so up to column densities of 1024 cm-2 and distances of 300 pc (Schleicher et al. 2010).

In our local neighborhood, the initial mass function (IMF) is observed to be a power-law, nicely following a Salpeter slope. The IMF is assumed to have a universal shape and spatially well resolved studies (Chabrier 2003; Sabbi et al. 2007; Elmegreen et al. 2008) confirm this at low and high masses. When we turn to observations of more radical environments like the Galaxy center, hints for different IMFs appear (Figer 2005a,b; Paumard et al. 2006; Espinoza et al. 2009; Elmegreen 2009; Bartko et al. 2010). However, the evidence is also claimed to be controversial (Bastian et al. 2010). Nuclei of active galaxies, e.g., ULIRGs like Arp 220 and Markarian 231, see van der Werf et al. (2010), enjoy conditions much more extreme than the Milky Way galactic center. Little is known on the IMF there, so theory and simulations are needed to guide understanding. The idea arose that massive stars form more easily in these extreme regions and that the IMF might be top-heavy (Brassington et al. 2007; Klessen et al. 2007; Dabringhausen et al. 2009). Other explanations are possible, still, it is worthwile to test this hypothesis numerically.

Until now, there have been many studies of feedback processes and external influences on star formation (Bonnell & Rice 2008; Wada 2008; Wada et al. 2009; Krumholz et al. 2010; Bate 2010). However, no numerical simulations have been performed yet of a collapsing molecular cloud under strong incident X-ray radiation. In this paper, we present the first 3D numerical simulations on the effects of X-rays on star forming molecular clouds that reside close to a massive black hole. In the next sections, we compare such a molecular cloud to an X-ray free environment and evaluate the star formation efficiency, stellar masses, and the resulting IMF.

2. Numerical model

2.1. The FLASH code

All simulations in this study have been performed with the Eulerian hydrodynamical and N-body code FLASH (Fryxell et al. 2000; Dubey et al. 2009). We used FLASH with the improved gravity solver that was added in version 3 (Ricker 2008). FLASH is a modular based, strongly-scaled parallel code that is specialized in adaptively refined meshes in which one can use particles in conjunction with the grid.

The hydrodynamic equations are solved using the piecewise parabolic method (PPM, Colella & Woodward 1984), which is an improved version of Godunov’s method (Godunov 1959). PPM is particularly well suited for flows involving discontinuities, such as shocks, that are strongly present in this study. FLASH is provided with many and extensively tested modules that encompass a broad range of physics. Many of these modules, including hydrodynamics, thermodynamics, (self-)gravity, particles, turbulence, and shocks are standard ingredients for interstellar physics and star formation and are thus incorporated in this work. Several additions were made to the code to follow the (radiation)physics in detail, such as sink particles, radiative transfer, multi-scale turbulence, and refinement criteria based on Jeans length and sink particles. The main additions are described in the following subsections.

2.2. Sink particles

Sink particles are point particles that can grow in mass by accreting gas and merge with other particles, however, they cannot lose mass or fragment. These particles represent compact objects, in our case protostars. Sink particles are a necessary ingredient if one wants to follow a density evolution over a large dynamic range, like a collapsing cloud. We created a sink particle module that takes care of the particle creation, determines their masses, handles Bondi-Hoyle type accretion, and tracks the particles if they are eligible for merging.

A sink particle is created when it passes several criteria for indefinite collapse as given by Federrath et al. (2010). The mass it obtains at creation is determined by a density threshold that is set by the maximum resolution. The Truelove et al. (1997) criterion states that the Jeans length should not be resolved by less then 4 cells if one wants to avoid artificial fragmentation. This criterion can be rewritten in the form of a density limit. Any excess density is taken away from the grid and added to the mass of the particle. The sink particle algorithm created for this study follows the method created by Krumholz et al. (2004). For any details of this method, we refer the reader to that paper.

2.3. Radiative transfer

In order to update the effects caused by X-rays on the temperature of the gas, we ported an X-ray dominated region chemical code (XDR code) created by Meijerink & Spaans (2005) into FLASH. This code incorporates all of the heating (photo-ionization, yielding non-thermal electrons) and cooling processes from atomic (fine-structure, semi-forbidden) and molecular transitions (CO, H2, H2O). Effects from internal UV, cosmic rays, and dust-gas coupling are treated as well. Given an X-ray flux, gas density, and column density along the line of sight to the source, the XDR code calculates the temperature and the chemical abundances. This output is fed into the simulation at every iteration. Most of the computation is spent finding the column densities for every cell. A ray-tracing algorithm, specifically created for this purpose, searches the grid and sums up the column densities of each cell lying along the line of sight from the source, the accreting black hole, to the target cell. The X-ray flux is an E-0.9 power law between 1 and 100 keV. X-ray scattering is not very important, but is nonetheless treated in the XDR-code. A uniform background of cosmic rays prevents the temperature from dropping below 10 K. For this, a cosmic ray ionization rate typical for the Milky Way ζ = 5 × 10-17 s-1, is assumed (Spaans & Silk 2005)

2.4. Initial conditions

We create a spherical gas cloud with solar abundances at a distance of 10 parsec from a 107 M black hole. We run two separate simulations and name them simulation A and B. Simulation A represents a molecular cloud near an active black hole under the impact of X-rays. Simulation B has a cloud near an inactive black hole and has isothermal conditions, with an equation of state of the form P ∝ ρ. The XDR code that updates the thermodynamics is only linked with simulation A where we do have X-rays. The 107 M black hole, at 10% of Eddington, yields a flux of 160 erg s-1 cm-2 with some extinction (Meijerink et al. 2007). Wada et al. (2009) show that column densities of 1022.5 cm-2 can typically exist in the central R ≃ 10 pc of an AGN, leading to a 1 keV optical depth of 2–3. Furthermore, Wada et al. (2009) show that column densities as large as 1024 cm-2 occur and persist in a statistical sense in the dynamically active inner 20 pc. A gaseous cloud can thus be shielded by this clumpy medium around the AGN and remain cold, while the temperature rises rapidly once the cloud is exposed to the radiation. Both our simulations start with the same initial conditions, but we expose simulation A to an X-ray source once the simulation starts. In this, heating is nearly immediate, since the timescale for heating is much shorter than the collapse time, theat ≪ tff, and of the order of 10-1 years. All other conditions are the same for both simulations.

The simulations are set up with an initial random, divergence-free turbulent velocity field and a characteristic FWHM of 5 km s-1 that agrees well with molecular clouds observed in active regions (Pérez-Beaupuits et al. 2009). These are supersonic flows with Mach numbers of up to 25, where the isothermal sound speed of the cloud is cs = 0.19 km s-1 (for T = 10 K) and can go up to a maximum of 5 km s-1 when the cloud is heated by X-rays (T ≃ 104 K). We do not drive the turbulence but follow its decay. The turbulence is applied over all scales with a power spectrum of P(k) ∝ k-4, following the empirical laws for compressible fluids (Larson 1981; Myers & Gammie 1999; Heyer & Brunt 2004). We start the simulations with a cloud that is in a stable Keplerian orbit around the black hole. Shear that is introduced by the black hole is taken into account. The maximum velocity difference imposed by the black hole,  △ vshear = 2.2 km s-1 across the cloud, is of the order of the applied initial turbulence. This process keeps the turbulence strong to large dynamical times, i.e.,  > 1 tff, with and 105 years1 throughout this work. The shearing time, tshear = Dcloud/ △ vshear, is almost 3 times larger than the cloud free-fall time and gravitationally bound (roughly) spherical clouds are likely to exist at densities of  ~105  cm-3.

The molecular cloud starts with a uniform number density of 105 cm-3 and has a size of 0.33 pc in radius. With a mean molecular weight of μ = 2.3, the total mass of the cloud amounts to 800 solar masses. The rest of the medium is filled with gas that has a uniform density of 100 cm-3. The simulation box, a cube of size 24 pc, has outflow boundaries and is isolated in terms of gravity. We increase the resolution where needed according to a self-developed Jeans criterion. The algorithm calculates the Jeans length at every grid cell, compares it against the Truelove criterion, and adds resolution when this is about to be violated. The maximum grid resolution that we allow for any simulation is 81923 cells. With the box size of 24 pc, the maximum spatial resolution becomes 8.8 × 1015 cm. If gas continues to collapse and needs higher resolution beyond the maximum refinement level, sink particles are created thereby taking density away from the gas such that mass and momentum is conserved and the resolution criterion is not violated. We give sink particles an accretion radius of effectively 3 cells, that is, 2.6 × 1016 cm  ≃  1760 AU. The temperature in both simulations is initialized at 10 K, but for the simulation with X-rays, the temperature is updated by the XDR code, and changes quickly (≲ 10-4 tff = 1 timestep) after initialization.

3. Results

3.1. X-rays versus no X-rays: IMF and Jeans mass

Both simulations, A and B, are followed for 2 × 105 years. This is approximately two free-fall times for a 105 cm-3 cloud. X-rays can heat the cloud to as high as 6000 K at low column densities (< 1021 cm-2), this is the case for the side that is directly exposed to the X-ray source, but also to as low as 10 K at high column densities (> 1024 cm-2). Figure 1 shows a column density plot to this effect.

The directional heating increases the pressure and causes the gas to expand and evaporate on the irradiated face of the cloud. The molecular cloud is compressed, loses mass, and an ionizing pressure flow travels inward. We see that this compression creates a density increase of about half an order of magnitude within a free-fall time as compared to simulation B. The conical compression front is disrupted where the turbulence creates sub-pc scale gaps (0.01–0.05 pc) and radiation is able to penetrate, as is illustrated in Fig. 2. Those regions are also heated up and pressurized. This causes the pressure front of the irradiated side of the cloud to break up. We see finger-like shapes forming, with a high density head, and the gas that is lying in its shadow is well shielded and very cold (about 10 K). The increased density induces star formation. We find that sink particles are created in the compressed cloud edge much earlier than in the shielded and colder parts. A phase diagram is plotted in Fig. 3 that shows the decline of temperature as the density increases. The secondary band with a steep decrease in temperature at low densities is the direct result of shielding. We note that the reason that fewer points appear in the plot in the shielded regions is merely due to the difference in resolution determined by the refinement criteria.

thumbnail Fig. 1

Column density (cm-2) slice of simulation A at t = tff. The colors represent the column density along the line of sight to the black hole, which is located at the upper left side (the black arrows indicate the direction of radiation). The density is shown in black contours. Contour levels are 1, 4, 16, 64, 256, 1024  × 104 cm-3. The white spheres display the location of the sink particles and the red arrow indicates the direction of the cloud’s orbital motion.

Open with DEXTER
thumbnail Fig. 2

Density morphology of simulation A at t = 0.4tff in 2D slices through the center of the cloud with axes in parsec. The color represents the number density (cm-3). On the left, an XY slice is shown with arrows representing the direction of radiation emanating from the black hole, which is located at the upper left side. On the right, an XZ slice is shown with radiation that is coming from the left and under an angle of 16 degrees with the slice.

Open with DEXTER

Simulation A starts forming protostars at around two-thirds of the initial free-fall time, which is later than simulation B, but at about three-quarters of the free-fall time, we see a sudden increase of star formation. Simulation B, on the other hand, has a gradual increase of protostars, starting at about a half tff, only to peak close to the free-fall time of 105 years. In the end, simulation B has created more sink particles, 153 against 118, but their masses are much lower. We plot the particle masses against the number of particles in Fig. 4.

The plots show that the X-ray dominated environment clearly has a higher characteristic mass as well as a higher minimum and maximum mass compared to simulation B. Higher sink particle masses are seen consistently throughout simulation A. Another difference is evident when looking at the slope of the IMF, which is defined as

(1)

with N the number of stars in a range of mass dM, α the power-law index, and Γ the slope above the characteristic mass of  ~0.3–0.5 M.

Simulation A has a much flatter slope than simulation B, which is ΓA =  −0.28 after 105 years of evolution, flatter than the Salpeter slope of –1.35 (Salpeter 1955). There is also a steep decrease in the number of stars above  ~3 M. This slope is very steep but uncertain and can depend strongly on late accretion and merging (Dib et al. 2010). Simulation B, has a slope of ΓB =  −1.37 which, in fact, is very similar to the Salpeter slope. Our simulations show that these trends for simulation B are maintained even at higher dynamical times, t = 2tff, in agreement with the time independent, competitive accretion driven results of Clark et al. (2008) for initial energies of |Egrav|  ⩾ Ekin. We attribute the difference in Γ partly to the higher gas temperatures (~50 K) caused by X-rays, which in turn increases the Jeans mass and the fragmentation mass scale. Most of the massive sink particles are created in the dense fragments and trace the warm X-ray heated (T > 50 K) molecular gas. The Jeans mass in case B is about 0.2 M after one free-fall time, whereas the Jeans mass in case A can be much higher, but typically ranges from 0.5 to 2 M in the  ≳ 106 cm-3 gas. The Jeans mass is approximately given by (Frieswijk et al. 2007);

(2)

In response to the high temperature in the XDR gas, only high density regions collapse creating more massive protostars and, in the end, fewer of them. Furthermore, the cloud is shielded from radiation with increasing column density, causing less X-ray penetration and thus less heating, lowering the Jeans mass. As a result, gas becomes more compressible and the density increases, starting a snowball effect. This causes a steep decrease of the Jeans mass, giving rise to a sharp increase in star formation. Another feature of the aforementioned effect is that almost all of the particles that are created in simulation A, are formed in the high density front, the finger head, of the cloud, where the snowball effect is initiated.

3.2. Sink particle behavior

We let the code check whether sink particles come too close to one another and if they should merge. This is the case if the velocity difference between two particles is less than the escape velocity and if the merging time is shorter than the simulation time step. Merging is not often seen. We define the merging time as half of a degenerate orbit time, , where m1 and m2 are the masses of two independent sink particles and R12 is the distance between the two.

Mergers can potentially increase the mass considerably, but most of the protostellar mass gain comes from accretion. To this end, we adopted a Bondi-Hoyle type of accretion (Krumholz et al. 2004). This type of accretion arises when a homogeneous flow of matter at infinity moves non-radially towards the accretor. Accretion increases with protostellar mass, but drops with decreasing ambient density and higher Mach numbers. We see that the ratio of the most massive to the second most massive sink particle in a local cloud region grows in time. With this, our expectation is confirmed that protostars compete with each other for material and that high mass stars as well as stars that lie in deep potential wells accrete more (Bonnell & Bate 2006; Bonnell 2008).

The total mass of the particles after one free-fall time for simulation B is 104 M, i.e., about 1/8th of the initial cloud mass. Quite interestingly, simulation A is able to convert more gas mass into stars, 230 M. This is a strong argument for efficient X-ray induced star formation in AGN. The star formation efficiencies mentioned here are upper limits, since these calculations do not contain feedback effects from young stars. Feedback effects, like outflows, should decrease the efficiency overall (Wang et al. 2010).

thumbnail Fig. 3

Temperature-density phase diagram of the X-ray irradiated simulation A after 1 dynamical evolution, t = tff.

Open with DEXTER
thumbnail Fig. 4

The initial mass function of the simulations at t = tff. Left panel shows the IMF for simulation A, and the right panel shows the IMF for simulation B. The green dashed line represents the Salpeter function. The blue dot-dashed line represents the Chabrier IMF, as fitted to simulation B. A power-law fit is applied to the data above the turn-over mass,  ~1.2 M left and  ~0.7 M right, and is shown as the solid purple line.

Open with DEXTER

4. Conclusions and discussion

We have performed two 3D simulations of similar molecular clouds, each at 10 pc distance from a supermassive black hole and followed their evolution. In one case A, we expose the cloud to an active black hole, producing a strong, 160 erg s-1 cm-2, X-ray flux. In the other case B the black hole was inactive and the molecular cloud had isothermal (10 K) conditions throughout the run. We saw clear differences between the simulations.

For the X-rays included run, we find that the molecular cloud is heated at the irradiated side and an ionizing pressure front is formed. This conic pressure front breaks up due to turbulent motions, forming finger-like structures that can be seen in column density plots, see Fig. 1. The density increases at the top of the compression front and the shielded gas cools. Star formation is initially delayed due to the higher temperatures, but the cloud continues to contract and after the critical density is reached sink particles start to from. Since the temperature in the shielded parts of the cloud decreases with increasing density, the Jeans mass drops at an amplified rate. Consequently, protostar formation increases sharply around 0.75 tff. Despite this, fewer sink particles are created in total with respect to simulation B, but they are more massive. The latter is a consequence of the high temperature (~50 K) and Jeans mass (≳ 1 M) in X-ray irradiated gas. In the end, case A ends up with more total stellar mass after one free-fall time. These protostars also accrete more material as the accretion rate scales with mass and density. This becomes increasingly more important for the mass growth due to the deeper potential well created and favored by massive stars. The two effects together cause that competitive accretion dominates the mass growth and strongly affects the shape of the mass function. The resulting IMF has a higher characteristic mass and a near-flat, non-Salpeter slope. We summarize the main points by stating how case A behaves with respect to B:

  • Protostars are created at a later stage, around0.65  t ff (0.45  t ff for B).

  • Fewer protostars are created but they have higher masses, M = 0.8–6 M (0.2–3.6 M for B).

  • Although protostars are created later, there is a burst mode around 3/4th of the collapse time of 10 5  years.

  • The total stellar mass formed from the gas is more than twice higher after 105 years and the efficiency is about 28% (13% for B).

  • Competitive accretion is more influential in shaping the IMF.

  • The IMF has a much flatter slope, Γ =   −0.28 (–1.37 for B).

  • The characteristic mass of the IMF is higher, Mchar ~ 1.2 M (~ 0.7 M for B).

All of these results tell us two main things. Firstly, star formation is induced by X-rays and massive stars can form in this way with high efficiency. Secondly, the resulting IMF differs from a Salpeter shape.

The pressure force exerted by radiation from the black hole is found to be modest. The bulk of the radiation pressure will results from the UV, since it dominates the bolometric luminosity L for very massive black holes. Typically, the change in momentum of the gas is proportional to τFIR L/c, for the speed of light c and mean FIR optical depth τFIR. The latter comes in because the dust, if it absorbs the bulk of the hard radiation, re-emits this in the FIR (100–300 μm). The ratio of radiation pressure over thermal pressure, for a typical τFIR of  ~0.1, is about 10% at 10 pc from the black hole.

We did not include any model for stellar feedback in this study. The total time of the simulations is not long enough for stellar evolution effects to play a major role, however, outflows from young stars can produce significant winds that disrupt gas build-up. This may affect the later time (> 1 tff) accretion onto protostars.

There is also strong gravitational shear, which helps to drive the turbulence but also compresses some regions in the center of the cloud along the direction perpendicular to the orbit. This increases the star formation efficiency overall (Bonnell & Rice 2008), but does not affect the conclusions when comparing the two cases with each other.

There are several other parameters that can be studied for their effectiveness, like the X-ray flux, turbulent strength, cloud size and initial gas density, distance to, and mass of, the black hole. In a follow-up paper we plan to present a detailed parameter study to this effect including the impact of UV and cosmic rays.


1

Given the initial condition of ρ = 3.84 × 10-19 g/cm-3.

Acknowledgments

We are very grateful to the anonymous referee for an insightful and constructive report 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. Part of the simulations have been run on the dedicated special purpose machines “Gemini” at the Kapteyn Astronomical Institute, University of Groningen. We are also grateful for the subsidy and time granted, by NCF, on the “Huygens” national supercomputer of SARA Amsterdam to run our high resolution simulations.

References

All Figures

thumbnail Fig. 1

Column density (cm-2) slice of simulation A at t = tff. The colors represent the column density along the line of sight to the black hole, which is located at the upper left side (the black arrows indicate the direction of radiation). The density is shown in black contours. Contour levels are 1, 4, 16, 64, 256, 1024  × 104 cm-3. The white spheres display the location of the sink particles and the red arrow indicates the direction of the cloud’s orbital motion.

Open with DEXTER
In the text
thumbnail Fig. 2

Density morphology of simulation A at t = 0.4tff in 2D slices through the center of the cloud with axes in parsec. The color represents the number density (cm-3). On the left, an XY slice is shown with arrows representing the direction of radiation emanating from the black hole, which is located at the upper left side. On the right, an XZ slice is shown with radiation that is coming from the left and under an angle of 16 degrees with the slice.

Open with DEXTER
In the text
thumbnail Fig. 3

Temperature-density phase diagram of the X-ray irradiated simulation A after 1 dynamical evolution, t = tff.

Open with DEXTER
In the text
thumbnail Fig. 4

The initial mass function of the simulations at t = tff. Left panel shows the IMF for simulation A, and the right panel shows the IMF for simulation B. The green dashed line represents the Salpeter function. The blue dot-dashed line represents the Chabrier IMF, as fitted to simulation B. A power-law fit is applied to the data above the turn-over mass,  ~1.2 M left and  ~0.7 M right, and is shown as the solid purple line.

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.