Very massive runaway stars from three-body encounters

Abstract

1 INTRODUCTION

There is growing observational evidence that most (and possibly all) massive stars are formed in clusters (Lada & Lada 2003; cf. de Wit et al. 2005; Schilbach & Röser 2008; Gvaramadze & Bomans 2008b) and that all (or most) O-type stars (either in clusters or in the field) are (or were) members of binary or multiple systems (Mason et al. 1998; Preibisch, Weigelt & Zinnecker 2001; García & Mermilliod 2001; Kobulnicky & Fryer 2007; Clark et al. 2008). Observations also show that the binary frequency increases with stellar mass (Larson 2001; Clark et al. 2008) and that the most massive binaries are usually the most short-period (tight) ones (e.g. Mermilliod & García 2001). Moreover, a high proportion of massive binaries have mass ratios close to unity (Clarke & Pringle 1992; Pinsonneault & Stanek 2006; Kobulnicky & Fryer 2007) and this proportion is larger in close binaries (Mason et al. 1998). The best examples of very massive and tight binaries with companions of comparable mass are the WN6ha binary systems WR 20a (83 M_⊙+ 82 M_⊙; Bonanos et al. 2004; Rauw et al. 2005) and NGC 3603-A1 (116 M_⊙+ 89 M_⊙; Schnurr et al. 2008a).

The WR 20a binary is also remarkable by its significant displacement (∼1 pc) from the centre of the parent cluster Westerlund 2. This displacement strongly suggests that the binary experienced a dynamical encounter with another massive star (either single or binary) and consequently recoiled (cf. Rauw et al. 2005; see also Section 6). The purpose of this paper is to study numerically the role of three-body dynamical encounters between tightly bound (hard; Aarseth & Hills 1972; Heggie 1975; Hills 1975) binaries and a third star in the production of very massive runaway stars (e.g. λ Cep, BD+43°3654, Sk −67°22, BI 237, 30 Dor 016; see Section 2 for a summary of these stars) whose high peculiar velocities can hardly be explained within the framework of the binary–supernova scenario (Section 3). In Section 4 we discuss the dynamical ejection scenario. In Section 5 we present the results of our numerical experiments. The obtained results are discussed in Section 6.

2 VERY MASSIVE RUNAWAY STARS

The stars with peculiar velocities exceeding 30 km s⁻¹ are called runaway stars (Blaauw 1961). The origin of these stars can be attributed to two basic processes: (i) disruption of a short-period binary system following the (asymmetric) supernova explosion of one of the binary components (Blaauw 1961; Stone 1991) and (ii) dynamical three- or four-body encounters in dense stellar systems (Poveda, Ruiz & Allen 1967; Gies & Bolton 1986). Observations show that the percentage of runaway stars is highest (∼25 per cent) among the O stars and steeply decreases to several per cent for the B stars and to even smaller values for the less massive stars (e.g. Gies 1987; Blaauw 1993; Zinnecker & Yorke 2007). This tendency is consistent with the fact that the massive stars prefer to reside in the cores of the parent clusters (either due to dynamical or due to primordial mass segregation) where the dynamical encounters between the cluster members are most frequent and energetic, and suggests that the production of runaway stars is dominated by the second process. On the other hand, in both processes the less massive stars could be accelerated to larger velocities, which is consistent with the observed anticorrelation between the mass and the velocity of runaway stars (Gies & Bolton 1986). The record-holder among the Galactic massive runaway stars is the early B-type (≃11 M_⊙) star HD 271791, whose peculiar velocity (≃530–920 km s⁻¹; Heber et al. 2008) is about an order of magnitude larger than that of BD+43° 3654 – the fastest known Galactic early O-type star (the linear momenta of both stars, however, are comparable).

BD+43°3654 is an O4If (Comerón & Pasquali 2007) runaway star ejected from the Cyg OB2 association about 1.8 Myr ago with a peculiar transverse velocity of ≃40 km s⁻¹ (Comerón & Pasquali 2007; cf. Gvaramadze & Bomans 2008a). Like many other runaway stars, BD+43°3654 generates a bow shock visible in the infrared; see van Buren & McCray (1988) and Comerón & Pasquali (2007) for the Infrared Astronomical Satellite (IRAS) and the Midcourse Space Experiment (MSX) satellite images of the bow shock, respectively. The heliocentric radial velocity of BD+43°3654 of −66 ± 9 km s⁻¹ (Kobulnicky, Gilbert & Kiminki 2010) is much larger than the mean systemic velocity of Cyg OB2 of −10.3 ± 0.3 km s⁻¹ (Kiminki et al. 2007), which also supports the runaway nature of the star. With an initial mass of ≃70 ± 15 M_⊙ (Comerón & Pasquali 2007) and a total peculiar velocity of ≃70 km s⁻¹, BD+43°3654 is the most massive known runaway star in the Galaxy.

Another example of very massive Galactic runaway stars is the O6I(n)f (Walborn 1973) star λ Cep (HD 210839). Like BD+43°3654, λ Cep produces a bow shock, originally discovered with IRAS by van Buren & McCray (1988). In Fig. 1 we present for the first time the Spitzer Space Telescope 24-μm image showing the fine structure of the bow shock. (The image was retrieved from the Spitzer archive using the Leopard software.) Using the parallax and the proper motion from the new reduction of the Hipparcos data (van Leeuwen 2007), one finds the peculiar velocity of λ Cep in Galactic coordinates: v_l≃−32 km s⁻¹, v_b≃−7 km s⁻¹[we used here the Galactic constants R₀ = 8.4 kpc and Θ₀ = 254 km s⁻¹ (Reid et al. 2009) and the solar peculiar motion (U_⊙, V_⊙, W_⊙) = (10.0, 11.0, 7.2) km s⁻¹ (McMillan & Binney 2010)]. The orientation of the bow shock and the direction of the transverse peculiar velocity are consistent with the possibility that λ Cep was ejected about 2.5 Myr ago from the Cep OB3 association, located at to the east from the star (cf. Hoogerwerf, de Bruijne & de Zeeuw 2001). The current mass of λ Cep is ≃45–60 M_⊙ (Repolust, Puls & Herrero 2004; Martins, Schaerer & Hillier 2005), which implies that the star was ejected very soon after its birth in the association. Adopting the heliocentric radial velocity of λ Cep of ≃−76 km s⁻¹ (Conti, Leep & Lorre 1977) and the systemic velocity of Cep OB3 of ≃−23 km s⁻¹ (Mel'nik & Dambis 2009), one finds a total peculiar velocity for the star of ≃60 km s⁻¹.

Figure 1

Spitzer 24 μm image of the bow shock associated with the O6 I(n)f runaway star λ Cep. The position of the star is indicated by a circle.

Open in new tabDownload slide

Finally, we note several very massive (O2-type) stars in the Large Magellanic Cloud (LMC) whose high radial velocities (∼40–150 km s⁻¹ greater than the LMC systemic velocity) were interpreted as an indication that these stars are runaways (Massey et al. 2005; Evans et al. 2006, 2010; cf. Nota et al. 1994; Danforth & Chu 2001; Schnurr et al. 2008b). A strong support for this interpretation comes from the discovery of bow shocks associated with one of these candidate runaway stars, BI 237, and several other OB stars in the field of the LMC (Gvaramadze, Kroupa & Pflamm-Altenburg 2010; see also Section 6). The most striking runaway in the LMC is the O2III-If* star 30 Dor 016 (Evans et al. 2010), which is located in the periphery of the star-forming complex 30 Doradus, at ≃8 arcsec (≃120 pc in projection) from R136 (the central cluster of 30 Doradus) or ≃5 arcsec (≃70 pc) from NGC 2060 (another cluster in 30 Doradus). The very large mass of 30 Dor 016 (∼90 M_⊙; Evans et al. 2010) makes this star the most massive known runaway.

The summary of the very massive runaway stars in the Galaxy and the LMC is given in Table 1. The last column gives the birthplaces of the stars.

Table 1

Open in new tab

Very massive runaway stars. For the first two (Galactic) stars we give their full (three-dimensional) peculiar velocities, while for the remaining four stars (located in the LMC) we give the peculiar radial velocities only.

Star	Sp. type	v (⁠⁠)	Association
BD+43°3654	O4If	≃70(1)	Cyg OB2
λ Cep	O6I(n)f	≃60(2)	Cep OB3
N11-026	O2.5 III(f*)	≃35(3)	LH 10 ?
Sk −67°22	O2 If*	≃150(4)	?
BI 237	O2 V((f*))	≃120(4)	LH 82
30 Dor 016	O2III-If*	≃85(5)	30 Doradus

Star	Sp. type	v (⁠⁠)	Association
BD+43°3654	O4If	≃70(1)	Cyg OB2
λ Cep	O6I(n)f	≃60(2)	Cep OB3
N11-026	O2.5 III(f*)	≃35(3)	LH 10 ?
Sk −67°22	O2 If*	≃150(4)	?
BI 237	O2 V((f*))	≃120(4)	LH 82
30 Dor 016	O2III-If*	≃85(5)	30 Doradus

(1) Comerón & Pasquali (2007), Kobulnicky et al. (2010); (2) van Leeuwen (2007), Conti et al. (1977); (3) Evans et al. (2006); (4) Massey et al. (2005); (5) Evans et al. (2010).

Table 1

Open in new tab

Very massive runaway stars. For the first two (Galactic) stars we give their full (three-dimensional) peculiar velocities, while for the remaining four stars (located in the LMC) we give the peculiar radial velocities only.

Star	Sp. type	v (⁠⁠)	Association
BD+43°3654	O4If	≃70(1)	Cyg OB2
λ Cep	O6I(n)f	≃60(2)	Cep OB3
N11-026	O2.5 III(f*)	≃35(3)	LH 10 ?
Sk −67°22	O2 If*	≃150(4)	?
BI 237	O2 V((f*))	≃120(4)	LH 82
30 Dor 016	O2III-If*	≃85(5)	30 Doradus

Star	Sp. type	v (⁠⁠)	Association
BD+43°3654	O4If	≃70(1)	Cyg OB2
λ Cep	O6I(n)f	≃60(2)	Cep OB3
N11-026	O2.5 III(f*)	≃35(3)	LH 10 ?
Sk −67°22	O2 If*	≃150(4)	?
BI 237	O2 V((f*))	≃120(4)	LH 82
30 Dor 016	O2III-If*	≃85(5)	30 Doradus

(1) Comerón & Pasquali (2007), Kobulnicky et al. (2010); (2) van Leeuwen (2007), Conti et al. (1977); (3) Evans et al. (2006); (4) Massey et al. (2005); (5) Evans et al. (2010).

3 BINARY–SUPERNOVA SCENARIO

We now consider whether or not the origin of very massive runaway stars can be explained within the framework of the binary–supernova scenario. According to this scenario, a runaway star attains its peculiar velocity in the process of disintegration of a binary system following the supernova explosion of the primary (initially more massive) component of the binary (Blaauw 1961). Since we are interested in the production of very massive (≥60–70 M_⊙) runaways, one should assume that the primary star was a very massive star as well, with an initial mass comparable to that of the secondary star. (Note that the secondary star can increase its mass due to Roche-lobe overflow of the primary star so that the mass of the runaway star could be larger than the initial mass of the primary star.)

Stellar evolutionary models suggest that the pre-supernova mass of stars with initial (zero-age main-sequence) masses, M_ZAMS, from 12 to 120 M_⊙ do not exceed ∼10–17 M_⊙ (Schaller et al. 1992; Vanbeveren, De Loore & Van Rensbergen 1998; Woosley, Heger & Weaver 2002; Meynet & Maeder 2003) and that it is maximum for stars with M_ZAMS≃ 20–25 M_⊙ and ≳ 80 M_⊙[see fig. 6 of Meynet & Maeder (2003)]. In the first case, the supernova explosion leaves behind a neutron star, while in the second one the stellar supernova remnant is a black hole of mass of ≃5–10 M_⊙ (e.g. Woosley et al. 2002; Eldridge & Tout 2004). It is clear that a very massive binary cannot be unbound by a symmetric supernova explosion since the system loses much less than a half of its pre-supernova mass (Boersma 1961). One should therefore assume that the stellar supernova remnant, a 5–10 M_⊙ black hole, received at birth a kick velocity exceeding the escape velocity from the system (Stone 1982; Tauris & Takens 1998). In this case, the peculiar velocity of the runaway star strongly depends on the magnitude and the orientation of the kick attained by the black hole (Tauris & Takens 1998; Gvaramadze 2006a, 2009). Moreover, the runaway star can achieve the highest velocity if the pre-supernova binary was as tight as possible, i.e. if the secondary star of radius r₂ was close to filling its Roche lobe, r₂∼r_L, where r_L is the radius of the Roche lobe, given by (Eggleton 1983)

1

where a is the binary semimajor axis; q=m₂/m₁; and m₁ and m₂ are the masses of the primary and the secondary stars, respectively. Assuming that the secondary is a 70 M_⊙ main-sequence star of radius

2

(Habets & Heintze 1981) and adopting the maximum pre-supernova mass of the exploding star of 15 M_⊙, one finds from equations (1) and (2) that a≃ 30 R_⊙ (≃2r₂). Then using equations (44)–(47) and (51)–(56) given in Tauris & Takens (1998), one finds that to produce a 70 M_⊙ runaway star with a velocity of ≃70 km s⁻¹ the black hole should attain a kick velocity of at least ≃250 km s⁻¹ (if the mass of the black hole is 5 M_⊙; see Fig. 2) or ≃400 km s⁻¹ (if the mass of the black hole is 10 M_⊙), while the angle, θ, between the kick vector and the direction of motion of the exploding star should be in a certain range (the smaller the kick the narrower the range of allowed angles). Although one cannot exclude a possibility that a 5–10 M_⊙ black hole can attain a kick of several hundreds of and of appropriate orientation (see Gualandris et al. 2005 for an example of such case), we consider the binary–supernova scenario as highly unlikely (cf. Gvaramadze 2007, 2009; Gvaramadze & Bomans 2008a).

Figure 2

The peculiar velocity of a 70 M_⊙ (secondary) star as a function of the angle, θ, between the kick vector and the direction of motion of the primary (exploding) 15 M_⊙ star and the magnitude of the kick, w, attained by the stellar supernova remnant (black hole) of mass M_BH. Upper panel: M_BH = 5 M_⊙, w = 250 km s⁻¹ (short–dashed line), 300 km s⁻¹ (middle–dashed line), 350 km s⁻¹ (long–dashed line). Bottom panel: M_BH = 10 M_⊙, w = 400 km s⁻¹ (short–dashed line), 500 km s⁻¹ (middle–dashed line), 600 km s⁻¹ (long–dashed line). The discontinuities in the curves correspond to a range of angles θ for which the system remains bound. The horizontal line indicates the peculiar velocity of the runaway star of 70 km s⁻¹. See text for details.

Open in new tabDownload slide

It is obvious that the binary–supernova scenario cannot be applied to the runaway stars ejected from young (≤2–3 Myr) clusters. The most massive stars in these clusters simply have no time to end their lives in supernova explosions. Similarly, the young (∼1–2 Myr) ages of O2-type runaway stars in the LMC are also inconsistent with the binary–supernova scenario (cf. Evans et al. 2010). Moreover, the large separations of the very massive runaway stars (listed in Table 1) from their plausible parent clusters and OB associations imply that the majority of these stars were ejected very soon after the birth, which also argues against the binary–supernova scenario (Gvaramadze et al. 2010).

4 DYNAMICAL EJECTION SCENARIO: MASSIVE RUNAWAY STARS FROM THREE-BODY ENCOUNTERS

An alternative to the binary–supernova scenario is the scenario based on three- and four-body dynamical encounters in dense stellar systems (Poveda et al. 1967; van Albada 1968; Aarseth 1974; Kroupa 1998; Gualandris, Portegies Zwart & Eggleton 2004; Pflamm-Altenburg & Kroupa 2006; Gvaramadze 2007, 2009; Gvaramadze & Bomans 2008a; Gvaramadze, Gualandris & Portegies Zwart 2008, 2009). The best-studied and possibly the most efficient process responsible for the origin of high-velocity stars is the close dynamical encounter between two hard binary stars (Mikkola 1983; Leonard & Duncan 1990). Numerical experiments performed by Leonard (1991) show that in the course of binary–binary encounters one of the binary components can occasionally be ejected with a velocity comparable to the escape velocity from the surface of the most massive star in the binaries, i.e. with a velocity of ∼1000 km s⁻¹.

Three-body encounters can also produce very high velocity runaway stars, provided that the mass ratio of the single star to the mass of the (hard) binary is either ≫ or ≪1 (Hills & Fullerton 1980). In the first case, one of the binary components could be replaced by the very massive star in a so-called exchange encounter, while the second component is ejected with a high velocity (e.g. Gvaramadze, Gualandris & Portegies Zwart 2009). For equal-mass binary components and zero-impact parameter, exchange encounters produce a typical ejection velocity of ∼1.8V_orb, where V_orb is the orbital velocity of the ejected star in the original binary (Hills & Fullerton 1980). In the second case, the low-mass star is scattered by the very massive binary and attains a typical velocity of ∼0.8V_orb in a so-called fly-by encounter (here V_orb is the orbital velocity in the very massive binary with equal-mass components). In both cases, the ejected star gains its kinetic energy at the expense of the increased binding energy of the post-encounter binary. In response to the encounter, the binary recoils with a fraction M₃/(M₁+M₂) of the velocity of the ejected star (M₁+M₂ and M₃ are the masses of the post-encounter binary and the runaway star).

In the following, we will concentrate on dynamical encounters between very massive binaries and single stars of mass comparable to that of the binary components. For illustrative purposes, we assume that the massive binaries have parameters similar to those of the most massive known binary systems in the Galaxy, WR 20a and NGC 3603-A1 (see Table 2). Both systems are very tight, with the semimajor axes of only ≃3 times larger than the radii of the primary stars. The binding energy of these binaries is comparable to the energy of a supernova explosion, ∼10⁵¹ erg. If a 70 M_⊙ intruding star extracts only 1 per cent of this energy, it will attain a peculiar velocity of ∼100 km s⁻¹, which is large enough to explain the observed peculiar velocities of the very massive runaway stars (see Table 1).

5 NUMERICAL EXPERIMENTS

In this section, we perform numerical simulations of three-body encounters in order to obtain the velocity distribution for runaway stars produced in the course of interactions between (massive) single stars and a hard very massive binary. The simulations are carried out with the sigma3 package included in the starlab software environment (McMillan & Hut 1996; Portegies Zwart et al. 2001; http://www.ids.ias.edu/~starlab). The stars are treated as point masses interacting gravitationally. However, the monitoring of the relative distances between pairs of stars combined with the information on the stellar radii allows us to identify collisions. The stellar radii are determined via the mass–radius relationship given by equation (2). If two stars come closer than the sum of their radii, the calculation is stopped and the encounter is classified as a collision (merger).

We consider a target binary composed of stars of mass M₁ and M₂; initial semimajor axis a and eccentricity e; and an intruding star of mass M₃ with an initial velocity V_rel relative to the centre of mass of the binary. V_rel is set to ⁠, in accordance with typical dispersion velocities in young massive clusters. The angles that define the spatial orientation of the binary with respect to the single star are randomized in a Monte Carlo fashion (see Hut & Bahcall 1983). The eccentricity of the binary is drawn from a thermal distribution P(e) = 2e (Heggie 1975), having set a maximum value in order to guarantee that the two binary components do not come into contact at the first pericenter passage. The impact parameter b is randomized according to an equal probability distribution for b² in the range [0 −b_max]. The maximum value b_max is determined automatically for each experiment [see Gualandris et al. (2004) for a description]. Energy conservation is usually better than one part in 10⁶ and, in case the error exceeds 10⁻⁵, the encounter is rejected. The accuracy in the integrator is chosen in such a way that at most 5 per cent of the encounters are rejected.

5.1 WR 20a-like binary

In the first set of simulations we focus on interactions in which a single star of mass M₃ (ranging from 3 to 80 M_⊙) encounters a WR 20a-like binary (M₁=M₂ = 80 M_⊙, a = 55 R_⊙≃ 0.25 au; cf. Table 2).

In Fig. 3 we present the probability of different outcomes (branching ratios) as a function of M₃. For each value of M₃ we perform a total of 3000 scattering experiments, which result either in a fly-by or a merger. Ionizations never take place as the binary is too hard to be dissociated by the incoming star. The small binary separation (comparable to the radii of the binary components) makes exchange encounters very rare so that their contribution to the production of runaway stars is negligible. The only encounters which can produce runaways are fly-by encounters. The probability of these encounters decreases with increasing M₃ and drops to ≃20 per cent for M₃ = 80 M_⊙. Correspondingly, the fraction of mergers increases to ≃80 per cent.

Figure 3

Branching ratio for the outcome of encounters between a (80, 80) M_⊙ binary and a single star as a function of its mass. The different outcomes are: merger (stars), fly-by (circles) and exchange (squares). The error bars represent the formal (1σ) Poissonian uncertainty of the measurement.

Open in new tabDownload slide

Fig. 4 shows the average velocity of escapers produced in fly-by encounters as a function of their mass, M₃. It can be seen that the velocity is a non-monotonic function of M₃: it rapidly grows to its maximum value at M₃ = 10 M_⊙ and gradually decreases afterwards. This behaviour is the result of the interplay between two competing factors: (i) the more massive the intruding star the greater its ability to affect the binding energy of the binary system and thereby to increase the kinetic energy of the system, and (ii) the more massive the intruding star the larger the recoil velocity of the binary and the smaller the velocity of the ejected star relative to the centre of mass. Fig. 4 also shows that the average velocity attained by the 70–80 M_⊙ stars is quite moderate, <30 km s⁻¹, so that they cannot be formally classified as runaways. On the other hand, in 10 per cent of the fly-by encounters the 70–80 M_⊙ stars attain velocities of >70 km s⁻¹, and occasionally (in 1 per cent of the fly-by encounters) can be accelerated to even larger (>150 km s⁻¹) velocities.

Figure 4

Average velocity of escapers as a function of mass. Circles represent the average velocity, diamonds indicate the velocity V_max for which 10 per cent of the encounters have V_∞ > V_max, and stars indicate the velocity V_max for which 1 per cent of the encounters have V_∞ > V_max. The error bars indicate the 1σ deviation from the mean. For clarity, we only show them for one data set.

Open in new tabDownload slide

In Fig. 5 we show the probability of fly-by encounters resulting in ejection velocities from 30 to 100 km s⁻¹ (top to bottom). For M₃≳ 70 M_⊙ about 6–8 per cent of all encounters produce escapers with peculiar velocities of >30 km s⁻¹, i.e. typical of runaway stars, while the higher velocities (>80 km s⁻¹) can be attained in ≃2–3 per cent of the encounters. The percentage of high-velocity (>100 km s⁻¹) runaway OB stars (M₃≥ 8 M_⊙) increases with decrease of M₃ and reaches a maximum (≃20 per cent) for M₃≃ 10 M_⊙.

Figure 5

The probability of fly-by encounters between a (80, 80) M_⊙ binary (with a = 0.25 au) and a single star (of mass M₃) resulting in different ejection velocities: >30 km s⁻¹ (squares), >50 km s⁻¹ (diamonds), >80 km s⁻¹ (stars), >100 km s⁻¹ (circles).

Open in new tabDownload slide

To study the effect of the initial semimajor axis of the very massive binary, we performed further scattering experiments. Fig. 6 shows the branching ratios as function of a for four different values of the mass of the intruding star, M₃ = 20, 40, 60 and 80 M_⊙. One can see that the larger the semimajor axis the larger the percentage of fly-by encounters. Fig. 6 also shows that the wider the binary and the more massive the intruding star the larger the percentage of exchange encounters, which for M₃ = 80 M_⊙ and a = 4 au reaches ≃20 per cent.

Figure 6

Branching ratios for the outcomes of encounters between a (80, 80) M_⊙ binary and a single star of mass M₃ = 20, 40, 60 and 80 M_⊙ as a function of the binary semimajor axis. The different outcomes are: merger (stars), fly-by (circles) and exchange (squares). The error bars represent the formal (1σ) Poissonian uncertainty of the measurement.

Open in new tabDownload slide

In Fig. 7 we show the average velocity of escapers as a function of a. One can see that the velocity increases with increasing a and then (for a≳ 0.4–0.6 au) gradually decreases. This counterintuitive growth can be understood if one takes into account that the percentage of fly-by encounters grows with a as well (Fig. 6), so that for a range of semimajor axes (up to several au) the increasing number of encounters producing runaway stars compensates and even overcomes the velocity decrease caused by the increase of the semimajor axis. Particularly, one can see that the 80 M_⊙ stars attain velocity of ∼100–140 km s⁻¹ in about 2–4 per cent of all encounters, if the binary separation is ≃0.3–3 au. Thus, three-body encounters involving binaries with the mass of WR 20a and the semimajor axes up to several au are quite efficient in producing very massive runaway stars.

Figure 7

Average velocity of escapers as a function of the initial binary semimajor axis in the interaction of a (80, 80) M_⊙ binary star with a single star of different mass: M₃ = 20 M_⊙ (circles), M₃ = 40 M_⊙ (squares), M₃ = 60 M_⊙ (triangles), M₃ = 80 M_⊙ (diamonds). Solid symbols represent the average velocity while the empty symbols indicate the velocity V_max for which 10 per cent of the encounters have V_∞ > V_max. The error bars indicate the 1σ deviation from the mean. For clarity, we only show them for one data set.

Open in new tabDownload slide

5.2 NGC 3603-A1-like binary

In the second set of simulations we consider three-body encounters involving a NGC 3603-A1-like binary (M₁ = 120 M_⊙, M₂ = 90 M_⊙, a = 60 R_⊙≃ 0.3 au; cf. Table 2). The branching ratios for these encounters (Fig. 8) are almost identical to those given in Fig. 3. The only difference is the somewhat larger percentage of fly-by encounters, which is due to the larger semimajor axis of the binary.

Figure 8

Branching ratio for the outcome of encounters between a (120, 90) M_⊙ binary and a single star as a function of its mass. The different outcomes are: merger (stars), fly-by (circles) and exchange (squares). The error bars represent the formal (1σ) Poissonian uncertainty of the measurement.

Open in new tabDownload slide

The average velocity of escapers is shown in Fig. 9. As expected, the larger mass of the binary results in a higher velocities of the escapers. Fig. 9 shows that in 10 per cent of the fly-by encounters (or in ≃2–3 per cent of all encounters) the 70–80 M_⊙ stars attain velocities exceeding 110–120 km s⁻¹. One can also see that the NGC 3603-A1-like binary is most efficient at accelerating the 20 M_⊙ stars, whose average velocity is ≃100 km s⁻¹, while in about 6 per cent of all encounters these stars attain velocities of ≥320 km s⁻¹.

Figure 9

Average velocity of escapers as a function of the mass. Circles represent the average velocity, diamonds indicate the velocity V_max for which 10 per cent of the encounters have V_∞ > V_max, and stars indicate the velocity V_max for which 1 per cent of the encounters have V_∞ > V_max. The error bars indicate the 1σ deviation from the mean. For clarity, we only show them for one data set.

Open in new tabDownload slide

Fig. 10 shows that the wider the very massive binary the more frequent are the fly-by and the exchange encounters. Correspondingly, the average velocity of escapers grows with a and after reaching the maximum value at a≃ 0.4–0.6 au it gradually decreases (see Fig. 11). Still, for binaries as wide as several au, the average velocity of escapers remains higher than that attained in the encounters with the NGC 3603-A1-like binary.

Figure 10

Branching ratios for the outcomes of encounters between a (120, 90) M_⊙ binary and a single star of mass M₃ = 20, 40, 60 and 80 M_⊙ as a function of the binary semimajor axis. The different outcomes are: merger (stars), fly-by (circles) and exchange (squares). The error bars represent the formal (1σ) Poissonian uncertainty of the measurement.

Open in new tabDownload slide

Figure 11

Average velocity of escapers as a function of the initial binary semimajor axis in the interaction of a (120, 90) M_⊙ binary star with a single star of different mass: M₃ = 20 M_⊙ (circles), M₃ = 40 M_⊙ (squares), M₃ = 60 M_⊙ (triangles), M₃ = 80 M_⊙ (diamonds). Solid symbols represent the average velocity while the empty symbols indicate the velocity V_max for which 10 per cent of the encounters have V_∞ > V_max. The error bars indicate the 1σ deviation from the mean. For clarity, we only show them for one data set.

Open in new tabDownload slide

6 DISCUSSION

We performed numerical simulations of dynamical encounters between very massive hard binaries and a single (massive) star in order to explore the possibility that this three-body process is responsible for the origin of very massive (≥60–70 M_⊙) runaway stars, whose high peculiar velocities (≥70 km s⁻¹) cannot be easily produced by the disintegration of a binary system following a supernova explosion. For illustrative purposes, we considered encounters with binaries whose parameters are similar to those of the two most massive known Galactic binaries, WR 20a and NGC 3603-A1 (see Table 2). Our simulations were motivated by the observational fact that one of these binaries (WR 20a) is significantly offset (∼1 pc) from the centre of the parent cluster (Westerlund 2), which strongly suggests that the binary was involved in a dynamical encounter with another massive star and thereby was kicked out of the cluster (cf. Rauw et al. 2005). We estimated the typical velocities produced in encounters between very massive binaries and single stars and found that ≃10 per cent of the fly-by encounters (or ≃2 per cent of all encounters) between the WR 20a-like binary and a 70–80 M_⊙ star produce escapers with velocities (>70 km s⁻¹) similar to those of the very massive Galactic runaway stars, λ Cep and BD+43°3654. We also found that in about 2 per cent of all encounters between the NGC 3603-A1-like binary and a 80 M_⊙ star the escaper attains a velocity of ≳120 km s⁻¹, which is comparable to that of the most massive (∼90 M_⊙) runaway star in the LMC, 30 Dor 016 (Evans et al. 2010; see also Table 1). The ejection velocities could be even higher if the semimajor axes of the very massive binaries were larger than those of WR 20a and NGC 3603-A1. In about 2–5 per cent of encounters involving binaries with semimajor axes in the range from 0.3 to ∼4 au, the ejection velocity of 80 M_⊙ stars is ≥100–160 km s⁻¹. We therefore argue that the origin of (at least) some very massive high-velocity runaway stars in the Galaxy and the LMC is associated with dynamical three-body encounters.

Production of high-velocity massive stars via three-body encounters is accompanied by recoil of the very massive binary in the opposite direction to the ejected star. If the very massive runaway star 30 Dor 016 were ejected in the field via the three-body encounter in the central cluster, R136, of the 30 Doradus nebula, then one would expect to find a very massive binary on the opposite side of the cluster. Interestingly, such a binary does indeed exist. The very massive binary R145 (HD 269928), whose mass is of the same order of magnitude as those of WR 20a and NGC 3603-A1 (Schnurr et al. 2009, also private communication), is located at (or ≃19 pc in projection) from R136, just on the opposite side of 30 Dor 016 (see Fig. 12).1 If one assumes that 30 Dor 016 and R145 were ejected from R136 owing to the same three-body encounter, then the conservation of the linear momentum implies that the mass of the binary should be ≃570 M_⊙, which is too large to be realistic (see Schnurr et al. 2009). From this it follows that either another very massive binary exists at a larger distance from R136 or 30 Dor 016 attained its peculiar velocity in the course of a binary–binary encounter (cf. Gvaramadze & Bomans 2008a). In the latter case, 30 Dor 016 and other stars involved in the encounter should not lie on the same line.

Figure 12

Spitzer 24 μm image of the 30 Doradus star-forming complex with position of its central cluster, R136, marked by the diamond point. The positions of the very massive binary R145 and two runaway stars, 30 Dor 016 and Sk −69°206, are marked by circles. See text for details.

Open in new tabDownload slide

Similarly, the large offset of R145 from R136 could be interpreted as an indication that the binary was involved in an energetic gravitational interaction in the parent cluster and that a massive runaway star was ejected in the opposite direction. Fig. 12 shows that the B2 star (Rousseau et al. 1978) Sk −69°206 could be such a runaway. This star, located ≃17 arcmin to the west of R136, was identified as a runaway via detection of its associated bow shock, whose orientation is consistent with the possibility that Sk −69°206 was ejected from 30 Doradus (Gvaramadze et al. 2010). Assuming that the mass of Sk −69°206 is ∼10–15 M_⊙,2 one finds that R145 should be as massive as 130–200 M_⊙, which is consistent with the mass estimate given in Schnurr et al. (2009, also private communication).

The presence of numerous very massive (O2-3 and WN6h) stars spread all around 30 Doradus suggests that despite the young age (1–2 Myr) of R136, the cluster has already experienced a violent dynamical evolution during which it lost a significant fraction of its massive (single and binary) stars (cf. Brandl et al. 2007; see also Pflamm-Altenburg & Kroupa 2006; Moeckel & Bate 2010). We therefore predict that some of the very massive stars in the 30 Doradus region are binary systems recoiled from R136 due to three-body encounters in the cluster's core. The spectroscopic monitoring of the brightest stars around 30 Doradus would allow us to reveal the radial velocity variability and thereby to identify binaries among them (e.g. Schnurr et al. 2008b), while the future proper motion measurements for these stars with the space astrometry mission Gaia will allow us to determine the timing of their ejection and thereby to link several stars to the same ejection event.

Our simulations also showed that the dynamical three-body encounters provide an efficient channel for production of high-velocity early B-type stars – the progenitors of the majority of neutron stars (pulsars). We found that ≃6–8 per cent of all encounters between 10–20 M_⊙ single stars and very massive binaries produce escapers with velocities (>200–350 km s⁻¹) typical of pulsars (e.g. Hobbs et al. 2005), and thereby could contribute to the origin of peculiar velocities of these objects (cf. Gvaramadze 2006b, 2007; Gvaramadze et al. 2008; Gvaramadze & Bomans 2008b).

We also found that in ≃1 per cent of encounters involving early B-type stars, the ejected star attains a velocity of ≥550–600 km s⁻¹. It is worth noting that the velocity of this order of magnitude was measured for the 11 ± 1 M_⊙ star HD 271791 – the fastest known massive runaway in the Galaxy (Heber et al. 2008). It was shown by Gvaramadze (2009, 2010) that the most likely origin of this extremely high velocity star is through dynamical interaction in the dense core of the parent star cluster. Three-body encounters discussed in the present paper could be one of the possible dynamical processes responsible for the origin of HD 271791.

A by-product of our simulations is the finding that about 1 per cent of all encounters involving low-mass (3–5 M_⊙) stars produces escapers with velocities of >500–600 km s⁻¹, typical of the so-called hypervelocity stars – the ordinary stars moving with velocities exceeding the Milky Way's escape velocity (Brown et al. 2005). The existence of the hypervelocity stars was foreseen by Hills (1988), who showed that a close encounter between a tight binary and the supermassive black hole in the Galactic Centre can produce escapers with a velocity of up to several 1000 km s⁻¹. An alternative explanation of the origin of hypervelocity stars is that they attain extremely high velocities via strong dynamical three- or four-body encounters in the dense cores of massive star clusters located in the Galactic disc (Gvaramadze 2006b, 2007, 2009; Gvaramadze et al. 2008, 2009) or in the Large Magellanic Cloud (Gualandris & Portegies Zwart 2007). The three-body encounters between single late B-type stars and a very massive binary provide an additional channel for production of the hypervelocity stars.

To conclude, we note that numerous uncertainties about the initial conditions and early dynamical evolution of young star clusters precludes us from making any estimates of the production rate of high-velocity runaway stars (cf. Gvaramadze et al. 2008, 2009).

We are grateful to P. Kroupa, N. Langer, J. Pflamm-Altenburg, O. Schnurr and T. M. Tauris for useful discussions, and to P. P. Eggleton (the referee) for useful suggestions on the manuscript. VVG acknowledges the Deutsche Forschungsgemeinschaft for financial support. This research has made use of the NASA/IPAC Infrared Science Archive, which is operated by the Jet Propulsion Laboratory, California Institute of Technology, under contract with the National Aeronautics and Space Administration.

REFERENCES