format: text

text:
# 1 Introduction

## 1 Introduction

Stars are important. They are the primary source of radiation (with competition from the 3K black body radiation of the cosmic microwave background and from accretion processes onto black holes in active galactic nuclei, which themselves are likely to have formed from stars), and of all chemical elements heavier than the H and He that made up the primordial gas. The Earth itself consists primarily of these heavier elements, called metals in astronomical terminology. Metals are produced by nuclear fusion in the interior of stars, with the heaviest elements produced during the passage of the final supernova shockwave through the most massive stars. To reach the chemical abundances observed today in our solar system, the material had to go through many cycles of stellar birth and death. In a literal sense, we are star dust.

Stars are also our primary source of astronomical information and, hence, are essential for our understanding of the universe and the physical processes that govern its evolution. At optical wavelengths almost all natural light we observe in the sky originates from stars. In daytime this is more than obvious, but it is also true at night. The Moon, the second brightest object in the sky, reflects light from our Sun, as do the planets, while virtually every other extraterrestrial source of visible light is a star or a collection of stars. Throughout the millenia, these objects have been the observational targets of traditional astronomy, and define the celestial landscape, the constellations. When we look at a dark night sky, we can also note dark patches of obscuration along the band of the Milky Way. These are clouds of dust and gas that block the light from stars further away.

Since about half a century ago we know that these clouds are associated with the birth of stars (for a historic account see Herbig 2002). The advent of new observational instruments and techniques gave access to astronomical information at wavelengths far shorter and longer that visible light. It is now possible to observe astronomical objects at wavelengths ranging high-energy $`\gamma `$-rays down to radio frequencies. Especially useful for studying these dark clouds are radio and sub-mm wavelengths, at which they are transparent. Observations now show that all star formation occurring in the Milky Way is associated with these dark clouds.

These clouds are dense enough, and well enough protected from dissociating UV radiation by self-shielding and dust scattering in their surface layers for hydrogen to be mostly in molecular form in their interior. The density and velocity structure of molecular clouds is extremely complex and follows hierarchical scaling relations that appear to be determined by supersonic turbulent motions (e.g. Williams, Blitz, & McKee 2000). Molecular clouds are large, and their masses exceed the threshold for gravitational collapse by far when taking only thermal pressure into account. Naively speaking, they should be contracting rapidly and form stars at very high rate. This is generally not observed. The star formation efficiency of molecular clouds in the solar neighborhood is estimated to be of order of a few percent (e.g. Elmegreen 1991, McKee 1999).

For many years it was thought that support by magnetic pressure against gravitational collapse offered the best explanation for the low rate of star formation. In this so called “standard theory of star formation”, developed by Shu (1977; and see Shu, Adams, & Lizano 1987), Mouschovias & Spitzer (1976), Nakano (1976), and others, interstellar magnetic fields prevent the collapse of gas clumps with insufficient mass to flux ratio, leaving dense cores in magnetohydrostatic equilibrium. The magnetic field couples only to electrically charged ions in the gas, though, so neutral atoms can only be supported by the field if they collide frequently with ions. The diffuse interstellar medium (ISM) with number densities $`n`$ of order unity remains ionized highly enough so that neutral-ion collisional coupling is very efficient (see Mouschovias 1991). In dense cores, where $`n>10^5`$ cm<sup>-3</sup>, ionization fractions drop below parts per million. Neutral-ion collisions no longer couple the neutrals tightly to the magnetic field, so the neutrals can diffuse through the field in a process known in astrophysics as ambipolar diffusion. This allows gravitational collapse to proceed in the face of magnetostatic support, but on a timescale as much as an order of magnitude longer than the free-fall time, drawing out the star formation process.

Recently, however, both observational and theoretical results have begun to cast doubt on the “standard theory” (for a recent compilation see Mac Low & Klessen 2003). While theoretical considerations point against singular isothermal spheres as starting conditions of protostellar collapse as postulated by the theory (see Whitworth et al. 1996, Nakano 1998, Desch & Mouschovias 2001), there is a series of observational findings that put other fundamental assumptions of the “standard theory” into question as well. For example, the observed magnetic field strengths in molecular cloud cores appear too weak to provide support against gravitational collapse (Crutcher 1999, Bourke et al. 2001). At the same time, the infall motions measured around star forming cores extend too broadly (e.g. Tafalla et al. 1998 or Williams et al. 1999 for L1544), while the central density profiles of cores are flatter than expected for isothermal spheres (e.g. Bacmann et al. 2000). Furthermore, the chemically derived ages of cloud cores are comparable to the free-fall time instead of the much longer ambipolar diffusion timescale (Bergin & Langer 1997). Observations of young stellar objects also appear discordant. Accretion rates appear to decrease rather than remain constant, far more embedded objects have been detected in cloud cores than predicted, and the spread of stellar ages in young clusters does not approach the ambipolar diffusion time (as discussed in the review by André et al. 2000).

These inconsistencies suggest to look beyond the standard theory, and we do so by seeking inspiration from the classical dynamical picture of star formation which we reconsider in the light of the recent progress in describing and understanding molecular cloud turbulence. Rather than relying on quasistatic evolution of magnetostatically supported objects, a new dynamical theory of star formation invokes supersonic interstellar turbulence to control the star formation process. We argue that this is both sufficient to explain star formation in Galactic molecular clouds and more consistent with observations.

Our line of reasoning leads us first to a general introduction of the concept of turbulence (Section 2), which is then followed by an analysis of its decay properties (Section 3). As our arguments rely to a large degree on results from numerical models we give a brief introduction into numerical simulations of supersonic turbulence (Section 4). We then discuss how local collapse can occur in globally stable interstellar gas clouds (Section 5) and investigate the physical processes that may prevent or promote this collapse (Section 6) leading to either more clustered or more isolated modes of star formation (Section 7). We deal with the timescales of star formation (Section 8) and discuss how the inclusion of magnetic fields may influence molecular cloud fragmentation (Section 9). We also discuss specific predictions of the new theory of turbulent star formation for protostellar mass accretion rates (Section 10) and for the resulting stellar mass spectra (Section 11). We then speculate about physical scales of interstellar turbulence in our Galaxy (Section 12), and ask what sets the overall efficiency of star formation (Section 13) and what terminates the process on scales of individual star forming regions (Section 14). At the end of this review (Section 15), we summarize our results and conclude that indeed the hypothesis that stellar birth is controlled by the complex interplay between supersonic turbulence and self-gravity offers an attractive pathway towards a consistent and comprehensive theory of star formation.

## 2 Turbulence

At this point, we should briefly discuss the concept of turbulence, and the differences between supersonic, compressible (and magnetized) turbulence, and the more commonly studied incompressible turbulence. We mean by turbulence, in the end, nothing more than the gas flow resulting from random motions at many scales. We furthermore will use in the discussion below only the very general properties and scaling relations of turbulent flows, focusing mainly on effects of compressibility. For a more detailed discussion of the complex statistical characteristics of turbulence, we refer the reader to the book by Lesieur (1997).

Most studies of turbulence treat incompressible turbulence, characteristic of most terrestrial applications. Root-mean-square (rms) velocities are subsonic, and density remains almost constant. Dissipation of energy occurs entirely in the centers of small vortices, where the dynamical scale $`\mathrm{}`$ is shorter than the length on which viscosity acts $`\mathrm{}_{\mathrm{visc}}`$. Kolmogorov (1941) described a heuristic theory based on dimensional analysis that captures the basic behavior of incompressible turbulence surprisingly well, although subsequent work has refined the details substantially. He assumed turbulence driven on a large scale $`L`$, forming eddies at that scale. These eddies interact to from slightly smaller eddies, transferring some of their energy to the smaller scale. The smaller eddies in turn form even smaller ones, until energy has cascaded all the way down to the dissipation scale $`\mathrm{}_{\mathrm{visc}}`$.

In order to maintain a steady state, equal amounts of energy must be transferred from each scale in the cascade to the next, and eventually dissipated, at a rate

$$\dot{E}=\eta v^3/L,$$
(1)
where $`\eta `$ is a constant determined empirically. This leads to a power-law distribution of kinetic energy $`Ev^2k^{10/3}`$, where $`k=2\pi /\mathrm{}`$ is the wavenumber, and density does not enter because of the assumption of incompressibility. Most of the energy remains near the driving scale, while energy drops off steeply below $`\mathrm{}_{\mathrm{visc}}`$. Because of the local nature of the cascade in wavenumber space, the viscosity only determines the behavior of the energy distribution at the bottom of the cascade below $`\mathrm{}_{\mathrm{visc}}`$, while the driving only determines the behavior near the top of the cascade at and above $`L`$. The region in between is known as the inertial range, in which energy transfers from one scale to the next without influence from driving or viscosity. The behavior of the flow in the inertial range can be studied regardless of the actual scale at which $`L`$ and $`\mathrm{}_{\mathrm{visc}}`$ lie, so long as they are well separated. The behavior of higher order structure functions $`S_p(\stackrel{}{r})=\{v(\stackrel{}{x})v(\stackrel{}{x}+\stackrel{}{r})\}^p`$ in incompressible turbulence has been successfully modeled by She & Leveque (1994) by assuming that dissipation occurs in the filamentary centers of vortex tubes.

Gas flows in the ISM vary from this idealized picture in a number of important ways. Most significantly, they are highly compressible, with Mach numbers $``$ ranging from order unity in the warm, diffuse ISM, up to as high as 50 in cold, dense molecular clouds. Furthermore, the equation of state of the gas is very soft due to radiative cooling, so that pressure $`P\rho ^\gamma `$ with the polytropic index falling in the range $`0.4<\gamma <1.2`$ (e.g. Scalo et al. 1998, Ballesteros-Paredes, Vázquez-Semadeni, & Scalo 1999b, Spaans & Silk 2000). Supersonic flows in highly compressible gas create strong density perturbations. Early attempts to understand turbulence in the ISM (von Weizsäcker 1943, 1951, Chandrasekhar 1949) were based on insights drawn from incompressible turbulence. Although the importance of compressibility was already understood, how to incorporate it into the theory remained unclear. Furthermore, compressible turbulence is only one physical process that may cause the strong density inhomogeneities observed in the ISM. Others are thermal phase transitions (Field, Goldsmith, & Habing 1969, McKee & Ostriker 1977, Wolfire et al. 1995) or gravitational collapse (e.g. Wada & Norman 1999).

In supersonic turbulence, shock waves offer additional possibilities for dissipation. Shock waves can transfer energy between widely separated scales, removing the local nature of the turbulent cascade typical of incompressible turbulence. The spectrum may shift only slightly, however, as the Fourier transform of a step function representative of a perfect shock wave is $`k^2`$, so the associated energy spectrum should be close to $`\rho v^2k^4`$, as was indeed found by Porter, Pouquet, & Woodward (1994). However, even in hypersonic turbulence, the shock waves do not dissipate all the energy, as rotational motions continue to contain a substantial fraction of the kinetic energy, which is then dissipated in small vortices. However, Boldyrev (2002) has proposed a theory of structure function scaling based on the work of She & Leveque (1994) using the assumption that dissipation in supersonic turbulence primarily occurs in sheet-like shocks, rather than linear filaments. This model has been extended to describe density structure functions and exhibits in general good agreement with numerical results (see Boldyrev, Nordlund, & Padoan 2002).

The driving of interstellar turbulence is neither uniform nor homogeneous. Controversy still reigns over the most important energy sources at different scales, but it appears likely that isolated and correlated supernovae dominate (Mac Low & Klessen 2003). However, it is not yet understood at what scales expanding, interacting blast waves contribute to turbulence. Analytic estimates have been made based on the radii of the blast waves at late times (Norman & Ferrara 1996), but never confirmed with numerical models (much less experiment).

Finally, interstellar gas is magnetized. Although magnetic field strengths are difficult to measure, with Zeeman line splitting being the best quantitative method, it appears that fields within an order of magnitude of equipartition with thermal pressure and turbulent motions are pervasive in the diffuse ISM, most likely maintained by a dynamo driven by the motions of the interstellar gas. A model for the distribution of energy and the scaling behavior of strongly magnetized, incompressible turbulence based on the interaction of shear Alfvén waves is given by Goldreich & Sridhar (1995, 1997) and Ng & Bhattacharjee (1996). The scaling properties of the structure functions of such turbulence was derived from the work of She & Leveque (1994) by Müller & Biskamp (2000; also see Biskamp & Müller 2000) by assuming that dissipation occurs in current sheets. A theory of very weakly compressible turbulence has been derived by using the Mach number $`1`$ as a perturbation parameter (Lithwick & Goldreich 2001), but no further progress has been made towards analytic models of strongly compressible magnetohydrodynamic (MHD) turbulence with $`1`$.

With the above in mind, we propose that stellar birth is regulated by interstellar turbulence and its interplay with gravity. Turbulence, even if strong enough to counterbalance gravity on global scales, will usually provoke local collapse on small scales. Supersonic turbulence establishes a complex network of interacting shocks, where converging flows generate regions of high density. This density enhancement can be sufficient for gravitational instability. Collapse sets in. However, the random flow that creates local density enhancements also may disperse them again. For local collapse to actually result in the formation of stars, collapse must be sufficiently fast for the region to ‘decouple’ from the flow, i.e. it must be shorter than the typical time interval between two successive shock passages. The shorter this interval, the less likely a contracting region is to survive. Hence, the efficiency of star formation depends strongly on the properties of the underlying turbulent velocity field, on its lengthscale and strength relative to gravitational attraction. This principle holds for star formation throughout all scales considered, ranging from small local star forming regions in the solar neighborhood up to galaxies as a whole (see Mac Low & Klessen 2003).

## 3 Decay and Maintenance of Supersonic Motions

We first consider the question of how to maintain the observed supersonic motions in molecular clouds. As described above, magnetohydrodynamic waves were generally thought to provide the means to prevent the dissipation of interstellar turbulence. However, numerical models have now shown that they probably do not. One-dimensional simulations of decaying, compressible, isothermal, magnetized turbulence by Gammie & Ostriker (1996) showed quick decay of kinetic energy $`K`$ in the absence of driving, but found that the quantitative decay rate depended strongly on initial and boundary conditions because of the low dimensionality. Mac Low et al. (1998), Stone, Ostriker & Gammie (1998), and Padoan & Nordlund (1999) measured the decay rate in direct numerical simulations in three dimensions, using a number of different numerical methods. They uniformly found rather faster decay, with Mac Low et al. (1998) characterizing it as $`Kt^\eta `$, with $`0.85<\eta <1.1`$. A resolution and algorithm study is shown in Figure 1. Magnetic fields with strengths ranging up to equipartition with the turbulent motions (ratio of thermal to magnetic pressures as low as $`\beta =0.025`$) do indeed reduce $`\eta `$ to the lower end of this range, but not below that, while unmagnetized supersonic turbulence shows values of $`\eta 11.1`$.

Stone et al. (1998) and Mac Low (1999) showed that supersonic turbulence decays in less than a free-fall time under molecular cloud conditions, regardless of whether it is magnetized or unmagnetized. The hydrodynamical result agrees with the high-resolution, transsonic, decaying models of Porter et al. (1994). Mac Low (1999) showed that the formal dissipation time $`\tau _\mathrm{d}=K/\dot{K}`$ scaled in units of the free fall time $`t_{\mathrm{ff}}`$ is

$$\tau _\mathrm{d}/\tau _{\mathrm{ff}}=\frac{1}{4\pi \xi }\left(\frac{32}{3}\right)^{1/2}\frac{\kappa }{_{\mathrm{rms}}}\mathrm{\hspace{0.17em}3.9}\frac{\kappa }{_{\mathrm{rms}}},$$
(2)
where $`\xi =0.21/\pi `$ is the energy-dissipation coefficient, $`_{\mathrm{rms}}=v_{\mathrm{rms}}/c_\mathrm{s}`$ is the rms Mach number of the turbulence, and $`\kappa `$ is the ratio of the driving wavelength to the Jeans wavelength $`\lambda _\mathrm{J}`$, which is the critical scale for gravitational collapse to set in (Jeans 1902). In molecular clouds, $`_{\mathrm{rms}}`$ is typically observed to be of order 10 or higher. If the ratio $`\kappa <1`$, as is probably required to maintain gravitational support (Léorat et al. 1990), then even strongly magnetized turbulence will decay long before the cloud collapses and not markedly retard the collapse.

Either observed supersonic motions must be continually driven, or molecular clouds must be less than a single free-fall time old. Observational evidence does suggest that clouds are a few free-fall times old, on average, though perhaps not more than two or three, so there is likely some continuing energy input into the clouds (Ballesteros-Paredes, Hartmann, & Vázquez-Semadeni 1999a, Fukui et al. 1999, Elmegreen 2000).

## 4 Modeling Turbulence in Self-Gravitating Gas

This leads to the question of what effects supersonic turbulence will have on self-gravitating clouds. Can turbulence alone delay gravitational collapse beyond a free-fall time? Most analytical approaches to that question are based on the assumption that the turbulent flow is close to incompressible, and are therefore not applicable to interstellar turbulence. However, some more recent models have made certain progress in recovering the velocity structure of compressible turbulence as well (Boldyrev 2002, Boldyrev, Nordlund, & Padoan 2002).

Numerical models of highly compressible, self-gravitating turbulence have shown the importance of density fluctuations generated by the turbulence to understanding support against gravity. Early models were done by Bonazzola et al. (1987), who used low resolution ($`32\times 32`$ collocation points) calculations with a two-dimensional spectral code to support their analytical results. The hydrodynamical studies by Passot et al. (1988), Léorat et al. (1990), Vázquez-Semadeni, Passot, & Pouquet (1995) and Ballesteros-Paredes et al. (1999b), were also restricted to two dimensions, and were focused on the interstellar medium at kiloparsec scales rather than molecular clouds, although they were performed with far higher resolution (up to $`800\times 800`$ points). Magnetic fields were introduced in these models by Passot, Vázquez-Semadeni, & Pouquet (1995), and extended to three dimensions with self-gravity (though at only $`64^3`$ resolution) by Vázquez-Semadeni, Passot, & Pouquet (1996). One-dimensional computations focused on molecular clouds, including both MHD and self-gravity, were presented by Gammie & Ostriker (1996) and Balsara, Crutcher & Pouquet (2001). Ostriker, Gammie, & Stone (1999) extended their work to 2.5 dimensions more recently.

These models at low resolution, low dimension, or both, suggested several important conclusions. First, gravitational collapse, even in the presence of magnetic fields, does not generate sufficient turbulence to markedly slow continuing collapse. Second, turbulent support against gravitational collapse may act at some scales, but not others.

More recently, three-dimensional high-resolution computations by Klessen (2000), Klessen, Heitsch, & Mac Low (2000) and Heitsch, Mac Low, & Klessen (2001a) have confirmed both of these results. These authors used two different numerical methods: ZEUS-3D (Stone & Norman 1992ab), an Eulerian MHD code; and an implementation of smoothed particle hydrodynamics (SPH; Benz 1990, Monaghan 1992), a Lagrangian hydrodynamics method using particles as an unstructured grid. Both codes were used to examine the gravitational stability of three-dimensional hydrodynamical turbulence at high resolution. The use of both Lagrangian and Eulerian methods to solve the equations of self-gravitating hydrodynamics in three dimensions (3D) allowed them to attempt to bracket reality by taking advantage of the strengths of each approach. This gives some protection against interpreting numerical artifacts as physical effects (for a detailed discussion see Klessen et al. 2000).

The computations discussed here were done on periodic cubes, with an isothermal equation of state, using up to $`256^3`$ zones (with one model at $`512^3`$ zones) or $`80^3`$ SPH particles. To generate turbulent flows Gaussian velocity fluctuations are introduced with power only in a narrow interval $`k1|\stackrel{}{k}|k`$, where $`k=L/\lambda _\mathrm{d}`$ counts the number of driving wavelengths $`\lambda _\mathrm{d}`$ in the box (Mac Low et al. 1998). This offers a simple approximation to driving by mechanisms acting on that scale. To drive the turbulence, this fixed pattern is normalized to maintain constant kinetic energy input rate $`\dot{E}_{\mathrm{in}}=\mathrm{\Delta }E/\mathrm{\Delta }t`$ (Mac Low 1999). Self-gravity is turned on only after a state of dynamical equilibrium has been reached.

## 5 Local versus Global Collapse

First we examine the question of whether gravitational collapse can generate enough turbulence to prevent further collapse. Hydrodynamical SPH models initialized at rest with Gaussian density perturbations show fast collapse, with the first collapsed objects forming in a single free-fall time (Klessen, Burkert, & Bate 1998; Klessen & Burkert 2000, 2001). Models set up with a freely decaying turbulent velocity field behaved similarly (Klessen 2000). Further accretion of gas onto collapsed objects then occurs over the next free-fall time, defining the predicted spread of stellar ages in a freely-collapsing system. The turbulence generated by the collapse (or virialization) does not prevent further collapse contrary to what sometimes has been suggested (e.g. by Elmegreen 1993). The presence of magnetic fields does not change that conclusion (Balsara et al. 2001) as accretion down filaments aligned with magnetic field lines onto cores can occur readily. This allows high mass-to-flux ratios to be maintained even at small scales, which is necessary for supercritical collapse to continue after fragmentation occurs.

Second, we examine whether continuously driven turbulence can provide support against gravitational collapse. The models of driven, self-gravitating turbulence by Klessen et al. (2000) and Heitsch et al. (2001a) show that local collapse occurs even when the turbulent velocity field carries enough energy to counterbalance gravitational contraction on global scales. An example of local collapse in a globally supported cloud is given in Figure 2. A hallmark of global turbulent support is isolated, inefficient, local collapse.

Thus, highly compressible turbulence does both, it promotes as well as prevents collapse. Its net effect is to inhibit collapse globally, while at the same time promoting it locally. The resolution to this apparent paradox lies in the requirement that any substantial turbulent support must come from supersonic flows, as otherwise pressure support would be at least equally important. Supersonic flows compress the gas in shocks. In isothermal gas with density $`\rho `$ the postshock gas has density $`\rho ^{}=^2\rho `$, where $``$ is the Mach number of the shock. The turbulent Jeans length $`\lambda _\mathrm{J}\rho ^{1/2}`$ in these density enhancements, so it drops by a factor of $``$ in isothermal shocks making shock compressed gas clumps more susceptible to gravitational collapse. On the other hand, if we consider the system on scales exceeding the lengthscale of turbulence (i.e. in the limit of microturbulence), we can follow the classical picture that treats turbulence as an additional pressure and define an effective sound speed $`c_{\mathrm{s},eff}^2=c_\mathrm{s}^2+v_{\mathrm{rms}}^2/3`$ (Chandrasekhar 1949). The critical mass for gravitational collapse, the Jeans mass $`M_\mathrm{J}\rho ^{1/2}c_\mathrm{s}^3`$, then strongly increases with the turbulent rms velocity dispersion $`v_{\mathrm{rms}}`$, so that for $`v_{\mathrm{rms}}c_\mathrm{s}`$ turbulence ultimately does inhibit collapse on global scales. Between these two scales, there is a broad intermediate region, especially for long wavelength driving, where local collapse can occur despite global support.

Klessen et al. (2000) demonstrated that turbulent support can completely prevent collapse only when it can support not just the average density, but also these high-density shocked regions, a point that was appreciated already by Elmegreen (1993) and Vázquez-Semadeni et al. (1995). Two criteria must be fulfilled: the rms velocity must be sufficiently high for the turbulent Jeans criterion to be met in these regions, and the driving wavelength $`\lambda _\mathrm{d}<\lambda _\mathrm{J}(\rho ^{})`$. If these two criteria are not fulfilled, the high-density regions collapse, although the surrounding flow remains turbulently supported. The efficiency of collapse depends on the properties of the supporting turbulence. Sufficiently strong driving on short enough scales can prevent local collapse for arbitrarily long periods of time, but such strong driving may be rather difficult to arrange in a real molecular cloud. Furthermore, if we assume that stellar driving sources have an effective wavelength close to their separation, then the condition that driving acts on scales smaller then the Jeans wavelength in ‘typical’ shock generated gas clumps requires the presence of an extraordinarily large number of stars evenly distributed throughout the cloud, with typical separation 0.1 pc in Taurus, or only 350 AU in Orion. This is not observed. Very small driving scales seem also to be at odds with the observed large-scale velocity fields at least in some molecular clouds (e.g. Ossenkopf & Mac Low 2002).

## 6 Promotion and Prevention of Local Collapse

The origin of local collapse can also be understood in terms of a timescale argument. Roughly speaking, the lifetime of a clump is determined by the interval between two successive passing shocks: the first creates it, while if the second is strong enough, it disrupts the clump again if it has not already collapsed (Klein, McKee & Colella 1994, Mac Low et al. 1994). If its lifetime is long enough, a Jeans unstable clump can contract to sufficiently high densities to effectively decouple from the ambient gas flow. It then becomes able to survive the encounter with further shock fronts (e.g. Krebs & Hillebrandt 1983), and continues to accrete from the surrounding gas, forming a dense core. The weaker the passing shocks, and the greater the separation between them, the more likely that collapse will occur. Equivalently, weak driving and long typical driving wavelengths enhance collapse. The influence of the driving wavelength is more pronounced, however, because individual shocks sweep up more mass when the typical wavelength is longer, so density enhancements resulting from the interaction of shocked layers will have larger masses, and so are more likely to exceed their local Jeans limit. Turbulent driving mechanisms that act on large scales will produce large coherent structures (filaments of compressed gas with embedded dense cores) on relatively short timescales compared to small-scale driving even if the total kinetic energy in the system is the same. Examples of the density structure of long and small-wavelength driving, respectively, are given in Figure 3, which can be directly compared to Figure 2b.

Further insight of how local collapse proceeds comes from examining the mass growth rates in each model. Figure 4 shows the accretion history for three sets of models from Klessen et al. (2000). The driving strength increases from $`𝒜`$ over $``$ to $`𝒞`$, but is held constant for each set of models with the effective driving wavelength $`\lambda _\mathrm{d}`$ being varied. All models show local collapse, except at the extreme end, when $`\lambda _\mathrm{d}<\lambda _\mathrm{J}(\rho ^{})`$ (model $`5`$).

The cessation of strong accretion onto cores occurs long before all gas has been accreted. This is because the time that dense cores spend in shock-compressed, high-density regions decreases with increasing driving wavenumber and increasing driving strength. In the case of long wavelength driving, cores form coherently in high-density regions associated with one or two large shock fronts that can accumulate a considerable fraction of the total mass of the system. The overall accretion rate is high and cores spend sufficient time in this environment to accrete a large fraction of the total mass in the region. Any further mass growth has to occur from chance encounters with other dense regions. In the case of short wavelength driving, the network of shocks is tightly knit. Cores form in shock generated clumps of small masses because individual shocks are not able to sweep up much matter. Furthermore, in this rapidly changing environment the time interval between the formation of clumps and their destruction is short. The period during which individual cores are located in high-density regions where they are able to accrete at high rate is short as well. So altogether, the global accretion rates are small and saturate at lower values of $`M_{}`$ as the driving wavelength is decreased.

## 7 Clustered versus Isolated Star Formation

Different star formation regions present different distributions of protostars and pre-main sequence stars. In some regions, such as the Taurus molecular cloud, stars form isolated from other stars, scattered throughout the cloud (Mizuno et al. 1995). In other regions, they form in clusters, as in L1630 in Orion (Lada 1992), or even more extremely in starburst regions such as 30 Doradus (Walborn et al. 1999; for a review see Zinnecker et al. 1993).

Numerical simulations of self-gravitating turbulent clouds demonstrate that the length scale and strength at which energy is inserted into the system determine the structure of the turbulent flow and therefore the locations at which stars are most likely to form. Large-scale driving leads to large coherent shock structures (see e.g. Figure 3a). Local collapse occurs predominantly in filaments and layers of shocked gas and is very efficient in converting gas into stars. This leads to what we can identify as ‘clustered’ mode of star formation: stars form in coherent aggregates and clusters. Even more so, this applies to regions of molecular gas that have become decoupled from energy input. As turbulence decays, these regions begin to contract and form dense clusters of stars with very high efficiency on about a free-fall time scale (Klessen et al. 1998, Klessen & Burkert 2000). The same holds for insufficient support, i.e. for regions where energy input is not strong enough to completely balance gravity. They too will contract to form dense stellar clusters.

The ‘isolated’ mode of star formation occurs in molecular cloud regions that are supported by driving sources that act on small scales and in an incoherent or stochastic manner. In this case, individual shock induced density fluctuations form at random locations and evolve more or less independently of each other. The resulting stellar population is widely dispersed throughout the cloud and, as collapsing clumps are exposed to frequent shock interaction, the overall star formation rate is low.

These points are illustrated in Figure 5, which shows the distribution of collapsed cores in several models with strong enough turbulence to formally support against collapse. Coherent, efficient local collapse occurs in model $`1`$, where the turbulence is driven strongly at long wavelengths (compare with Figure 3). Incoherent, inefficient collapse occurs in model $`3`$, on the other hand, where turbulence is driven at small scales. Individual cores form independently of each other at random locations and random times, are widely distributed throughout the entire volume, and exhibit considerable age spread. In the decaying turbulence model, once the kinetic energy level has decreased sufficiently, all spatial modes of the system contract gravitationally, including the global ones (Klessen 2000). As in the case of large-scale shock compression, stars form more or less coevally in a limited volume with high efficiency.

Despite the fact that both turbulence driven on large scales and freely decaying turbulence lead to star formation in aggregates and clusters, Figure 6 suggests a possible way to distinguish between them. Decaying turbulence typically leads to the formation of a bound stellar cluster, while aggregates associated with large-scale, coherent, shock fronts often have higher velocity dispersions that result in their complete dispersal. Note, however, that at the late stages of dynamical evolution shown in Figure 6, the model becomes less appropriate, as feedback from newly formed stars is not included. Ionization and outflows from the stars formed first will likely retard or even prevent the accretion of the remaining gas onto protostars, possibly preventing a bound cluster from forming even in the case of freely decaying turbulence.

The control of star formation by supersonic turbulence gives rise to a continuous but articulated picture. There may not be physically distinct modes of star formation, but qualitatively different behaviors do appear over the range of possible turbulent flows. The apparent dichotomy between a clustered mode of star formation and an isolated one, as discussed by Lada (1992) for L1630 and Strom, Strom, & Merrill (1993) for L1641, disappears, if a different balance between turbulent strength and gravity holds at the relevant length scales in these different clouds.

Turbulent flows tend to have hierarchical structure (Lesieur 1997) which may explain the hierarchical distribution of stars in star forming regions shown by statistical studies of the distribution of neighboring stars in young stellar clusters (e.g. in Taurus, see Larson 1995). Hierarchical clustering seems to be a common feature of all star forming regions (e.g. Efremov & Elmegreen 1998). It is a natural outcome of turbulent fragmentation.

## 8 Timescales of Star Formation

Turbulent control of star formation predicts that stellar clusters form predominantly in regions that are insufficiently supported by turbulence or where only large-scale driving is active. In the absence of driving, molecular cloud turbulence decays more quickly than the free-fall timescale $`\tau _{\mathrm{ff}}`$ (Eq. 2), so dense stellar clusters will form on the free-fall timescale. Even in the presence of support from large-scale driving, substantial collapse still occurs within a few free-fall timescales, see Figure and 7a. If the dense cores followed in these models continue to collapse on a short timescale to build up stellar objects in their centers, then this directly implies the star formation timescale. Therefore the age distribution will be roughly $`\tau _{\mathrm{ff}}`$ for stellar clusters that form coherently with high star formation efficiency. When scaled to low densities, say $`n(\mathrm{H}_2)10^2\mathrm{cm}^3`$ and $`T10`$K, the global free-fall timescale in the models is $`\tau _{\mathrm{ff}}=3.3\times 10^6`$years. If star forming clouds such as Taurus indeed have ages of order $`\tau _{\mathrm{ff}}`$, as suggested by Ballesteros-Paredes et al. (1999), then the long star formation time computed here is quite consistent with the very low star formation efficiencies seen in Taurus (e.g. Leisawitz et al. 1989, Palla & Stahler 2000, Hartmann 2001), as the cloud simply has not had time to form many stars. In the case of high-density regions, $`n(\mathrm{H}_2)10^5\mathrm{cm}^3`$ and $`T10`$K, the dynamical evolution proceeds much faster and the corresponding free-fall times drops to $`\tau _{\mathrm{ff}}=10^5`$years. These values are indeed supported by observational data such as the formation timescale of the Orion Trapezium cluster. It is inferred to stem from gas of density $`n(\mathrm{H}_2)\stackrel{<}{}\mathrm{\hspace{0.25em}10}^5\mathrm{cm}^3`$, and is estimated to be less than $`10^6`$ years old (Hillenbrand & Hartmann 1998). The age spread in the models increases with increasing driving wavenumber $`k`$ and increasing $`M_\mathrm{J}_{\mathrm{turb}}`$, as shown in Figure 7. Long periods of core formation for globally supported clouds appear consistent with the low efficiencies of star-formation in regions of isolated star formation, such as Taurus, even if they are rather young objects with ages of order $`\tau _{\mathrm{ff}}`$.

## 9 Effects of Magnetic Fields

So far, we concentrated on the effects of purely hydrodynamic turbulence. How does the picture discussed here change, if we consider the presence of magnetic fields? Magnetic fields may alter the dynamical state of a molecular cloud sufficiently to prevent gravitationally unstable regions from collapsing (McKee 1999). They have been hypothesized to support molecular clouds either magnetostatically or dynamically through MHD waves.

Mouschovias & Spitzer (1976) derived an expression for the critical mass-to-flux ratio in the center of a cloud for magnetostatic support. Assuming ideal MHD, a self-gravitating cloud of mass $`M`$ permeated by a uniform flux $`\mathrm{\Phi }`$ is stable if the mass-to-flux ratio

$$\frac{M}{\mathrm{\Phi }}<\left(\frac{M}{\mathrm{\Phi }}\right)_{\mathrm{cr}}$$
(3)
with $`(M/\mathrm{\Phi })_{\mathrm{cr}}=c_\mathrm{\Phi }G^{1/2}`$. The exact value depends on the geometry and the field and density distribution of the cloud. A cloud is termed subcritical if it is magnetostatically stable and supercritical if it is not. Mouschovias & Spitzer (1976) determined that $`c_\mathrm{\Phi }=0.13`$ for spherical clouds. Assuming a constant mass-to-flux ratio in a cylindrical region results in $`c_\mathrm{\Phi }=1/(2\pi )0.16`$ (Nakano & Nakamura 1978). Without any other mechanism of support such as turbulence acting along the field lines, a magnetostatically supported cloud will collapse to a sheet which then will be supported against further collapse. Fiege & Pudritz (2000) discussed a sophisticated version of this magnetostatic support mechanism, in which poloidal and toroidal fields aligned in the right configuration could prevent a cloud filament from fragmenting and collapsing.

Investigation of the second alternative, support by MHD waves, concentrates mostly on the effect of Alfvén waves, as they (1) are not as subject to damping as magnetosonic waves and (2) can exert a force along the mean field, as shown by Dewar (1970) and Shu et al. (1987). This is because Alfvén waves are transverse waves, so they cause perturbations $`\delta \stackrel{}{B}`$ perpendicular to the mean magnetic field $`\stackrel{}{B}`$. McKee & Zweibel (1995) argue that Alfvén waves can even lead to an isotropic pressure, assuming that the waves are neither damped nor driven. However, in order to support a region against self-gravity, the waves would have to propagate outwardly, rather than inwardly, which would only further compress the cloud. Thus, as Shu et al. (1987) comment, this mechanism requires a negative radial gradient in wave sources in the cloud.

It can be demonstrated (e.g. Heitsch et al. 2001a) that supersonic turbulence does not cause a magnetostatically supported region to collapse, and vice versa, that in the absence of magnetostatic support, MHD waves cannot completely prevent collapse, although they can retard it to some degree. In the case of a subcritical region with $`M<M_{cr}`$ sheets of high density gas form perpendicular to the field lines. Turbulence can shift the sheets along the field lines without changing the mass-to-flux ratio, but collapse does not occur, because the shock waves cannot sweep gas across field lines and the entire region is initially supported magnetostatically.

A supercritical cloud with $`M>M_{cr}`$ could only be stabilized by MHD wave pressure. This is insufficient to completely prevent gravitational collapse, as shown in Figure 8. The effect of the magnetic field on the morphology of the cloud is week, and collapse occurs in all models of unmagnetized and magnetized turbulence regardless of the numerical resolution and magnetic field strength as long as the system is magnetically supercritical. This is shown more quantitatively in Figure 9. Effects of numerical resolution make themselves felt in different ways in hydrodynamical and MHD models. In the hydrodynamical case, higher resolution results in thinner shocks and thus higher peak densities. These higher density peaks form cores with deeper potential wells that accrete more mass and are more stable against disruption. Higher resolution in the MHD models, on the other hand, better resolves short-wavelength MHD waves, which apparently can delay collapse, but not prevent it. This result extends to models with $`512^3`$ zones (Heitsch et al. 2001b).

## 10 Mass Growth of Protostellar Cores

Supersonic turbulence is able to produce star forming regions that vary enormously in scale. The most likely outcome of turbulent molecular cloud fragmentation in the Milky Way are stellar aggregates or clusters (Adams & Myers 2001). The number density of protostars and protostellar cores in the extreme cases the can be high enough for mutual dynamical interaction to become important. This has important consequences for the mass growth history of individual stars and the subsequent dynamical evolution of the nascent stellar cluster, because this introduces a further degree of stochasticity to the star formation process in addition to the statistical chaos associated with turbulence and turbulent fragmentation in the first place.

Klessen (2001a) considers the formation of a nascent star cluster for the case where turbulence is decayed and has left behind random Gaussian fluctuations in the density structure. As the system contracts gravitationally, a dense cluster of protostellar cores builds up on a timescale of about two to three free-fall times. The protostellar accretion rates in this environment are strongly time variable, as illustrated in Figure 10, which is a direct result of the mutual dynamical interaction and competition between protostellar cores. While gas clumps collapse to build up protostars, they may merge as they follow the flow pattern towards the cluster potential minimum. The timescales for both processes are comparable. The density and velocity structure of merged gas clumps generally differs significantly from their progenitor clumps, and the predictions for isolated cores are no longer valid. More importantly, these new larger clumps contain multiple protostars, which subsequently compete with each other for the accretion from a common gas reservoir. The most massive protostar in a clump is hereby able to accrete more matter than its competitors (also Bonnell et al. 1997, Klessen & Burkert 2000, Bonnell et al. 2001). Its accretion rate is enhanced through the clump merger, whereas the accretion rate of low-mass cores typically decreases. Temporary accretion peaks in the wake of clump mergers are visible in abundance in Figure 10. Furthermore, the small aggregates of cores that build up are dynamically unstable and low-mass cores may be ejected. As they leave the high-density environment, accretion terminates and their final mass is reached.

The typical density profiles of gas clumps that give birth to protostars exhibit a flat inner core, followed by a density fall-off $`\rho r^2`$, and are truncated at some finite radius (see Figure 13 in Klessen & Burkert 2000), which in the dense centers of clusters often is due to tidal interaction with neighboring cores. As result, a short-lived initial phase of strong accretion occurs when the flat inner part of the pre-stellar clump collapses. This corresponds to the class 0 phase of protostellar evolution (André et al. 2000). If these cores were to remain isolated and unperturbed, the mass growth rate would gradually decline in time as the outer envelope accretes onto the center. This is the class I phase. Once the truncation radius is reached, accretion fades and the object enters the class II phase. This behavior is expected from analytical models (e.g. Henriksen et al. 1997) and agrees with other numerical studies (e.g. Foster & Chevalier 1993). However, collapse does not start from rest for the density fluctuations considered here, and the accretion rates exceed the theoretically predicted values even for the most isolated objects in the simulation.

The most massive protostars begin to form first and continue to accrete at high rate throughout the entire cluster evolution. As the most massive gas clumps tend to have the largest density contrast, they are the first to collapse and constitute the center of the nascent cluster. These protostars are fed at high rate and gain mass very quickly. As their parental clumps merge with others, more gas is fed into their ‘sphere of influence’. They are able to maintain or even increase the accretion rate when competing with lower-mass objects (e.g. core 1 and 8 in Figure 10). Low-mass stars, on average, tend to form somewhat later in the dynamical evolution of the system (as indicated by the absolute formation times in Figure 10), and typically have only short periods of high accretion.

As high-mass stars are associated with large core masses, while low-mass stars come from low-mass cores, the stellar population in clusters is predicted to be mass segregated right from the beginning. High-mass stars form in the center, lower-mass stars tend to form towards the cluster outskirts. This is in agreement with recent observational findings for the cluster NGC330 in the Small Magellanic Cloud (Sirianni et al. 2002). Dynamical effects during the embedded phase of star cluster evolution will enhance this initial segregation even further.

## 11 Mass Spectra from Turbulent Fragmentation

As discussed before, a full understanding of turbulent molecular cloud fragmentation should in principle allow for a prediction of the distribution of stellar masses (e.g. Larson 1981, Fleck 1982, Elmegreen 1993, Padoan 1995, Padoan & Nordlund). However, a complete theory of compressible interstellar turbulence is still out of reach, and we have to resort to numerical modeling instead to make some progress. To illustrate this point we examine the mass spectra of gas clumps and collapsed cores from models of self-gravitating, isothermal, supersonic turbulence driven with different wavelengths (Klessen 2001b). In the absence of magnetic fields and more accurate equations of state, these models can only be illustrative, not definitive, but nevertheless they offer insight into the processes acting to form the initial stellar mass function (IMF; for a review see Kroupa 2002). Figure 11 plots for four different models the mass distribution of gas clumps, of the subset of gravitationally unstable clumps, and of collapsed cores, at four different evolutionary phases. In the initial phase, before local collapse begins to occur, the clump mass spectrum is not well described by a single power law. During subsequent evolution, as clumps merge and grow bigger, the mass spectrum extends towards larger masses, approaching a power law with slope $`\alpha 1.5`$. Local collapse sets in and results in the formation of dense cores most quickly in the freely collapsing model. The influence of gravity on the clump mass distribution weakens when turbulence dominates over gravitational contraction on the global scale, as in the other three models. The more the turbulent energy dominates over gravity, the more the spectrum resembles the initial case of pure hydrodynamic turbulence. This suggests that the clump mass spectrum in molecular clouds will be shallower in regions where gravity dominates over turbulent energy. This may explain the observed range of slopes for the clump mass spectrum in different molecular cloud regions (e.g. Kramer et al. 1998).

Like the distribution of Jeans-unstable clumps, the mass spectrum of dense protostellar cores resembles a log-normal in the model without turbulent support and in the one with long-wavelength turbulent driving, with a peak at roughly the average thermal Jeans mass $`m_\mathrm{J}`$ of the system. However, models supported at shorter wavelength have mass spectra much flatter than observed, suggesting that clump merging and competitive accretion are important factors leading to a log-normal mass spectrum. The protostellar clusters discussed here only contain between 50 and 100 cores. This allows for comparison with the IMF only around the characteristic mass scale, typically about 1 M, since the numbers are too small to study the very low- and high-mass end of the distribution. Focusing on low-mass star formation, however, Bate, Bonnell, & Bromm (2002) demonstrate that brown dwarfs are a natural and frequent outcome of turbulent fragmentation. In this model, brown dwarfs form when dense molecular gas fragments into unstable multiple systems that eject their smallest members from the dense gas before they have been able to accrete to stellar masses. Numerical models with sufficient dynamic range to treat the full range of stellar masses remain yet to be done.

## 12 Scales of Interstellar Turbulence

Turbulence has self-similar properties only within a certain range of scales. The upper end typically is associated with the global extent of the system or with the scale at which energy is inserted. The lower scale is the energy dissipation scale of the system where turbulent kinetic energy is converted into random motion, into heat. In purely hydrodynamic systems this is the scale where molecular viscosity becomes important.

In interstellar clouds this situation may be different. It was first shown by Zweibel & Josafatsson (1983) that ambipolar diffusion would be the most important dissipation mechanism in typical molecular clouds with very low ionization fractions $`x=\rho _i/\rho _n`$, where $`\rho _i`$ is the density of ions, $`\rho _n`$ is the density of neutrals, and $`\rho =\rho _i+\rho _n`$. The ambipolar diffusion strength is defined as

$$\lambda _{\mathrm{AD}}=v_\mathrm{A}^2/\nu _{ni},$$
(4)
where $`v_\mathrm{A}^2=B^2/4\pi \rho _n`$ approximates the effective Alfvén speed for the coupled neutrals and ions if $`\rho _n\rho _i`$, and $`\nu _{ni}=\gamma \rho _i`$ is the rate at which each neutral is hit by ions. The coupling constant depends on the cross-section for ion-neutral interaction, and for typical molecular cloud conditions has a value of $`\gamma 9.2\times 10^{13}`$ cm$`^3`$s<sup>-1</sup>g<sup>-1</sup> (e.g. Smith & Mac Low 1997). Zweibel & Brandenburg (1997) define an ambipolar diffusion Reynolds number as

$$R_{\mathrm{AD}}=\stackrel{~}{L}\stackrel{~}{V}/\lambda _{\mathrm{AD}}=_\mathrm{A}\stackrel{~}{L}\nu _{ni}/v_\mathrm{A},$$
(5)
which must fall below unity for ambipolar diffusion to be important, where $`\stackrel{~}{L}`$ and $`\stackrel{~}{V}`$ are the characteristic length and velocity scales, and $`_\mathrm{A}=\stackrel{~}{V}/v_\mathrm{A}`$ is the characteristic Alfvén Mach number. In our situation we again can take the rms velocity as typical value for $`\stackrel{~}{V}`$. By setting $`R_{\mathrm{AD}}=1`$, we can derive a critical lengthscale below which ambipolar diffusion is important

$$\stackrel{~}{L}_{cr}=\frac{v_\mathrm{A}}{_\mathrm{A}\nu _{ni}}(0.041\text{ pc})\left(\frac{B}{10\mu \mathrm{G}}\right)_{\mathrm{A}}^{}{}_{}{}^{1}\left(\frac{x}{10^6}\right)^1\left(\frac{n_n}{10^3\mathrm{cm}^3}\right)^{3/2},$$
(6)
with the magnetic field strength $`B`$, the ionization fraction $`x`$, the neutral number density $`n_n`$, and where we have taken $`\rho _n=\mu n_n`$, with $`\mu =2.36m_\mathrm{H}`$. This is consistent with typical sizes of protostellar cores (e.g. Bacmann et al. 2000), if we assume that ionization and magnetic field both depend on the density of the region and follow the empirical laws $`n_i=3\times 10^3\mathrm{cm}^3(n_n/10^5\mathrm{cm}^3)^{1/2}`$ (e.g. Mouschovias 1991) and $`B30\mu \mathrm{G}(n_n/10^3\mathrm{cm}^3)^{1/2}`$ (e.g. Crutcher 1999).

On large scales, an maximum upper limit to the turbulent cascade in the Milky Way is given by the extent and thickness of the Galactic disk. This is indeed the true upper scale, if molecular clouds are created by converging large-scale flows generated by the collective influence of recurring supernovae explosions. For individual molecular clouds this means, that turbulent energy is fed in at scales beyond the size of the considered cloud. The bulk of its turbulent energy content is then generated in the event of cloud assembly, which then dissipates rapidly resulting in short cloud life times. The same compressional motions that created the cloud in the first place, however, may also act as continuing source of kinetic energy during some initial period (e.g. Walder & Folini 2000, Hartmann, Ballesteros-Paredes, & Bergin 2001), thus extending the overall cloud lifetime to a few crossing times. This energy cascades down to supply turbulence on smaller scales within the cloud. This picture of molecular cloud turbulence being driven by large-scale, external sources is strongly support by analysis of velocity structure which is always is dominated by the large-scale modes in all clouds observed (Ossenkopf & Mac Low 2002).

## 13 Efficiency of Star Formation

The global star formation efficiency in normal molecular clouds is usually estimated to be of the order of a few per cent. Their life times may be on the order of a few crossing times, i.e. a few $`10^6`$ years (Ballesteros-Paredes et al. 1999a, Fukui et al. 1999, Elmegreen 2000). In this case nearly all models of interstellar turbulence discussed above are consistent with the observed overall efficiencies. If molecular clouds survive for several tens of their free-fall time $`\tau _{\mathrm{ff}}`$ (i.e. a few $`10^7`$ years as proposed by Blitz & Shu 1980), turbulence models are more strongly constrained. However, even in this case models with parameters reasonable for Galactic molecular clouds can maintain global efficiencies below $`M_{}=5`$% for $`10\tau _{\mathrm{ff}}`$ (Klessen et al. 2000). Furthermore, it needs to be noted that the local star formation efficiency in molecular clouds can reach very high values. For example, the Trapezium star cluster in Orion is likely to have formed with an efficiency of about 50% (Hillenbrand & Hartmann 1998), compared to a value of 5% proposed for Taurus-Aurigae.

## 14 Termination of Local Star Formation

It remains quite unclear what terminates stellar birth on scales of individual star forming regions, and even whether these processes are the primary factor determining the overall efficiency of star formation in a molecular cloud. Three main possibilities exist. First, feedback from the stars themselves in the form of ionizing radiation and stellar outflows may heat and stir surrounding gas up sufficiently to prevent further collapse and accretion. Second, accretion might peter out either when all the high density, gravitationally unstable gas in the region has been accreted in individual stars, or after a more dynamical period of competitive accretion, leaving any remaining gas to be dispersed by the background turbulent flow. Third, background flows may sweep through, destroying the cloud, perhaps in the same way that it was created in the first place. Most likely the astrophysical truth lies in some combination of all three possibilities.

If a stellar cluster formed in a molecular cloud contains OB stars, then the radiation field and stellar wind from these high-mass stars strongly influence the surrounding cloud material. The UV flux ionizes gas out beyond the local star forming region. Ionization heats the gas, raising its Jeans mass, and possibly preventing further protostellar mass growth or new star formation. The termination of accretion by stellar feedback has been suggested at least since the calculations of ionization by Oort & Spitzer (1955). Whitworth (1979) and Yorke et al. (1989) computed the destructive effects of individual blister Hii regions on molecular clouds, while Franco et al. (1994) and Diaz-Miller et al. (1998) concluded that indeed the ionization from massive stars may limit the overall star forming capacity of molecular clouds to about 5%. Matzner (2002) analytically modeled the effects of ionization on molecular clouds, concluding as well that turbulence driven by Hii regions could support and eventually destroy molecular clouds. The key question facing these models is whether Hii region expansion couples efficiently to clumpy, inhomogeneous molecular clouds, a question probably best addressed with numerical simulations.

Bipolar outflows are a different manifestation of protostellar feedback, and may also strongly modify the properties of star forming regions (Norman & Silk 1980, Adams & Fatuzzo 1996). Recently Matzner & McKee (2000) modeled the ability of bipolar outflows to terminate low-mass star formation, finding that they can limit star formation efficiencies to 30–50%, although they are ineffective in more massive regions. How important these processes are compared to simple exhaustion of available reservoirs of dense gas (Klessen et al. 2000) remains an important question.

The models relying on exhaustion of the reservoir of dense gas argue that only dense gas will actually collapse, and that only a small fraction of the total available gas reaches sufficiently high densities, due to cooling (Schaye 2002), gravitational collapse and turbulent triggering (Elmegreen 2002), or both (Wada, Meurer, & Norman 2002). This of course pushes the question of local star formation efficiency up to larger scales, which may indeed be the correct place to ask it.

Other models focus on competitive accretion in local star formation, showing that the distribution of masses in a single group or cluster can be well explained by assuming that star formation is fairly efficient in the dense core, but that stars that randomly start out slightly heavier tend to fall towards the center of the core and accrete disproportionately more gas (Bonnell et al. 1997, 2001). These models have recently been called into question by the observation that the stars in lower density young groups in Serpens simply have not had the time to engage in competitive accretion, but still have a normal IMF (Olmi & Testi 2002).

Finally, star formation in dense clouds created by turbulent flows may be terminated by the same flows that created them. Ballesteros-Paredes et al. (1999a) suggested that the coordination of star formation over large molecular clouds, and the lack of post-T Tauri stars with ages greater than about 10 Myr tightly associated with those clouds, could be explained by their formation in a larger-scale turbulent flow. Hartmann et al. (2001) make the detailed argument that these flows may disrupt the clouds after a relatively short time, limiting their star formation efficiency that way. It can be argued that field supernovae are the most likely driver for this background turbulence in spiral galaxies like the Milky Way (Mac Low & Klessen 2003).

## 15 Summary: The Control of Star Formation by Supersonic Turbulence

In this review we have proposed that star formation is regulated by interstellar turbulence and its interplay with gravity. We have discussed that this new approach can explain the same observations successfully described by the so called “standard theory”, while also addressing (and resolving!) its inconsistencies with other observed properties of Galactic star forming regions.

The key point to this new understanding of star formation in Galactic molecular clouds lies in the properties of interstellar turbulence. Turbulence is observed in the interstellar medium almost ubiquitously and is typically supersonic as well as super-Alfvénic. It is energetic enough to counterbalance gravity on global scales, but at the same time it may provoke local collapse on small scales. This apparent paradox can be resolved when considering that supersonic turbulence establishes a complex network of interacting shocks, where converging flows generate regions of high density. This density enhancement can be sufficiently large for gravitational instability to set in. The same random flow that creates density enhancements, however, may disperse them again. For local collapse to result in stellar birth, it must progress sufficiently fast for the region to ‘decouple’ from the flow. Typical collapse timescales are hereby of the same order as the lifetimes of shock-generated density fluctuations in the turbulent gas. This makes the outcome highly unpredictable. As stars are born through a sequence of stochastic events, any theory of star formation is in essence a statistical one with quantitative predictions only possible for an ensemble of stars.

In the new picture, the efficiency of protostellar core formation, the growth rates and final masses of the protostars, essentially all properties of nascent star clusters depend on the intricate interplay between gravity on the one hand side and the turbulent velocity field in the cloud on the other. The star formation rate is regulated not just at the scale of individual star-forming cores through ambipolar diffusion balancing magnetostatic support, but rather at all scales (Elmegreen 2002), via the dynamical processes that determine whether regions of gas become unstable to prompt gravitational collapse. The presence of magnetic fields does not alter that picture significantly, as long as they are too weak for magnetostatic support, which is indicated by observations (Crutcher 1999, Bourke et al. 2001). In particular, magnetic fields cannot prevent the decay of interstellar turbulence, which in turn needs to be continuously driven or else stars form quickly and with high efficiency

Inefficient, isolated star formation will occur in regions which are supported by turbulence carrying most of its energy on very small scales. This typically requires an unrealistically large number of driving sources and appears at odds with the measured velocity structure in molecular clouds which in almost all cases is dominated by large-scale modes. The dominant pathway to star formation therefore seems to involve cloud regions large enough to give birth to aggregates or clusters of stars. This is backed up by careful stellar population analysis indicating that most stars in the Milky Way formed in open clusters with a few hundred member stars (Adams & Myers 2001).

Clusters of stars build up in molecular cloud regions where self-gravity overwhelms turbulence, either because such regions are compressed by a large-scale shock, or because interstellar turbulence is not replenished and decays on short timescales. Then, many gas clumps become gravitationally unstable synchronously and start to collapse. If the number density is high, collapsing gas clumps may merge to produce new clumps which now contain multiple protostars. Mutual dynamical interactions become common, with close encounters drastically altering the protostellar trajectories, thus changing the mass accretion rates. This has important consequences for the IMF. Already in their infancy, i.e. already in the deeply embedded phase, very dense stellar clusters are expected to be strongly influenced by collisional dynamics.

### Acknowledgments

This review would not have been possible without long-term collaboration and exchange of ideas with M.-M. Mac Low. Special thanks also to P. Bodenheimer and D. Lin for many vivid scientific discussions and for their warm hospitality at UC Santa Cruz; and thanks to J. Ballesteros-Paredes, F. Heitsch, P. Kroupa, E. Vázquez-Semadeni, and H. Zinnecker.

I want to express my gratitudes to the members of the Astronomische Gesellschaft for awarding the Ludwig Biermann Preis to me; in particular, I want to thank G. Hensler and A. Burkert in this context. I furthermore acknowledge support by the Emmy Noether Program of the Deutsche Forschungsgemeinschaft (DFG: KL1358/1).

### References

Adams, F. C., and M. Fatuzzo, 1996, ApJ, 464, 256

Adams, F. C., and P. C. Myers, 2001, ApJ, 553, 744

André, P., D. Ward-Thompson, and M. Barsony, 2000, in Protostars and Planets IV, edited by V. Mannings, A. P. Boss, and S. S. Russell (University of Arizona Press, Tucson), p. 59

Bacmann, A., P. André, J. -.L Puget, A. Abergel, S. Bontemps, and D. Ward-Thompson, 2000, A&A, 361, 555

Ballesteros-Paredes, J., L. Hartmann, and E. Vázquez-Semadeni, 1999a, ApJ, 527, 285

Ballesteros-Paredes, J., E. Vázquez-Semadeni, and J. Scalo, 1999b, ApJ, 515, 286

Balsara, D. S., R. M. Crutcher, and A. Pouquet, 2001, ApJ, 557, 451

Bate, M. R., I. A. Bonnell, and V. Bromm, 2002, MNRAS, 332, L65

Benz, W., 1990, in The Numerical Modelling of Nonlinear Stellar Pulsations, edited by J. R. Buchler (Kluwer, Dordrecht), 269

Bergin, E. A., and W. D. Langer, 1997, ApJ, 486, 316

Biskamp, D., and W.-C. Müller, 2000, Phys. Plasmas, 7, 4889

Blitz, L., and F. H. Shu, 1980, ApJ, 238, 148

Boldyrev, S., 2002, ApJ, 569, 841

Boldyrev, S., Å. Nordlund, and P. Padoan, 2002, ApJ, 573, 678

Bonazzola, S., E. Falgarone, J. Heyvaerts, M. Perault, and J. L. Puget, 1987, A&A, 172, 293

Bonnell, I. A., M. R. Bate, C. J. Clarke, and J. E. Pringle, 1997, MNRAS, 285, 201

Bonnell, I. A., M. R. Bate, C. J. Clarke, and J. E. Pringle, 2001, MNRAS, 323, 785

Bourke, T. L., P. C. Myers, G. Robinson, and A. R. Hyland, 2001, ApJ, 554, 916

Chandrasekhar, S., 1949, ApJ, 110, 329

Crutcher, R. M., 1999, ApJ, 520, 706

Desch, S. J. and T. C. Mouschovias, 2001, ApJ, 550, 314

Dewar, R. L., 1970, Phys. Fluids, 13, 2710

Diaz-Miller, R. I., J. Franco, and S. N. Shore, 1998, ApJ, 501, 192.

Efremov, Y. N., and B. G. Elmegreen, 1998, MNRAS, 299, 588

Elmegreen, B. G., 1991, in NATO ASIC Proc. 342: The Physics of Star Formation and Early Stellar Evolution, edited by C. J. Lada and N. D. Kylafis (Kluwer Academic Publishers), p. 35

Elmegreen, B. G., 1993, ApJ, 419, L29

Elmegreen, B. G., 2000, MNRAS, 311, L5

Elmegreen, B. G., 2002, ApJ, 577, 206

Fiege, J. D., and R. E. Pudritz, 2000, MNRAS, 311, 85

Field, G. B., D. W. Goldsmith, and H. J. Habing, 1969, ApJ, 155, L49

Fleck, R. C., 1982, MNRAS, 201, 551

Foster, P. N., and R. A. Chevalier, 1993, ApJ, 416, 303

Franco, J., S. N. Shore, and G. Tenorio-Tagle, 1994, ApJ, 436, 795

Fukui, Y. et al., 1999, PASJ, 51, 745

Gammie, C. F., and E. C. Ostriker, 1996, ApJ, 466, 814

Goldreich, P. and S. Sridhar, 1995, ApJ, 438, 763

Goldreich, P. and S. Sridhar, 1997, ApJ, 485, 680

Hartmann, L., 2001, AJ, 121, 1030

Hartmann, L., J. Ballesteros-Paredes, and E. A. Bergin, 2001, ApJ, 562, 852

Heitsch, F., M. Mac Low, and R. S. Klessen, 2001a, ApJ, 547, 280

Heitsch, F., E. G. Zweibel, M.-M. Mac Low, P. Li, and M. L. Norman, 2001b, ApJ, 561, 800

Herbig, G, 2002, in Saas Fee Advanced Course 29: Physics of Star Formation in Galaxies, edited by F. Palla, H. Zinnecker, A. Maeder, and G. Meynet (Springer, Berlin, Heidelberg), p. 1

Hillenbrand, L. A., and L. W. Hartmann, 1998, ApJ, 492, 540

Hendriksen, R. N., P. André, and S. Bontemps, 1997, A&A, 323, 549

Jeans, J. H., 1902, Phil. Trans. A., 199, 1

Klein, R. I., C. F. McKee, and P. Colella, 1994, ApJ, 420, 213

Klessen, R. S., 2000, ApJ, 535, 869

Klessen, R. S., and A. Burkert, 2000, ApJS, 128, 287

Klessen, R. S., and A. Burkert, 2001, ApJ, 549, 386

Klessen, R. S., A. Burkert, and M. R. Bate, 1998, ApJ, 501, L205

Klessen, R. S., F. Heitsch, and M.-M. Mac Low, 2000, ApJ, 535, 887

Kolmogorov, A. N., 1941, Dokl. Akad. Nauk SSSR, 30, 301 (reprinted in Proc. R. Soc. Lond. A, 434, 9-13 )

Kramer, C., J. Stutzki, R. Rohrig, U. Corneliussen, 1998, A&A, 329, 249

Krebs, J., and W. Hillebrandt, 1983, A&A, 128, 411

Kroupa, P., 2002, Science 295, 82

Lada, E. A., 1992, ApJ, 393, L25

Larson, R. B., 1981, MNRAS, 194, 809

Larson, R. B., 1995, MNRAS, 272, 213

Leisawitz, D., F. N. Bash, and P. Thaddeus, 1989, ApJS, 70, 731

Léorat, J., T. Passot, and A. Pouquet, 1990, MNRAS, 243, 293

Lesieur, M., 1997, Turbulence in Fluids, 3rd ed. (Kluwer, Dordrecht), p. 245

Lithwick, Y., and P. Goldreich, 2001, ApJ, 562, 279

Mac Low, M.-M., 1999, ApJ, 524, 169

Mac Low, M.-M., Klessen, R. S., 2003, Rev. Mod. Phys., in press (astro-ph/0301093)

Mac Low, M.-M., R. S. Klessen, A. Burkert, and M. D. Smith, 1998, Phys. Rev. Lett., 80, 2754

Mac Low, M.-M.,, C. F. McKee, R. I. Klein, J. M. Stone, and M. L. Norman, 1994, ApJ, 433, 757

Matzner, C. D., 2002, ApJ, 566, 302

Matzner, C. D., and C. F. McKee, 2000, ApJ, 545, 364

McKee, C. F., 1999, in NATO ASIC Proc. 540: The Origin of Stars and Planetary Systems, edited by C. J. Lada and N. D. Kylafis (Kluwer Academic Publishers), p. 29

McKee, C. F., and J. P. Ostriker, 1977, ApJ, 218, 148

McKee, C. F., and E. G. Zweibel, 1995, ApJ, 440, 686

Mizuno, A., T. Onishi, Y. Yonekura, T. Nagahama, H. Ogawa, and Y. Fukui, 1995, ApJ, 445, L161

Monaghan, J. J., 1992, ARAA, 30, 543

Mouschovias, T. C., 1991, in The Physics of Star Formation and Early Stellar Evolution, edited by C. J. Lada and N. D. Kylafis (Kluwer, Dordrecht), p. 449

Mouschovias, T. C., and L. Spitzer, Jr., 1976, ApJ, 210, 326

Müller, W.-C., and D. Biskamp, 2000, Phys. Rev. Lett., 84, 475

Nakano, T., 1976, PASJ, 28, 355

Nakano, T., 1998, ApJ, 494, 587

Nakano, T., and T. Nakamura, 1978, PASJ, 30, 681

Ng, C. S., and A. Bhattacharjee, 1996, ApJ, 465, 845

Norman, C. A., and A. Ferrara, 1996, ApJ, 467, 280

Norman, C. A., and J. Silk, 1980, ApJ, 239, 968

Olmi, L., and L. Testi, 2002, A&A, 392, 1053

Oort, J. H., & L. Spitzer, Jr., 1955, ApJ, 121, 6

Ostriker, E. C., C. F. Gammie, and J. M. Stone, 1999, ApJ, 513, 259

Ossenkopf V., and M.-M. Mac Low, 2002, A&A, 390, 307

Padoan, P., 1995, MNRAS, 277, 377

Padoan, P., and Å. Nordlund, 1999, ApJ, 526, 279

Padoan, P., and Å. Nordlund, 2002, ApJ, 576, 870

Palla, F., and S. W. Stahler, 2000, ApJ, 540, 255

Passot, T., A. Pouquet, and P. R. Woodward, 1988, A&A, 197, 392

Passot, T., E. Vázquez-Semadeni, and A. Pouquet, 1995, ApJ, 455, 536

Porter, D. H., A. Pouquet, and P. R. Woodward, 1994, Phys. Fluids, 6, 2133

Salpeter, E. E., 1955, ApJ, 121, 161

Scalo, J. M., E. Vázquez-Semadeni, D. Chappell, T. Passot, 1998, ApJ, 504, 835

Schaye, J., 2002, ApJ, submitted (astro-ph/0205125)

She, Z., and E. Leveque, 1994, Phys. Rev. Lett., 72, 336

Shu, F. H., 1977, ApJ, 214, 488

Shu, F. H., F. C. Adams, and S. Lizano, 1987, ARAA, 25, 23

Sirianni, M., A. Nota, G. De Marchi, C. Leitherer, and M. Clampin, 2002, ApJ, 579, 275

Smith, M. D., and M.-M. Mac Low, 1997, A&A, 326, 801

Spaans, M., and J., Silk, 2000, ApJ, 538, 115

Stone, J. M., and M. L. Norman, 1992a, ApJS, 80, 753

Stone, J. M., and M. L. Norman, 1992b, ApJS, 80, 791

Stone, J. M., E. C. Ostriker, and C. F. Gammie, 1998, ApJ, 508, L99

Strom, K. M., S. E. Strom, and K. M. Merrill, 1993, ApJ, 412, 233

Tafalla, M., D. Mardones, P. C. Myers, P. Caselli, R. Bachiller, and P. J. Benson, 1998, ApJ, 504, 900

Vázquez-Semadeni, E., T. Passot, and A. Pouquet, 1995, ApJ, 441, 702

Vázquez-Semadeni, E., T. Passot, and A. Pouquet, 1996 ApJ, 473, 881

von Weizsäcker, C. F., 1943, Z. Astrophys., 22, 319

von Weizsäcker, C. F., 1951, ApJ, 114, 165

Wada, K., and C. A. Norman, 1999, ApJ, 516, L13

Wada, K., G. Meurer, and C. A. Norman, 2002, ApJ, 577, 197

Walborn, N. R., R. H. Barbá, W. Brandner, M. ;. Rubio, E. K. Grebel, and R. G. Probst, 1999, AJ, 117, 225

Walder, R. and D. Folini, 2000, ApSS, 274, 343

Whitworth, A. P., 1979, MNRAS, 186, 59

Whitworth, A. P., A. S. Bhattal, N. Francis, and S. J. Watkins, 1996, MNRAS, 283, 1061

Williams, J. P., L. Blitz, and C. F. McKee, 2000, in Protostars and Planets IV, edited by V. Mannings, A. P. Boss, and S. S. Russell (University of Arizona Press, Tucson), p. 97

Williams, J. P., P. C. Myers, D. J. Wilner, and J. di Francesco, 1999, ApJ, 513, L61

Wolfire, M. G., D. Hollenbach, C. F. McKee, A. G. G. M. Tielens, and E. L. O. Bakes, 1995, ApJ, 443, 152

Yorke, H. W., Tenorio-Tagle, G., Bodenheimer, P., and M. Różyczka, 1989, A&A, 216, 207

Zinnecker, H., M. J. McCaughrean, and B. A. Wilking. 1993, in Protostars and Planets III, edited by E. H. Levy and J. I. Lunine (University of Arizona Press, Tucson), p. 429

Zweibel, E. G., and A. Brandenburg, 1997, ApJ, 478, 563

Zweibel, E. G., and K. Josafatsson, 1983, ApJ, 270, 511


text:
# Non-equilibrium distributions at finite noise intensities

## Abstract

We analyse the non-equilibrium distribution in dissipative dynamical systems at finite noise intensities. The effect of finite noise is described in terms of topological changes in the pattern of optimal paths. Theoretical predictions are in good agreement with the results of numerical solution of the Fokker-Planck equation and Monte Carlo simulations.

The analysis of the probability distribution in systems away from thermal equilibrium is one of the main themes in statistical physics Onsager:53 ; Haken:75  which has recently attracted renewed interest in the context of stochastic resonance Nicolis:82  and stochastic ratchets Dykman:97a . Despite its importance, the theory of the probability distributions in non-equilibrium systems is still in its infancy.

One of the keys to the solution of this problem is an understanding of how a non-equilibrium system deviates far away from a steady state Onsager:53 . The probability of such deviations is small and can be approximated by a Boltzmann-like distribution of the form Ludwig:75a ; Freidlin:84a

$$\rho (x)=z(x)\mathrm{exp}\left(S(x)/ϵ\right),ϵ0.$$
(1)
Here the action function $`S(x)`$ plays the role of a non-equilibrium potential Graham:84a , $`z(x)`$ is a prefactor and $`ϵ`$ is the noise intensity. The majority of earlier studies have been restricted to the analysis of zero noise limit. In this limit properties of $`\rho (x)`$ are determined by the properties of optimal paths that are deterministic trajectories along which the system moves in the course of large deviations Freidlin:84a ; Ludwig:75a . Accordingly, the pattern of optimal paths plays a fundamental role in the zero-noise-limit theory.

In “real life”, however, the noise intensity is always finite. Calculations of $`\rho (x)`$ at finite noise intensity, that allow for a direct comparison with an experiment, require analysis of the all-important prefactor. Such analysis has recently been initiated in a number of publications Smelyanskiy:99a ; Hanggi:00a ; Maier:01a  dealing with the closely related problem of noise-induced escape over an unstable limit cycle.

The research has shown that one has “to leave room for the possibility that, even for small noise strengths $`ϵ`$, more than one” optimal path contributes significantly to the escape Hanggi:00a ; Maier:01a . However, this plausible assumption contradicts the fundamental properties of the pattern of optimal paths including the existence of the singularities and generic uniqueness of the most probable escape path (see e.g. Graham:84a ; Jauslin:87a ; Dykman:94a ). This inconsistency does not allow the use of the methods of Smelyanskiy:99a ; Hanggi:00a ; Maier:01a  to compute noise-induced variations of the non-equilibrium distribution in the whole of phase space. As mentioned in Maier:01a  the problem is “intricate” due to the subtle properties of the pattern of optimal paths near the basin boundary.

In this Letter, the $`\rho (x)`$ is calculated over the whole of phase space as a function of noise intensity. Specifically, we address the intricate problem of the noise-induced changes in the pattern of optimal paths near the singularities and show that noise-induced variations of $`\rho (x)`$ are related to the shift of the singularities and of the most probable escape paths. The theoretical predictions are verified by Monte Carlo simulations and the numerical integration of the Fokker-Planck equation (FPE). In the boundary region we find the position of singularities analytically, and provide a more accurate calculation of the oscillating probability distribution discussed in the recent publications Smelyanskiy:99a ; Hanggi:00a ; Maier:01a , and suggest an intuitively clear physical picture of these oscillations.

To see why the computation of the prefactor calls for the analysis of the noise-induced topological changes in the pattern of optimal paths, let us consider a Langevin equation of the form

$$\dot{x}_i=K_i(x)+\sqrt{ϵ}\sigma _{ij}\xi _j(t)i,j=1,2.$$
(2)
Here $`x_i`$ are coordinates of a 2-dimensional dissipative system, $`\xi _j\left(t\right)`$ are Gaussian zero-mean uncorrelated sources of noise linearly mixed by a given noise matrix $`\sigma _{ij}`$, and $`ϵ`$ is a small parameter of the theory.

The evolution of the probability density $`\rho (x,t)`$ for the system (2) is described by the FPE

$`{\displaystyle \frac{\rho }{t}}={\displaystyle \frac{}{x_i}}\left(\rho K_i+{\displaystyle \frac{ϵ}{2}}{\displaystyle \frac{}{x_j}}\left[Q_{ij}\rho \right]\right),`$ (3)
where $`Q=\sigma ^T\sigma `$ is a positively defined diffusion matrix.

To the leading order of approximation in $`ϵ`$, $`S`$ satisfies Ludwig:75a ; Freidlin:84a  a Hamilton-Jacobi equation for a classical action in the form $`H(x,S)=0,`$ where $`H(x,p)`$ is the Hamiltonian function equals $`\frac{1}{2}Q_{ij}p_ip_j+K_i(x)p_i`$. The pattern of extreme trajectories emanating from the stationary state of the system (2) is found by integrating Hamiltonian equations with appropriate initial conditions (see e.g. Ludwig:75a ), and the action is found as an integral along these trajectories

$`\dot{x}_i=K_i+Q_{ij}p_j,\dot{p}_i={\displaystyle \frac{K_j}{x_i}}p_j,\dot{S}={\displaystyle \frac{1}{2}}Q_{ij}p_ip_j.`$ (4)
In general more than one trajectory can arrive at the same point $`x`$ Dykman:94a , but in the limit of $`ϵ0`$ the only trajectories contributing to the $`\rho (x)`$ are those that provide the global minimum of the action and form the pattern of optimal paths.

As a result of global minimization the state space of 2-dimensional system (2) is separated by so-called switching lines into areas to which the system arrives along optimal paths of topologically different types Dykman:94a . In particular, the escape from a metastable state of the system (2) in the $`ϵ0`$ limit is governed by the properties of the so-called most probable escape path (MPEP) that emanates (at $`t\mathrm{}`$) from a stationary state of (2) and arrives (at $`t+\mathrm{}`$) at the boundary of its basin of attraction. The pattern of optimal paths, switching line, and the most probable escape paths are shown in Fig. 1(a) for a periodically driven overdamped motion in a bistable Duffing (PDOD) system: $`K_1=xx^3+A\mathrm{sin}\left(\mathrm{\Omega }t\right)`$, $`K_2=\mathrm{\Omega }=2\pi /T`$, $`\sigma _{ij}=\sqrt{2}\delta _{i1}\delta _{j1}`$ ($`x_1=x`$ and $`x_2=\mathrm{\Omega }t`$).

The pattern is periodic in time and the behavior of the switching line and of the MPEP near the basin boundary is described by the scaling properties of the solution of the linearized Hamiltonian equations (c.f. Smelyanskiy:99a ; Hanggi:00a ; Maier:01a ; Smelyanskiy:97b ) in a Poincarè section (see Fig. 2) in the form

$`p_n=y_0^{}W_0a^n,`$ (5)
$`y_n=y_0^{}a^n+p_nW_0,`$ (6)
where $`y_n=x(t)x_u(t)`$ is the distance to the limit cycle at time $`t=t_0+nT`$, $`a=e^{\lambda T}`$, $`\lambda =1/T_t^{t+T}_{x_1}K_1dt^{}`$, and $`W_0=(1a^2)\left[2_t^{t+T}a^2(t^{},t)𝑑t^{}\right]^1`$ is the second derivative of the action in the $`x`$ direction ($`W_0=_y_yS`$ see eq. (8) below).

The first term in the RHS equation (6) describes the divergence of the zero-momentum trajectory from $`x_u(t)`$ corresponding to the motion on the unstable manifold of this cycle $`_u^{(u)}`$. The second term corresponds to the motion on the stable manifold $`_u^{(s)}`$ of $`x_u(t)`$ with non-zero momentum. We note that the initial position of the most probable escape path $`y_0^{}=x(t_0)x_u(t_0)`$ is the only global parameter of the linear problem. It defines the intersection of $`_u^{(s)}`$ and the unstable manifold of the stable limit cycle $`_s^{(u)}`$ (see Fig. 2(a)). The latter is a set of horizontal lines. The physically observable part of $`_s^{(u)}`$ formed by the optimal paths is a set of segments (see Fig. 2(a)) limited by the switching line.

The position of the switching line is found by noting that the cross-section of the action $`S(y)`$, as a function of the distance $`y`$ from the cycle, is a piecewise linear approximation of an ideal parabola $`S_{esc}\frac{1}{2}y^2W_0`$ Graham:84a ; Maier:01a ; Smelyanskiy:97b , where $`S_{esc}`$ is an escape action in zero-noise limit. The MPEP corresponds to the points where $`S(y)`$ touches the parabola. The switching line corresponds to the points of singularities of $`S(y)`$.

It follows from this geometrical picture that at $`ϵ=0`$: (i) MPEP is the only optimal path that approaches the unstable limit cycle at $`t\mathrm{}`$, (ii) the position of the switching line scales as $`y_n^{sl}=(a+1)y_0^{}a^n/2`$ and (iii) the switching line also approaches the unstable limit cycle at $`t\mathrm{}`$ in agreement with the qualitative picture shown in Fig. 1(a). A similar picture can be obtained for the inverted van der Pol oscillator (see Maier:01a (b)) and is generic for 2D systems with unstable limit cycles.

In the presence of the finite noise intensity, this picture breaks down. As we have mentioned above, one has to assume that, even for small noise intensity, more than one trajectory contributes significantly to the probability of crossing the boundary Hanggi:00a . The choice of a right set of trajectories that cross the boundary is a subtle problem Maier:01a  owing to the fact Graham:84a  that extreme trajectories intersect one another wildly near the unstable limit cycle. For example, in Smelyanskiy:99a  a contribution to the escape from all the trajectories that emanate from the steady state is estimated analytically. In Hanggi:00a ; Maier:01a  it was suggested that the main contribution to the escape comes from an infinite discrete set of trajectories that are small perturbations of the MPEP, which is found numerically. Neither of these choices is consistent with the existence of the switching line, and the methods of Smelyanskiy:99a ; Hanggi:00a ; Maier:01a  cannot be applied to compute noise-induced variation of the $`\rho (x)`$ over the whole of phase space (cf. Bandrivskyy:02 ).

We therefore suggest that, to analyze the distribution in the whole phase space, it is necessary to investigate the noise-induced changes in the pattern of optimal paths. Such an analysis can be conveniently performed in the next-to-leading order of approximation in $`ϵ`$ of the solution of the FPE, yielding Ludwig:75a  an equation for the prefactor in (1) and for the Hessian matrix $`S_{ij}_i_jS`$ Maier:97a

$`{\displaystyle \frac{dz}{dt}}=z\left(_iK_i+{\displaystyle \frac{1}{2}}Q_{ij}S_{ij}\right),`$ (7)
$`\dot{S}_{ij}=p_m_i_jK_mS_{im}_jK_m`$
$`S_{jm}_iK_mS_{jm}S_{ik}Q_{km}.`$ (8)
The second derivative of the action vanishes rapidly near $`x_u(t)`$ on $`_s^{(u)}`$, corresponding to the fact that $`_s^{(u)}`$ is a set of horizontal lines near $`x_u(t)`$. Therefore (for the PDOD) the noise induced correction to the action $`ϵ\mathrm{ln}(z(t))`$ over one period can be estimated as $`ϵ\lambda T`$. As a result, a straightforward application of the equations (4), (7), (8) leads to a discontinuity in the probability distribution at the position of the switching line of size $`ϵ\lambda T`$.

This discontinuity is unphysical. It can be removed if one defines the position of a generalized switching line as a set of points: $`z_1\mathrm{exp}(S_1/ϵ)=z_2\mathrm{exp}(S_2/ϵ)`$ for which probabilities of arrival along topologically different optimal paths are equal. We note that this is a natural generalization of the original definition Dykman:94a  and the latter is recovered when $`ϵ`$ goes to zero. It has physical significance of being the locus of experimentally observable points where the curvature of the probability distribution undergoes sharp qualitative changes (cf. Luchinsky:97c ). The noise-induced changes in the physically observable part of the $`_s^{(u)}`$ and in the $`ϵ\mathrm{ln}\rho (x)`$ are shown in Fig. 2. Note the overall shift of the $`ϵ\mathrm{ln}\rho (x)`$ due to the prefactor, and a stepwise shifts of the linear segments by $`ϵ\lambda T`$ without the change of their slope. The geometrical analysis of these changes gives the position of the modified switching line in the form $`y_n^{sl}=(a+1)y_n/2ϵTW_0/((1a)y_n)`$.

The change in the position of the switching line has profound consequences for the analysis of the non-equilibrium distribution and the escape problem for the finite noise intensities. Firstly, $`\rho (x,t)`$ calculated in the next to the leading order of WKB approximation is now continuous everywhere in phase space. Secondly, the modified switching line has to shift “up” (compare Fig. 1 and Fig. 2), and for the linear regime it crosses the boundary at the point where the distance from the original MPEP to the boundary is $`y_n^{(cr)}=\sqrt{2ϵTW_0/(1a^2)}`$. Since the new switching line crosses the boundary, infinitely many optimal paths can now cross the boundary. Moreover, the original MPEP ceases to be the most probable escape path for a finite noise intensity, since it is now “cut off” by the intersection with the new switching line. Its role is taken by a different path that crosses the boundary. Furthermore, the dominant contribution to the escape from a metastable state is now governed by the continuous set of trajectories which arrive to the boundary through a “corridor”between the modified switching line and the original MPEP. The modified pattern of optimal paths, including the modified MPEP and the modified switching line, are shown in Fig. 1(b). For small noise intensities, the set of optimal trajectories that cross the boundary can be treated as small perturbations of the original MPEP.

We note that the values of the prefactor $`z`$ and the Hessian computed along the paths diverge as a trajectory approaches a caustic. The approximation (7), (8) breaks down here, and the method described in this paper cannot be applied. The size of the breakdown region scales as $`O(ϵ^{2/3})`$ near caustic and as $`O(ϵ^{3/4})`$ near cusp (see e.g. Schulman:81 ; Dykman:94a ; Maier:97a ). Near the singularities $`\rho (x)`$ can be calculated using the Maslov-WKB approximation (see e.g. Dykman:94a ; Maier:97a ). The results of these calculations should merge in the “far filed”the approximation found by the method described above. However, in practice, when calculating $`Sϵ\mathrm{ln}(z)`$ using eqs. (4), (7), (8), the problem of divergence can be avoided by partitioning the state space into boxes with the size that scales as $`ϵ^{2/3}`$ and by assigning to each box the probability which is maximized over the set of trajectories that pass through this box.

To verify theoretical predictions we have performed Monte Carlo simulations and numerical integration of the FPE for two model systems: the PDOD and the inverted van der Pol oscillator $`K_1=x_2`$, $`K_2=x_12\eta x_2\left(1x_1^2\right)`$, $`\sigma _{ij}=\sqrt{4\eta }\delta _{i2}\delta _{j2}`$. Qualitatively similar behavior was obtained for both systems. In Fig. 3(a,b) we compare theoretical predictions for the function $`Sϵ\mathrm{ln}(z)`$ with the results of numerical calculations for the PDOD in different cross-sections. Both numerical techniques produce the same results and they are in good quantitative agreement with the theory for a wide range of noise intensity.

The behavior of $`ϵ\mathrm{ln}(\rho (x))`$ at the boundary on the metastable state is shown in Fig. 3(b) for various values of the noise intensity. It can be seen that the amplitude of the finite-noise induced oscillations increases with $`ϵ`$. The maxima of these oscillations correspond to the intersection of the boundary with the modified switching line and are the nondifferentiability points of the modified action surface $`Sϵ\mathrm{ln}(z)`$. The minima of the distribution correspond to the intersection of the boundary by the modified MPEP, and in the linear regime are shifted with respect to the maxima by $`\mathrm{\Delta }t\mathrm{ln}[2\lambda T/(1a_0^2)]/2\lambda `$. For arbitrary values of parameters the position of the modified switching line and MPEP found in simulations coincide with those predicted theoretically. As the noise intensity $`ϵ`$ decreases, the pattern of optimal paths and the modified action approach their zero-noise limits.

In conclusion, we have shown that the structure of singularities of non-equilibrium potential in the zero-noise limit can be modified in a self consistent way to include the effect of finite noise intensity. As a result, a simple numerical scheme can be used to calculate $`\rho (x)`$ over the whole of phase space in good quantitative agreement with the results of simulations and numerical integration of the FPE for two systems: inverted van der Pol and periodically driven overdamped Duffing oscillators.

In the boundary layer we find analytically the position of the modified switching line and modified MPEP, and give an intuitively clear physical interpretation to the oscillations of the probability at the boundary discussed in recent publications Smelyanskiy:99a ; Hanggi:00a ; Maier:01a . In particular, these results suggest that the oscillations of the escape rate as a function of the noise intensity predicted in Maier:01a  cannot be observed. The method suggested in this Letter can readily be extended to higher dimensional systems.

The work was supported by the Engineering and Physical Sciences Research Council (UK), the Joy Welch Trust (UK), the Russian Foundation for Fundamental Science, and INTAS.