Free Access
Issue
A&A
Volume 563, March 2014
Article Number A96
Number of page(s) 9
Section Interstellar and circumstellar matter
DOI https://doi.org/10.1051/0004-6361/201322308
Published online 17 March 2014

© ESO, 2014

1. Introduction

Runaway massive stars have large spatial velocities (v > 30 km s-1) (e.g., Gies & Bolton 1986; Tetzlaff et al. 2011). They move supersonically through the interstellar medium (ISM) forming a bowshock pointing in the direction of the star velocity (e.g., Van Buren et al. 1995). The shocked material is heated by the stellar radiation field, and the swept dust re-emits at infrared (IR) wavelengths (e.g., Van Buren & McCray 1988; Kobulnicky et al. 2010).

The high proper velocities of runaway stars can be produced by the supernova explosion of a presumed binary companion (Blaauw 1961), or by dynamical ejection (Leonard & Duncan 1988). Recently, Fujii & Zwart (2011) argued that the velocity originates from strong gravitational interactions between single stars and binary systems in the centres of stellar clusters. Runaway massive stars are usually present around young star clusters.

Bowshocks of early-type massive stars might produce non-thermal radiation (see del Valle & Romero 2012), a fact supported by several observations. Radio non-thermal emission was detected from the bowshock of the runaway star BD+43°3654 (Benaglia et. al 2010). Recently, the bowshock of the star AE Aurigae was detected at X-ray energies, and the emission is well described by a power-law spectrum that can be modelled as inverse Compton (IC) up-scattering of IR photons (López-Santiago et al. 2012). Finally, the bowshock of the well-known massive star HD 195519 has been associated with a Fermi source (see del Valle et al. 2013).

Runaway stars can eventually travel through the molecular cloud (MC) where they were formed, and interact with dense structures. The bowshockmedium interactions produce variable non-thermal emission as the star moves through the MC. In this work, we propose that runaway early-type stars moving within MCs can be variable gamma-ray sources. These objects might be counterpart of some of the unidentified variable gamma-ray sources, concentrated towards the galactic plane. We adopt the model developed by del Valle et al. (2012) to compute spectral energy distributions and light curves for different types of stars. We also calculate the photon absorption along the whole spectrum.

In the next section, we briefly describe bowshocks of runaway stars. In Sect. 3, we briefly discuss the molecular cloud structure we adopted. In the next section, we present the non-thermal emission and absorption calculations we implemented for O4I and O9I stars. The main results are also given in this section. Finally, in Sect. 5, we discuss the variability of the emission and we offer our conclusions.

2. Bowshocks of massive runaway stars

Significant research has been done on bowshock modelling (e.g., Van Buren & McCray 1988; Van Buren et al. 1990; Bandiera 1993; Brighenti & Dércole 1995; Chen et al. 1996; Wilkin 1996, 2000, Comerón 1997; Chen & Huang 1997; Comerón & Kaper 1998; Wareing et al. 2007). The collision of a stellar wind, of mass-loss rate w, density ρw and terminal velocity Vw, with the ISM, of density ρa, around a runaway star results in a system of two shocks. The ram pressure of the wind and the ISM balances at some distance from the star, i.e., , where ρw = w/4πR2Vw. Here R is the radial distance from the star. The value of R where this occurs is defined as the standoff radius R0: (1)As can be seen from Eq. (1), the standoff radius decreases for a denser ambient medium. This is because the ISM ram pressure becomes stronger.

Bowshocks of runaway stars are imaged in the IR because of the emission produced by the heated gas and dust that they sweep. The heating of the material can be produced by the UV emission of the runaway star or with the radiation of the shocked gas in the post-shock region. A simple energetic analysis shows that the former dominates by at least an order of magnitude (e.g., Van Buren & McCray 1988). The kinetic power of the stellar wind is erg s-1. The available power for heating the gas is a fraction ξ of the wind power, i.e., ξ Lw. On the other hand, the power from the star luminosity is 3.84 × 1033 (L/L) erg s-1. Again, the available power for heating the gas and dust with this radiation is a fraction of L. The power shock/radiation ratio for the case of an O4I star is Lw/LUV ~ 1.52 × 1038/2.69 × 1039 < 1 (see the stellar parameters in Table 1; the factor ξ is the same in both cases). It is clear then that the dominant heating mechanism is stellar radiation. In comparison, the radiation from the shocked gas plays a minor role.

In the case of a concrete star, such as ζ Oph (a well known runaway star, e.g., Peri et al. 2012), Lw~3.14 × 1035(10-7)(1.5 × 103)2 erg s-1~7 × 1034 erg s-1. The observed IR luminosity from the bowshock reaches values of 5.5 × 1035 erg s-1, i.e., it is higher than even the mechanical power of the wind, showing that the stellar luminosity is the main heating source (e.g., Povich et al. 2008).

The characteristic temperature of dust emission, TIR, can be estimated using a simplified dust model proposed by Draine & Lee (1984). Given a radiation field and the grain absorption efficiency, the dust temperature TIR can be computed by equating the dust heating by absorption with the dust cooling by emission. For the predominant UV radiation field, and using a dust emissivity law of the form j = λ-2B(T), where λ is the wavelength and B(T) is Planck’s emission law (Van Buren & McCray 1988), we find: (2)Here aμm~0.2μm is the dust grain radius, R0 pc is R0 in pc, and L 38 is the star luminosity in units 1038 erg s-1. More complex dust emission models can be found in Draine & Li (2007) and Draine (2011).

A fraction of the star’s bolometric luminosity is re-emitted in the IR by the dust grains. The re-processed luminosity, LIR, can be roughly estimated as that of a black body at T = TIR.

From Eqs. (1) and (2), it follows that , and . Therefore, the shape, luminosity and temperature of the bowshock depend on the ambient density.

3. Molecular cloud structure

Molecular clouds are the site of practically all star formation in the Galaxy. Typical densities are ~102 − 103 cm-3 (Crutcher et al. 2010). Here we adopt a value of nMC ~ 102 cm-3. These clouds have varied structures on different length scales. The clouds collapse to form dense cores1 through a combination of gravity and turbulence: gravoturbulent fragmentation (e.g., Klessen 2011). The higher-density structures have typical temperatures of the order of ~10 K and densities ~104 − 105 cm-3 (Bodenheimer 2011).

We consider a region in the MC with a plane-parallel density gradient of size Zc. The density profile is expected to be a power law (e.g., Smith et a. 2009; Donkov et al. 2011). We adopt a density profile of the form (see Fig. 1): (3)with δ ~ 3/2 (e.g., Smith, Clark & Bonnelli 2009) and Zcore ~ 10-2Zc. This value of Zcore ensures that n(Zc) →nMC. We adopt n0 = 105 cm-3.

thumbnail Fig. 1

Density profile of the molecular cloud.

Open with DEXTER

4. Radiative process in bowshock–medium interactions

We consider a bowshock of a runaway star that travels through a density gradient in an MC. As the star travels through the inhomogeneous medium, the emission produced in the bowshock varies (see Fig. 2 for a sketch of the situation).

thumbnail Fig. 2

Simplified schematic of a runaway star moving through a denser region in a molecular cloud (not to scale).

Open with DEXTER

Table 1

Parameters for the stars considered in this work, along with assumptions related to relativistic particles accelerated at the reverse shock.

As mentioned before, the collision of the supersonic stellar wind with the ISM results in a system of two shocks (e.g., Wilkin 2000). Following the model developed by del Valle & Romero (2012), we assume that relativistic particles are accelerated via the first-order Fermi mechanism in the reverse adiabatic shock. This shock propagates in the opposite direction of the stellar motion, with velocity vs~VW. The stellar wind can be considered as a continuous power source, therefore both shocks, the forward and reverse shock, reach a steady state. We assume that the bowshock reaches a steady state almost immediately in its way through the density gradient, so a steady-state system can be considered for each value of z.

We perform calculations for two types of massive stars: an O4I and an O9I star, as representative examples of a very powerful and a more modest case. Their adopted parameters are listed in Table 1. Many of the parameters that define the particle energy losses and the non-thermal emission change with n. In Tables 2 and 3, the model parameters as a function of z (i.e., n) are listed for both bowshocks. The shape of the bowshock surface also changes with n, but our model is not sensitive to these changes. We consider a one-zone homogeneous cap region where particles are accelerated and emit radiation. This region is located near the apex of the bowshock, where the shock is nearly planar.

Table 2

Model parameters as a function of z for an O4I star.

Table 3

Model parameters as a function of z for an O9I star.

The acceleration timescale as a function of the energy E, for a charged particle being accelerated in a magnetic field B, is given by (e.g., Gaisser 1990; Aharonian 2004; Bosch-Ramon 2009; Romero & Paredes 2011): (4)where rL is the Larmor radius rL = E/eB, and η is a phenomenological parameter related to the efficiency of the acceleration process involved. For a non-relativistic diffusive shock acceleration, in a plane shock in the test particle approximation η can be approximated by (Drury 1983): (5)where D is the diffusion coefficient and rg = E/(eB) is the particle gyro-radius. In the Bohm limit, D = DB and DB = rgc/3 so η becomes η ~ 20/3(c/vs)2.

We estimate the magnetic field considering that the magnetic energy density is in sub-equipartition with respect to the kinetic energy , by a 0.1 factor, i.e., (6)where A is the area of a sphere of radius R0. This guarantees that the plasma remains compressible and shocks are not suppressed by the magnetic fields. The magnetic field in the acceleration region is B~4Bsub due to the compression by the shock.

4.1. Energy losses

thumbnail Fig. 3

Electron losses for z = 0, 10-2Zc, 10-1Zc, and Zc, for an O4I star. The dotted line corresponds to slow convection time.

Open with DEXTER

thumbnail Fig. 4

Electron losses for z = 0, 10-2Zc, 10-1Zc, and Zc, for an O9I star.

Open with DEXTER

We calculate the energy losses at different positions z = 0, 10-2Zc, 10-1Zc, and Zc, as the star moves through the density gradient. It is clear that this scenario is symmetric: after reaching the maximum density at z = 0, the density decreases to nMC. Particles lose energy when interacting with the magnetic fields, radiation, and matter. The electrons lose energy mainly by IC scattering, synchrotron radiation, and relativistic Bremsstrahlung. Protons cool through proton-proton in-elastic collisions with the ambient gas, but they escape from the radiation region convected away by the stellar wind. The convection time might be longer than tconv ~ Δ/VW (here Δ is the width of the shocked wind) due to turbulence in the flow driven by instabilities. Therefore, for the O4I star, we consider a regular convection time, case a, and we also consider a longer convection time, case b, which is one order of magnitude longer than tconv. For the mathematical expressions of the losses, see del Valle & Romero (2012), and references therein.

Since convection imposes the upper limit to the energy of protons for all z, the maximum energy that these particles can reach almost remains unchanged with z. For the O4I star, the highest energies protons reach are ~70 TeV and 700 TeV for case a and case b, respectively. For the O9I star, protons are accelerated up to ~1 TeV. Most of the protons escape without losing much of their energy, and they might produce non-thermal radiation further away in the cloud (del Valle et al., in prep.). The dominant loss for electrons changes with z, therefore, the maximum energy that electrons can achieve varies for each z, in general, decreasing as z decreases. In Figs. 3 and 4, we show the most relevant energy losses, at z = 0, 10-2Zc, 10-1Zc, and Zc, for the O4I and O9I stars, respectively.

4.2. Particle energy distributions

We compute the particle energy distribution for both species of particles, solving the steady state transport equation in the homogeneous approximation (Ginzburg & Syrovatskii 1964): (7)Here Q(E) is the power-law injection function, tesc is the wind convection time, and (dE/dt)loss are the radiative losses (see del Valle & Romero 2012, for further details).

The solution of Eq. (7) is a broken power law in the particle energy. The ratio of relativistic proton power to relativistic electron power, a, is unknown. We consider a = 1 (equal energy density in both particle species) for the O4I star in case a, and for the O9I star; and a = 100 (as observed in galactic cosmic rays, Ginzburg & Syrovatskii 1964), for the O4I star in case b. Figures 57 show the computed particle distributions as a function of energy per unit volume for both electrons and protons, respectively. The normalisation constant of the injection function changes with z, and because of this, the particle distributions change.

thumbnail Fig. 5

Particle distribution for electrons (top) and protons (bottom), for the O4I star in case a, at z = 0, 10-2Zc, 10-1Zc, and Zc.

Open with DEXTER

thumbnail Fig. 6

Particle distribution for electrons (top) and protons (bottom), for the O4I star in case b, at z = 0, 10-2Zc, 10-1Zc, and Zc.

Open with DEXTER

thumbnail Fig. 7

Particle distribution for electrons (top) and protons (bottom), for the O9I star, at z = 0, 10-2Zc, 10-1Zc, and Zc.

Open with DEXTER

4.3. Spectral energy distributions

We compute the non-thermal luminosity for different values of z as the star moves through the density gradient. Figures 10 and 11 show the spectral energy distributions (SEDs) at the different locations of the stars. In the case of star O4I in case a, the synchrotron and IC of IR photons are strong and dominate the SEDs for energies E < 1 TeV. The IC cut-off decreases with z, and as n increases, the p − p component gets stronger, dominating the SEDs for 1 < E < 102 TeV. For case b, all leptonic contributions are weak; for z = Zc, the synchrotron and IC of IR photons dominate the SED up to the IC cut-off; as the star moves further in the MC, the hadronic contribution becomes stronger at high energies, and dominates the SEDs for E > MeV. In the case of the O9I system, the SEDs are dominated by leptonic contributions at high energies, while the synchrotron emission is weak. In the range z = 10-2 − 0, the IC of stellar photons dominates the SEDs for E > 1 MeV because as n grows R0 decreases and the emission region gets closer to the strong stellar photon field. In all cases, the synchrotron’s lower cut-off shifts to higher energies due to the synchrotron’s self absorption; this effect becomes important as the emission region becomes more compact.

4.4. Absorption

Photons can be absorbed by different mechanisms inside the emitting region, called internal absorption and in their way to the observer, external absorption. The internal absorption is produced via photon-photon pair production; the absorbing photon fields are: the IR, non-thermal, and stellar fields (see del Valle & Romero 2012). The external absorption is produced by the matter fields and the star photon field. This latter contribution depends on the inclination angle i with the line of sight; the closest distance to the star is given by R0sini, i.e., the absorption produced by this component would be non-negligible only in the particular case when i 0 (e.g., Romero et al. 2010).

For energies <10 keV, the absorption in the ambient material is important and catastrophic for a wide range of energies. The photons must travel through the dense MC before reaching the observer. They can be absorbed by photo-ionization for Eγ< 13.6 eV, and scattered by dust for lower energies, in the IR up to the ultra-violet (UV), (e.g., Ryter 1996; Reynoso et al. 2011). We estimate this opacity as: (8)where σγN is the interaction cross section (see, e.g., Reynoso et al. 2011) and NH is the column density of the MC. We adopt here a typical value NH = 1.5 × 1022 cm-2 (e.g., Solomon et al. 1987). The denser region in the MC does not produce significant absorption in γ-rays because its maximum column density is ~ncZc < 1022 cm-2, for zc < 0.1 pc. The same situation happens with the shocked ISM: it can be very dense, but it is confined to a very thin region.

thumbnail Fig. 8

SEDs for an O41-type star at z = 0, 10-2Zc, 10-1Zc, and Zc, case a (top) and case b (bottom).

Open with DEXTER

thumbnail Fig. 9

SEDs for an O91-type star at z = 0, 10-2Zc, 10-1Zc, and Zc.

Open with DEXTER

4.5. Results

In Figs. 10 and 11, the non-thermal luminosity vs. energy curves corrected by absorption are shown. In these figures, we also indicate the 1-yr Fermi sensitivity curves for different distances. The non-thermal emission increases with n at radio and X ray energies for both stars. The gamma-ray emission also increases, except for the star O4 in case a. For the O4 star, both the radio and the gamma-ray contributions are comparable. In the case of the O9 star, the high-energy emission dominates the SED. This is because the magnetic field B is relatively low, and the IR photon field is strong (this field is the main target of the IC scattering). The emission in the energy range between the near IR to the soft X-rays is fully suppressed. Internal absorption is only non-negligible at the tail of the high-energy spectrum.

5. Discussion: variability

In the scenario presented in this paper, the bowshock itself (manifested mainly by the IR signal) might not be detected because of the obscuring MC. Soft non-thermal X-rays might also be difficult to detect due to high absorption and contamination by stronger thermal radiation. However, the radio and gamma-ray emission are not affected by the MC absorption and might be detectable. Radio emission is strong for the O4I star in case a; and gamma-ray emission dominates the energy output in both types of stars. The variation on the gamma emission at TeV energies for the O4 star, in case a and case b, is significant due to the gradual increase of the hadronic contribution to the SED.

As the star moves through the density gradient, the emission varies with a timescale ~Zc/V. For Zc ~ 10-4 pc, and the values adopted for V, the variability timescale for the O4I star is ~1 yr, and ~3 yr for the case of the O9I, with the parameters given in Table 1. Figures 12 and 13 show the integrated luminosity for the energy ranges of radio (1 GHz), X-rays (1 − 10 keV), and gamma-rays (3 × 10-2 − 100 GeV).

Nearby giant MCs are located at distances d ≥ 100 pc. The O4I system, for case a, can be a variable gamma-ray source, detectable by Fermi at every z for a wide range of distances. In case b, the detection can occur for d ≤ 1 kpc (see Fig. 10). The O9I system might be detectable at every z for d ≤ 300 pc (see Fig. 11). In this case, the system would be a variable source, but with a longer variability timescale because the star moves slower. At larger distances the source should not be detectable by Fermi, although it might appear as a weak source for the future Cherenkov Telescope Array (CTA, see Actis et al. 2011).

The gamma-ray and/or radio emission, in some cases, can be detectable only at the maximum of the light curve. Years later the emission might reappear as the stars travel through other denser regions in the MC. These sources might turn on and off within years. This situation might occur for the O4I system, in case b, when d ~ 1.4 kpc (see Fig. 10), and for the O9I system when d ~ 300 pc (see Fig. 11).

The weak or undetectable sources for Fermi might be detected by the future CTA, since it is expected to reach higher sensitivities of almost one order of magnitude better than Fermi at E ~ 100 GeV. In particular, CTA might be able to detect runaway O9I stars, which are more common in the solar neighbourhood.

If the bowshock of a well-identified massive runaway star moving through an MC could be detected as a gamma-ray source, variations in its emission will allow the study of the fine structure of the MC, which is not possible at other energy bands where the emission is highly absorbed.

thumbnail Fig. 10

Total luminosity curves at z = 0, 10-2Zc, 10-1Zc, and Zc, for an O4I star, case a (top) and case b (bottom). Fermi sensitivities at 0.8 and 1.4 kpc are also shown.

Open with DEXTER

thumbnail Fig. 11

Total luminosity curves at z = 0, 10-2Zc, 10-1Zc, and Zc, for an O9I star. Fermi sensitivities at 100 pc and 300 pc are also shown.

Open with DEXTER

thumbnail Fig. 12

Variability curves at z = 0, 10-2Zc, 10-1Zc, and Zc, for an O4I star, case a (top) and case b (bottom).

Open with DEXTER

thumbnail Fig. 13

Variability curves at z = 0, 10-2Zc, 10-1Zc, and Zc, for an O9I star.

Open with DEXTER

The existence of a population of galactic variable gamma-ray sources is suspected since the epoch of the Energetic Gamma Ray Telescope (EGRET). Moreover, a statistically very significant positional correlation was found between gamma-ray sources in the third EGRET catalogue for OB star associations (e.g., Romero et al. 1999; Torres et al. 2001).

In the second Fermi catalogue (Nolan et al. 2012), 352 sources previously listed in the first catalogue (Abdo et al. 2010) did not show up. Most of these sources are concentrated along the Galactic plane. Although some sources may have disappeared due to improvements in the model of the diffuse background, which is most intense at low latitudes, some other sources might intrinsically vary their fluxes from one catalogue to another. This could be the case of massive runaway stars moving through molecular clouds.

We conclude that under some assumptions, bowshocks of massive runaway stars travelling through an MC might be variable gamma-ray sources; their time variability scale depends on the size scales of density inhomogeneities and the stellar velocities. Under some conditions, these sources might be galactic gamma-ray sources turning on and off over years. In this work, then, we propose a putative new class of galactic variable gamma-ray source.


1

A condensate structure of higher density than the average (>3σ) density in the MC.

Acknowledgments

We thank the anonymous referee and Dr M. Reynoso for valuable comments. This work is supported by PIP 0078 (CONICET) and PICT 2007-00848/2012-00878, Préstamo BID (ANPCyT).

References

All Tables

Table 1

Parameters for the stars considered in this work, along with assumptions related to relativistic particles accelerated at the reverse shock.

Table 2

Model parameters as a function of z for an O4I star.

Table 3

Model parameters as a function of z for an O9I star.

All Figures

thumbnail Fig. 1

Density profile of the molecular cloud.

Open with DEXTER
In the text
thumbnail Fig. 2

Simplified schematic of a runaway star moving through a denser region in a molecular cloud (not to scale).

Open with DEXTER
In the text
thumbnail Fig. 3

Electron losses for z = 0, 10-2Zc, 10-1Zc, and Zc, for an O4I star. The dotted line corresponds to slow convection time.

Open with DEXTER
In the text
thumbnail Fig. 4

Electron losses for z = 0, 10-2Zc, 10-1Zc, and Zc, for an O9I star.

Open with DEXTER
In the text
thumbnail Fig. 5

Particle distribution for electrons (top) and protons (bottom), for the O4I star in case a, at z = 0, 10-2Zc, 10-1Zc, and Zc.

Open with DEXTER
In the text
thumbnail Fig. 6

Particle distribution for electrons (top) and protons (bottom), for the O4I star in case b, at z = 0, 10-2Zc, 10-1Zc, and Zc.

Open with DEXTER
In the text
thumbnail Fig. 7

Particle distribution for electrons (top) and protons (bottom), for the O9I star, at z = 0, 10-2Zc, 10-1Zc, and Zc.

Open with DEXTER
In the text
thumbnail Fig. 8

SEDs for an O41-type star at z = 0, 10-2Zc, 10-1Zc, and Zc, case a (top) and case b (bottom).

Open with DEXTER
In the text
thumbnail Fig. 9

SEDs for an O91-type star at z = 0, 10-2Zc, 10-1Zc, and Zc.

Open with DEXTER
In the text
thumbnail Fig. 10

Total luminosity curves at z = 0, 10-2Zc, 10-1Zc, and Zc, for an O4I star, case a (top) and case b (bottom). Fermi sensitivities at 0.8 and 1.4 kpc are also shown.

Open with DEXTER
In the text
thumbnail Fig. 11

Total luminosity curves at z = 0, 10-2Zc, 10-1Zc, and Zc, for an O9I star. Fermi sensitivities at 100 pc and 300 pc are also shown.

Open with DEXTER
In the text
thumbnail Fig. 12

Variability curves at z = 0, 10-2Zc, 10-1Zc, and Zc, for an O4I star, case a (top) and case b (bottom).

Open with DEXTER
In the text
thumbnail Fig. 13

Variability curves at z = 0, 10-2Zc, 10-1Zc, and Zc, for an O9I star.

Open with DEXTER
In the text

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

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

Initial download of the metrics may take a while.