Gamma rays from red giant wind bubbles entering the jets of elliptical host blazars

TL;DR

This paper investigates how red giant wind bubbles entering blazar jets in elliptical galaxies produce gamma-ray emission, analyzing their dynamics and radiation to assess detectability and duty cycles.

Contribution

It provides a combined analytical and numerical study of red giant wind bubble interactions with blazar jets, estimating gamma-ray emission and detection prospects.

Findings

Analytical model effectively describes bubble-jet interactions.

Gamma-ray emission could be detectable within the local universe.

Duty cycle of gamma-ray production can be up to 1.

Abstract

Blazars in elliptical hosts have a population of red giants surrounding their jet. These stars can carry large wind-blown bubbles into the jets, leading to gamma-ray emission through bubble-jet interactions. We study the interaction dynamics and the gamma-ray emission produced when the bubbles formed by red giant winds penetrate the jet of a blazar in an elliptical galaxy. First, we characterized the masses and penetration rates of the red giant wind bubbles that enter the jet. Then, the dynamical evolution of these bubbles under the jet impact was analysed analytically and numerically, and the radiation losses of the particles accelerated in the interaction were characterised. Finally, the synchrotron and the inverse Compton contributions above 100 MeV were estimated under different jet magnetic fields, powers, and Lorentz factors. We find that an analytical dynamical model is a reasonable approximation for the red giant wind bubble-jet interaction. The radiation produced by these wind bubbles interacting with a jet can have a duty cycle of up to 1. For realistic magnetic fields, gamma rays could be detectable from sources within the local universe, preferentially from those with high Lorentz factors (of about 10), and this could be a relatively common phenomenon for these sources. For magnetic fields in equipartition with the jet power, and high acceleration rates, synchrotron gamma rays may be detectable even for modest Lorentz factors (of about 3), but with a much lower duty cycle.

Tables

Table 1. Table 1: Results.

$L_j$ [erg s-1]	${10}^{43}$	${10}^{44}$	${10}^{45}$	${10}^{43}$	${10}^{44}$	${10}^{45}$
$\Gamma$	3	3	3	10	10	10
Inverse Compton: $E_e^{\prime }=E_{\mathrm{IC}}^{\prime }$ $B=0.1B_{\mathrm{eq}}$
$L_{\mathrm{st}}$ [erg s-1]	$6.5\times {10}^{39}$	$1.7\times {10}^{40}$	$1.7\times {10}^{40}$	$2.5\times {10}^{39}$	$6.9\times {10}^{39}$	$1.8\times {10}^{40}$
$L_{\mathrm{pop}}$ [erg s-1]	$1.1\times {10}^{38}$	$6.5\times {10}^{38}$	$8.4\times {10}^{38}$	$3.7\times {10}^{40}$	$4.7\times {10}^{40}$	$3.4\times {10}^{40}$
$L_{\mathrm{dom}}$ [erg s-1]	$1.4\times {10}^{38}$	$8.8\times {10}^{37}$	$5.5\times {10}^{37}$	$3.5\times {10}^{40}$	$3.1\times {10}^{40}$	$2.5\times {10}^{40}$
[$M_{\odot }$]	$1.6\times {10}^{-6}$	$1.7\times {10}^{-6}$	$1.5\times {10}^{-6}$	$1.7\times {10}^{-6}$	$1.7\times {10}^{-6}$	$1.5\times {10}^{-6}$
[kpc]	0.33	1.0	1.6	0.67	1.1	1.6
$PR$ [yr-1]	$8.3\times {10}^{-3}$	$2.6\times {10}^{-2}$	$4.0\times {10}^{-2}$	$8.6\times {10}^{-3}$	$2.7\times {10}^{-2}$	$4.0\times {10}^{-2}$
$t_{\mathrm{obs}}$ [yr]	$54$	$170$	$270$	$40$	$33$	$24$
Synchrotron: $E_e^{\prime }=E_{\mathrm{Sy}}^{\prime }$, $B=B_{\mathrm{eq}}$
$L_{\mathrm{st}}$ [erg s-1]	$7.5\times {10}^{41}$	$1.9\times {10}^{42}$	$1.9\times {10}^{42}$	$1.3\times {10}^{41}$	$7.4\times {10}^{41}$	$2.5\times {10}^{42}$
$L_{\mathrm{pop}}$ [erg s-1]	$1.4\times {10}^{40}$	$4.0\times {10}^{40}$	$6.7\times {10}^{40}$	$1.1\times {10}^{42}$	$2.0\times {10}^{42}$	$2.5\times {10}^{42}$
$L_{\mathrm{peak}}$ [erg s-1]	$5.8\times {10}^{41}$	$5.8\times {10}^{42}$	$5.8\times {10}^{43}$	$6.3\times {10}^{42}$	$6.3\times {10}^{43}$	$6.4\times {10}^{44}$
[$M_{\odot }$]	$1.5\times {10}^{-6}$	$1.7\times {10}^{-6}$	$1.6\times {10}^{-6}$	$1.7\times {10}^{-6}$	$1.7\times {10}^{-6}$	$1.6\times {10}^{-6}$
[kpc]	0.26	0.88	1.9	0.91	1.3	1.9
$PR$ [yr-1]	$5.7\times {10}^{-3}$	$2.0\times {10}^{-2}$	$3.8\times {10}^{-2}$	$7.4\times {10}^{-3}$	$2.3\times {10}^{-2}$	$3.4\times {10}^{-2}$
$t_{\mathrm{peak}}$ [yr]	1.3	0.14	0.014	4.8	0.48	0.046
Synchrotron: $E_e^{\prime }=E_{\mathrm{Sy}}^{\prime }$, $B=B_{\min }$
$B_{\min }$ [$B_{\mathrm{eq}}$]	$0.02$	$0.01$	$0.01$	$0.04$	$0.03$	$0.02$
$L_{\mathrm{pop}}$ [erg s-1]	$2.1\times {10}^{40}$	$6.0\times {10}^{40}$	$1.4\times {10}^{41}$	$1.7\times {10}^{42}$	$1.8\times {10}^{42}$	$1.4\times {10}^{42}$
$L_{\mathrm{dom}}$ [erg s-1]	$6.9\times {10}^{40}$	$2.1\times {10}^{40}$	$4.3\times {10}^{40}$	$1.8\times {10}^{42}$	$1.1\times {10}^{42}$	$1.9\times {10}^{42}$
$t_{\mathrm{obs}}$ [yr]	$13$	$45$	$32$	$38$	$65$	$10$

Equations

$$P_{\mathrm{w}}(\dot{M})={\rho }_{\mathrm{w}}(\dot{M})v_{\mathrm{w}}^2=P_{\mathrm{ext}}$$

$$P_{\mathrm{orb}}(z)={\rho }_{\mathrm{ISM}}{v}_{\mathrm{orb}}(z{)}^2,$$

$$P_{\mathrm{lat}}(z)=L_{\mathrm{j}}/c\pi z^2{\Gamma }_{\mathrm{j}}^2,$$

$$P_{\mathrm{ext}}=\max (P_{\mathrm{lat}},P_{\mathrm{ISM}}+P_{\mathrm{orb}}).$$

$$R_{\mathrm{out}}^2=\frac{\dot{M}{v}_{\mathrm{w}}}{4\pi p_{\mathrm{ext}}},M_{\mathrm{out}}=4\pi R_{\mathrm{out}}^3{\rho }_{\mathrm{w}}(R_{\mathrm{out}}).$$

$$c_s(R_{\mathrm{in}})<v_{\mathrm{orb}},R_{\mathrm{in}}^2=\frac{v_{\mathrm{orb}}^2\dot{M}}{4\pi v_{\mathrm{w}}{P}_{\mathrm{j}}},$$

$$R_{\mathrm{j}}(z)=\{\begin{array}{ll}\theta z,& \text{mbox}ifz<z_{\mathrm{eq}}\\ Const,& \text{mbox}ifz\ge z_{\mathrm{eq}}\end{array},$$

$$t_{\mathrm{acc}}\simeq \{\begin{array}{ll}{z}_0/\beta c,& \text{mbox}ifD<1\\ z_0/D\beta c,& \text{mbox}ifD>1\end{array},D\equiv \frac{P_{\mathrm{j},0}\pi R_{\mathrm{b}}^2{z}_0}{4c^2{M}_{\mathrm{b}}{\Gamma }_{\mathrm{j}}^3},$$

$$\frac{B_0^2}{4\pi }=\frac{L_{\mathrm{j}}}{\pi R_{\mathrm{j}}(z_0{)}^{2}c},$$

$$B(z)=B_0{\left(\frac{z_0}{z}\right)}^2,$$

$$<A>=\frac{\int A(z,t_{\mathrm{RG}})L_{\mathrm{pop}}(z,t_{\mathrm{RG}})\mathrm{d}z\mathrm{d}{t}_{\mathrm{RG}}}{L_{\mathrm{pop}}(z,t_{\mathrm{RG}})\mathrm{d}z\mathrm{d}{t}_{\mathrm{RG}}},$$

$$\frac{\mathrm{d}{L}_{\mathrm{st}}(z)}{\mathrm{d}V}=\eta L_{\mathrm{j}}{F}_{\mathrm{rad}}(z)\frac{{\delta }_{\mathrm{j}}^3}{{\Gamma }_{\mathrm{j}}^3}\int ⟨\frac{S_{\mathrm{s}}(z,t)}{S_{\mathrm{j}}(z)}⟩n_{\mathrm{s}}(z,t)\mathrm{d}t,$$

$$E_{\mathrm{b}}^{\prime }(M_{\mathrm{b}},z_0)={\int }_{z_0}^{z_{\max }}\eta \frac{M_{\mathrm{b}}{c}^2}{t_{\mathrm{IC}}^{\prime }(z)}\frac{\mathrm{d}z}{\Gamma \beta c};$$

$$L_{\mathrm{pop}}={\delta }_{\mathrm{j}}^3{\int }_{z_{\mathrm{BH}}}^H{E}_{\mathrm{b}}^{\prime }(m,z)PR(m,z)\mathrm{d}m\mathrm{d}z,$$

$$E_{\mathrm{b}}(M_{\mathrm{b}},z)=\eta F_{\mathrm{rad}}\overline{F_{\mathrm{e}}}{M}_{\mathrm{b}}{c}^2{\delta }_j^3,$$

$$L_{\mathrm{peak}}(M_{\mathrm{b}},z)=\eta F_{\mathrm{rad}}c{F}_{\mathrm{e},\max }{P}_0{\delta }_{\mathrm{j}}^2\pi r_{\mathrm{b}}^2,$$

Full text

*Aims. *We study the interaction dynamics and the gamma-ray emission produced when the bubbles formed by red giant winds penetrate the jet of a blazar in an elliptical galaxy.

*Methods. *First, we characterized the masses and penetration rates of the red giant wind bubbles that enter the jet. Then, the dynamical evolution of these bubbles under the jet impact was analysed analytically and numerically, and the radiation losses of the particles accelerated in the interaction were characterised. Finally, the synchrotron and the inverse Compton contributions above $\sim 100$ MeV were estimated under different jet magnetic fields, powers, and Lorentz factors.

*Results. *We find that an analytical dynamical model is a reasonable approximation for the red giant wind bubble-jet interaction. The radiation produced by these wind bubbles interacting with a jet can have a duty cycle of up to $\sim 1$. For realistic magnetic fields, gamma rays could be detectable from sources within the local universe, preferentially from those with high Lorentz factors ($\sim 10$), and this could be a relatively common phenomenon for these sources. For magnetic fields in equipartition with the jet power, and high acceleration rates, synchrotron gamma rays may be detectable even for modest Lorentz factors ($\sim 3$), but with a much lower duty cycle.

*Conclusions. *Blazars in elliptical galaxies within the local universe can produce detectable transient or persistent gamma-ray emission from red giant wind bubbles entering their jets.

††offprints: N. Torres-Albà

1 Introduction

Supermassive black holes, present in the innermost regions of galaxies, may accrete the material surrounding them, becoming active galactic nuclei (AGN). Some AGN produce collimated relativistic outflows, or jets (e.g. Begelman et al., 1984), which propagate through the host galaxy. This propagation will inevitably lead to the jet interacting with a variety of obstacles including stars, gas, and dense clouds. These interactions may affect the jet dynamically (e.g. Blandford & Koenigl, 1979; Wang et al., 2000; Sutherland & Bicknell, 2007). In particular, stars with high mass-loss rates may load the jet with enough matter to result in deceleration (e.g. Komissarov, 1994; Bowman et al., 1996; Hubbard & Blackman, 2006; Bosch-Ramon et al., 2012; Perucho et al., 2014, 2017).

The winds of stars interacting with AGN jets produce a double bow-shock structure in which particles can be accelerated to relativistic energies, possibly contributing to the jet’s total non-thermal emission. Several works have explored this interaction and the resulting emission, both in the case of steady radiation and transient events (e.g. Bednarek & Protheroe, 1997; Barkov et al., 2010, 2012; Khangulyan et al., 2013; Araudo et al., 2013; Bednarek & Banasiński, 2015; Wykes et al., 2015; Bosch-Ramon, 2015; de la Cita et al., 2016; Vieyro et al., 2017), with results that often point towards possible detectability for nearby sources.

Previous works have considered persistent emission being generated from a whole population of stars, though in all cases they focus on the interaction that occurs within the jet once the star has already penetrated (Araudo et al., 2013; Wykes et al., 2015; Bosch-Ramon, 2015; Vieyro et al., 2017). We refer to this stage of interaction as “steady state”.

In this work, we focus on the possible mass-loading and emission generated at the moment when stars penetrate the jet, and the latter interacts with large “bubbles” of material formed by the collision between the stellar wind and the interstellar medium. Perucho et al. (2017) inferred possible significant non-thermal emission and mass-loading from this early stage in jet-star interaction, using 2D and 3D simulations of one single star with heavy mass loss. We use semi-analytical prescriptions to estimate if this is the case for a whole population of stars within a galaxy.

We focus here on the study of blazar sources, as Doppler boosting is an important factor in enhancing the resulting emission. We consider only low-luminosity sources (i.e. $L_{\rm j}=10^{43}-10^{45}$ erg s*-1*), more abundant in the local universe. As recent studies show that the preferred hosts of blazars are late-type galaxies (e.g. Scarpa et al., 2000; Nilsson et al., 2003; Falomo et al., 2014; Olguín-Iglesias et al., 2016), we model the red giant population within an elliptical bulge. A mostly phenomenological approach is adopted (with the exception of an illustrative numerical simulation), based on specific source knowledge and reference parameter values, as a first simplified step to explore the outcome of the wind bubble and jet interaction.

The paper is organized as follows. A description of the prescriptions used to characterize the stellar population is given in Sect. 2. In Sect. 3, the properties of the stellar bubbles outside and immediately after penetrating the jet are described. The bubble evolution within the jet is described through analytical estimates and compared with simulation results in Sect. 4. Then, the non-thermal emission generated by bubble-jet interactions is estimated in Sect. 5. Finally, the discussion is presented in Sect. 6.

2 Characterization of the stellar population in an elliptical galaxy

Elliptical galaxies contain large populations of red giants, which can have high mass-loss rates, in the range of $\dot{M}\sim 10^{-10}-10^{-5}\,M_{\odot}$ yr*-1*(Reimers, 1975). We model the red giant population of any elliptical galaxy by taking as reference values those of M87, as its proximity allows for a precise study.

2.1 Stellar number density, mass-loss rate, and wind speed

Assuming a Kroupa initial mass function (IMF), and normalizing it to the total mass of stars ($M_{\rm T}$), we can estimate the number of stars within the bulge. We take as index for the IMF $x_{1}=-1.3$ for $0.1$ M${}_{\odot}<M<0.5$ M⊙ and $x_{2}=-2.3$ for $0.5$ M${}_{\odot}<M<m_{2}$ (Kroupa, 2001), where $m_{2}$ is the mass of the stars exiting the red giant phase (i.e. the most massive stars present; see below). We estimate the total stellar mass adopting that contained within the galactic bulge in M87 (Gebhardt & Thomas, 2009), knowing that the radius is $\sim 40^{\prime\prime}$ (Harris et al., 1999), corresponding to a bulge radius of $R_{\rm b}\sim 3.1$ kpc, which we consider spherical.

From this point on, we follow the calculations in Vieyro et al. (2017) to derive the mass of the red giants in the bulge of the galaxy assuming that all stars formed at the birth of the galaxy (i.e. no star formation extended in time), which yields $\sim 0.83\,M_{\odot}$, and their number, which is $N_{\rm T}\sim 1.3\times 10^{9}$. From Gebhardt & Thomas (2009), we derive that the decay of the density with radial distance/jet height ($z$) can be approximated by $n_{\rm s}(z)\propto N_{\rm T}/z$ in the considered inner $\sim 40^{\prime\prime}$.

We derive the mass-loss rate and wind speed of the red giant population exactly as in Vieyro et al. (2017) for the particular case of M87. The result is a mass-loss rate that increases rapidly with the age of the red giant. Stars are thus modelled as a distribution that depends both on height and red-giant age (i.e. how deep into the red giant phase the star is), $n_{\rm s}(z,t_{\rm RG})$. For simplicity, we consider a stellar wind of $v_{\rm w}=10^{7}$ cm s*-1*, and consider it constant during its evolution.

2.2 Orbital velocities and penetration rate

In order to study the collective emission and mass-loading generated by the whole population of stars as they penetrate the jet, we need to determine the frequency at which these events occur, that is, the penetration rate ($PR$). Knowing the orbital velocities of stars, one can estimate the penetration rate into the jet for the distance interval $(z,z+{\rm d}z)$ as ${\rm d}P(z,t_{\rm RG})\simeq n_{\rm s}(z,t_{\rm RG})v_{\rm orb}(z)R_{\rm j}(z){\rm d}z$ (Khangulyan et al., 2013), where $R_{\rm j}$ is the jet radius.

At low $z$ the stellar orbital movement is dominated by the central supermassive black hole, and thus stars orbit it following a Keplerian motion, with $v_{\rm orb}=\sqrt{GM_{\rm BH}/z}$. The gravitational influence of the black hole is dominant within a radius, or jet height, $z_{\rm g}={GM_{\rm BH}}/{\sigma^{2}}$, where $\sigma$ is the stellar velocity dispersion of the bulge. At larger $z$, we consider the stars to move within the bulge at a constant velocity $\sigma$.

For the black-hole mass, we use the value derived by Gebhardt et al., 2011, $6.6\pm 0.4\times 10^{9}$ M*⊙. Dispersion velocity measurements decrease from $\sigma\sim 480$ km s-1* near the nucleus (where the supermassive black hole is dominant) to $\sigma\sim 320$ km s*-1* at $\sim 40^{\prime\prime}$ (Gebhardt et al., 2011). We take the average value of $\sigma(z)$, $\sigma=360$ km s*-1*, as the constant velocity for stars within the bulge, which gives $z_{\rm g}=220$ pc ($\sim 2.5^{\prime\prime}$).

3 Interaction with the jet

3.1 Stars outside the jet

As a star moves outside the jet, the ram pressure generated by its stellar wind is in equilibrium with all external pressures, meaning

[TABLE]

The external pressures are given by the interstellar medium (ISM) thermal pressure $P_{\rm ISM}\approx 10^{-12}$ erg cm*-2*, and the orbital motion of the star through this medium, or “orbital (ram) pressure”, is given by

[TABLE]

where $\rho_{\rm ISM}$ is the ISM density. We fix the ISM density taking one hydrogen atom per cm3, a typical value in the central regions of elliptical galaxies (Tan et al., 2008), throughout the entire bulge.

As the star approaches the jet, the jet lateral pressure $P_{\rm lat}$ may be larger than the pressures generated by the movement within the ISM:

[TABLE]

where $L_{\rm j}$ and $\Gamma_{\rm j}$ are the jet luminosity and Lorentz factor, respectively (Bosch-Ramon & Barkov, 2016). Therefore, in close proximity with the jet,

[TABLE]

We call $R_{\rm out}$ the distance from a star at which pressure equilibrium is reached outside the jet. At this distance, the colliding pressures generate a double bow shock in which both stellar-wind material and interstellar material are accumulated. This shocked layer surrounds the star in a bullet-like shape, with material gathering in the direction of the stellar movement toward the jet, and eventually escaping in the opposite direction.

We approximate this shocked region as a sphere of radius $R_{\rm out}$, in which an amount of material of mass $M_{\rm out}$ is contained, which can be estimated as

[TABLE]

We do not consider the mass accumulated within the shocked layer, which would be of the order of $\left(v_{\rm w}/v_{\rm orb}\right)^{2}M_{\rm out}\lesssim 0.1\,M_{\rm out}$.

The dominant pressure component for most of the considered jet range is $P_{\rm orb}$, with $P_{\rm lat}$ being of comparable value at the base of the jet.

3.2 Stars penetrating the jet

As a star begins to penetrate the jet, the most external layers of the bubble it carries are hit by the jet ram pressure. This results in a shock that starts to propagate through those layers with a speed $c_{s}\approx\sqrt{{L_{j}}/{S_{j}c\rho_{w}(R)}}$. If the speed at which the shock propagates is sufficiently slow (i.e. the bubble is dense enough), part of the bubble material will penetrate the jet along with the star before the shock reaches the stagnation radius, $R_{\rm s}$. This will happen for all layers for which

[TABLE]

where $P_{\rm j}$ is the jet pressure. This means layers with $R>R_{\rm in}$ will be expelled before the star fully penetrates the jet, and lost at the jet contact discontinuity (CD), while all layers within $R_{\rm in}$, with mass $M_{\rm in}=4\pi R_{\rm in}^{3}\rho_{\rm w}(R_{\rm in})$, will manage to penetrate.

Once inside the jet, the shock will continue to propagate in the wind until it reaches the stagnation radius $R_{\rm s}$, the distance from the star where the wind and jet ram pressures are equal. There, a double bow shock is formed in which particles can be accelerated up to relativistic energies. In this work, we refer to this emission as steady state emission (e.g. Vieyro et al., 2017), as the jet-wind interaction process is continuous while the stars are inside the jet.111Instabilities produced at the jet-wind interaction region may actually lead to individual star-jet interaction variability (de la Cita et al., 2016).

Stellar material contained in the range $R_{\rm s}-R_{\rm in}$ is expelled within the jet as a blob (Bosch-Ramon et al., 2012; Perucho et al., 2017). The typical initial sizes of the bubbles, compared with the jet radius, are shown in Fig. 1 for a red giant of average mass-loss rate ($\dot{M}=5.7\times 10^{-10}$ M*⊙* yr*-1*), and for the red giant mass-loss rate that generates the dominant events ($\dot{M}_{\rm dom}=1.4\times 10^{-8}$ M*⊙* yr*-1*, as calculated in Sect. 5.2).

For the particular set of parameters used to plot Fig. 1, above $z\sim 400$ pc the size of the bubble inside the jet is shown as larger than the size of the bubble outside the jet. In such a case, the material introduced into the jet would be that contained within $R_{\rm s}-R_{\rm out}$. Another possibility, for a different set of parameters, is that $R_{\rm s}$ is larger than $R_{\rm in}$, but smaller than $R_{\rm out}$. In such a case, the star would lose any outer layers in the CD, and once inside the jet, the stellar wind termination region would expand up to the stagnation radius.

Figure 1 shows the mass contained in the bubble outside and right after penetrating the jet, for the average and the event-dominant mass-loss rates. The jet is loaded with the external matter brought by the bubbles expelled by stars at penetration, though this mass-load rate is much lower than the jet $\dot{M}_{\rm jet}=L_{\rm j}/\Gamma_{\rm j}c^{2}$, and therefore unlikely to result in a dynamical effect on the jet.

As has been studied for example by Wykes et al. (2015) and Vieyro et al. (2017), a population of high-mass stars in starburst galaxies can interact strongly with the jet. Young OB stars have stronger winds than red giants, with higher mass-loss rates and speeds. However, these stars have such fast winds that $v_{\rm w}>v_{\rm orb}$. Following Eqs. 5, when that condition is met, $R_{\rm s}>R_{\rm in}$. Therefore, the star penetrates the jet and the size of the interaction region surrounding the star actually increases, with no significant external material introduced into the jet.

4 Bubble evolution within the jet

We consider a jet that initially expands with a conical geometry, launched close to the supermassive black hole in the centre of the galaxy. We consider that the jet recollimates, which we model as the jet becoming cylindrical:

[TABLE]

when its pressure becomes equal to that of the ISM.

After the star penetrates the jet and a bubble of stellar material is expelled, this bubble evolves as a result of the interaction with the jet. We have adopted the analytical modelling of the evolution of a blob impacted by a jet developed by Barkov et al. (2010, 2012) and Khangulyan et al. (2013) (for previous numerical simulations of this process, see e.g. Bosch-Ramon et al., 2012; Perucho et al., 2017).

The shock produced by the impact of the jet ram pressure causes the material of the bubble to heat up, quickly expand, and accelerate, resulting in the bubble reaching relativistic speeds. This acceleration occurs on timescales of

[TABLE]

where $D$ is a dimensionless parameter related to the jet luminosity and Lorentz factor and the mass of the bubble, $R_{\rm b}$ is the radius of the bubble after acceleration (see e.g. Barkov et al., 2012), and any subindex “0” refers to the value at $z_{\rm 0}$, the height at which the star penetrates. The value $D$ actually gives a comparison between bubble and jet mass on scales of $z_{0}$; that is, if $D>1$ ($<1$), the jet will (not) effectively accelerate the blob before it covers a distance $\sim z_{0}$. After the blob is accelerated, it is carried downstream (likely at least partially disrupted) until jet termination, at a height $H$ and at a speed $\sim\beta c$, thus on a timescale $t_{\rm esc}\sim(H-z_{\rm 0})/\beta c$.

This analytical approximation to the bubble evolution within the jet, which we use to calculate the results presented in Sect. 5, can be compared to results obtained through numerical simulations. For that, we simulate here the bubble-jet interaction solving numerically the equations of relativistic hydrodynamics (RHD) assuming axisymmetry, a gas with constant adiabatic index $\bar{\gamma}=$4/3, and a dynamically negligible magnetic field. The RHD equations are solved using the Marquina flux formula (Donat & Marquina, 1996; Donat et al., 1998); further details on the code can be found in de la Cita et al. (2016).

The resolution of the calculations is 1000 cells in the vertical direction, the z-axis, and 200 in the radial direction, the r-axis. Those numbers of cells correspond to a physical range of $z_{\rm grid}=1000-1300$ pc and $r_{\rm grid}=0-65$ pc. The number of cells was chosen such that the results did not change significantly when going to higher resolution. The boundary conditions were set to inflow at $z_{\rm grid}=0$ with parallel stream lines (the jet is recollimated), reflection at $r_{\rm grid}=0$, and outflow in the remaining boundaries.

Figure 2 shows three density maps, illustrative snapshots of a simulation of the jet-bubble interaction at $\sim$0.4, 400, and 950 yr after bubble penetration inside the jet. Initial properties are: $L_{\rm j}=10^{44}$ erg $s^{-1}$, $\Gamma_{\rm j}=10$, $M_{\rm b}=2\times 10^{-6}$ M*⊙*, and $r_{\rm b,0}=1.5$ pc. The simulation starts once the cloud has already expanded significantly, moving with a Lorentz factor of 2, in the relativistic regime (e.g. Barkov et al., 2012; Khangulyan et al., 2013). We focus only on the relativistic regime, as starting from its actual initial radius would require $10^{2}-100^{2}$ times more cells and thus huge computing costs. As seen in the figure, the shocked bubble evolves upstream of the jet; most of the bubble expansion occurs at early times, and despite partial disruption the structure keeps its integrity.

In Fig. 2 we also plot the bubble radius and Lorentz factor as a function of jet height obtained from the simulation, and compare it to the results of the analytical calculation. The simulated evolution is slightly slower than analytically predicted. We note, however, that while the Lorentz factor of the bubble in the simulation reaches a value of $\sim 8$, it should eventually reach the jet Lorentz factor at higher $z$. The chosen duration of the simulation of $\sim 1000$ yr was adopted as a trade-off between a moderated computational time and illustrative effectiveness.

A noticeable difference in Fig. 2 between simulated and analytical results is bubble radius, computed as the mass-averaged cylindrical radius to facilitate comparison. The radius evolves more slowly in the simulation, reaching a final value of a factor $\lesssim 2$ smaller than the analytical estimate. We expect nevertheless more convergence at larger $z$ values.

We show thus that the analytical first-order approximation to the bubble evolution used in this work, adopted to derive the results presented in Table 1, is reasonably accurate in the relativistic regime. While instabilities might be important, once the shocked bubble reaches relativistic speeds, its evolution is similar to the analytic predictions, even if disrupted as a filamentary fragmented structure mixed with shocked jet material. The reason is that the transversal expansion of the shocked bubble material is slowed down in the laboratory frame by flow-frame time dilation, and in the longitudinal direction by very similar bottom and the top speeds of the shocked structure in the laboratory frame.

Therefore, from the simulation results, we find it reasonable at this stage to adopt a final constant bubble radius, and a final bubble speed close to that of the jet if $D>1$. For the bubbles that dominate the non-thermal emission (see Sect. 5.2), this condition is fulfilled.

A more complete simulation, including earlier and later stages of the bubble evolution, is planned for future work.

5 Radiated non-thermal power

Particles are assumed to be efficiently accelerated at the interaction between the jet and the bubble. The non-thermal energy contained in these particles is uncertain, although we focus here on the case in which the energy budget is significant and detectable radiation may be produced. The specific acceleration mechanism is not considered, and it is just assumed that some fraction of the shocked jet luminosity goes to non-thermal particles.

5.1 Energy losses

We assume that electrons (and positrons) are the dominant non-thermal emitting particles. Accelerated protons may also be present in the jet, although on the scales of interest significant radiation losses are not expected for these particles (see however, e.g. Aharonian, 2000). For electrons, we considered the radiative losses via inverse Compton (IC) and synchrotron emission. Regarding non-radiative losses, we considered adiabatic losses for the conical jet, and none after recollimation (but escape when the bubble reaches $H$).

We computed the IC losses as in Bosch-Ramon & Khangulyan (2009) (see also Khangulyan et al., 2014). Their approximation is valid for a Planck distribution of target photons of temperature $T$, and must be renormalized to the energy density of the considered target photon field. We consider the photon field in an elliptical galaxy as derived by Vieyro et al. (2017).

We also considered synchrotron losses (e.g. Longair, 1981) when taking $B=B_{\rm eq}$, a magnetic field of equivalent energy density to that of the jet, meaning

[TABLE]

with the field depending on height as

[TABLE]

where $B_{0}=B(z_{0})$. We discuss the effects on the synchrotron radiation of considering a lower magnetic field in Sect. 5.2.3.

5.2 Radiation

For all combinations of $L_{\rm j}=10^{43},\ 10^{44},\ 10^{45}$, and $\Gamma_{\rm j}=3,\ 10$, we computed the average luminosity of the whole population of bubbles penetrating the jet (as seen by the observer), $L_{\rm pop}$. Results are presented in Table 1. All parameters except for jet luminosity, $L_{\rm j}$, and Lorentz factor, $\Gamma_{\rm j}$, are fixed as described in Sects. 2 and 3.

We evaluated the emitted non-thermal radiation at two characteristic electron energies: first, for IC emission, at the Thompson-Klein Nishina transition energy, at $E^{\prime}_{\rm IC}=(m_{\rm e}c^{2})^{2}/kT\Gamma$, where the maximum of IC emission is expected for reasonable electron energy distributions, falling in the gigaelectronvolt-Teraelectronvolt range; secondly, for synchrotron emission, for the electrons that generate synchrotron 100 MeV photons as seen by the observer, at $E^{\prime}_{\rm sy}$.

We assumed that at most a fraction $\eta$ of the energy acquired by jet acceleration can be radiated222The bubble acquires in the laboratory frame an energy $\sim\Gamma_{\rm j}M_{\rm b}c^{2}$ due to jet acceleration., fixed to 0.1 throughout this work. From first principles, it is not possible to derive the value of $\eta$, but the adopted choice maximizes the predicted gamma-ray emission without assuming a full non-thermal conversion of the available energy. All the predicted luminosities thus scale with $\eta/0.1\leq 10$.

We estimated the typical properties of the bubbles that contribute the most to the overall luminosity as

[TABLE]

where $A$ denotes the quantity we are interested in evaluating (i.e. $M_{\rm b}$ and $z_{\rm 0}$). We refer to events with these properties as the “dominant” events, and list their characteristics in Table 1. Their penetration rate (i.e. rate at which events occur) along with peak luminosity of the event, $L_{\rm dom}$ or $L_{\rm peak}$, and typical duration of the peak, $t_{\rm obs}$ or $t_{\rm peak}$, are also listed. We note that $L_{\rm dom}$ is the average luminosity of a dominant event along its evolution within the jet; in reality, the luminosity is a function of jet height (see Eq. 13), larger at $z_{0}$, where the star penetrates, and becomes progressively dimmer as the radiative cooling efficiency diminishes.

We also compared the emission generated through this interaction mechanism to that generated in steady state, $L_{\rm st}$, by the same population of red giants (see Sect. 3). The apparent non-thermal luminosity per unit volume at a height $z$ due to jet-star interactions in steady state is

[TABLE]

where $<{S_{\rm s}(m,t)}/{S_{\rm j}}>$ is the time-averaged fraction of jet area intercepted by one stellar interaction.

5.2.1 Inverse Compton emission

At electron energies $E^{\prime}_{\rm IC}$, for all jet parameters considered, $t^{\prime}_{\rm IC}\sim 10^{12}-10^{13}$ s. This is much larger than the characteristic acceleration time of any considered bubble, typically in the range $t^{\prime}_{\rm acc}\sim 10^{6}-10^{9}$ s. Therefore, the IC emission at the energy of interest will last until long after the bubble has been accelerated.

As the bubble propagates downstream in the jet, the former radiates an IC luminosity $L^{\prime}_{\rm b,IC}\sim E^{\prime}_{\rm b,max}/t^{\prime}_{\rm IC}(z)$ in the flow frame, where $E^{\prime}_{\rm b}$ is the total energy emitted by one single bubble in the same frame during its evolution

[TABLE]

$z_{\rm max}$ is the maximum height the bubble reaches before losing all energy, or the jet maximum height, $H$. This approximation implicitly assumes that the energy distribution of the non-thermal electrons is $\propto E^{-2}$ (see Sect. 6.2).

Characteristic times for energy losses evaluated at $E^{\prime}_{\rm IC}$, as well as $t_{\rm esc}$, are plotted in Fig. 3. The energy losses are dominated by adiabatic expansion up to the jet recollimation point. At higher $z$, the shortest timescale is $t^{\prime}_{\rm esc}$, meaning the bubbles exit the jet before emitting (mostly through IC) all of their available energy.

We note that, contrary to the study of synchrotron gamma-ray emission below, IC energy losses are computed taking $B=0.1B_{\rm eq}$, which is a reasonable value for a jet with power dominated by its kinetic energy. However, since the energy losses at $E^{\prime}_{\rm IC}$ are dominated by non-radiative processes, considering more (or less) intense magnetic fields does not significantly affect the results.

We can estimate the luminosity detected by an observer, generated by the whole bubble population, averaged in time as

[TABLE]

where $\delta_{\rm{j}}$ is the Doppler boosting factor. This is a valid description as long as there is more than one dominant bubble simultaneously emitting333The Doppler-boosted luminosity of one single bubble is $L^{app}_{\rm b}=\delta^{4}L^{\prime}_{\rm b}$. If more than one star is within the jet emitting at different heights, the photons emitted by the distribution of sources at the same time in $K^{\prime}$ are not observed simultaneously in $K$, leading to $\delta^{4}\rightarrow\delta^{3}/\Gamma$ (see Sikora et al., 1997). Thus, if insufficient events take place at the same time, the apparent luminosity of the whole population can be reduced to that of one single bubble., meaning that the source duty cycle is larger than 1. In addition, the emission generated by the whole population, $L_{\rm pop}$, as seen by the observer will be larger than the emission produced by one single, dominant event, $L_{\rm dom}$.

The predicted apparent luminosities of the IC emission are a few times $10^{40}$ erg s*-1*. These values are of the order of the steady state emission generated by the whole stellar population within the jet.

5.2.2 Synchrotron emission

If we evaluate the energy losses at $E^{\prime}_{\rm Sy}$, $t^{\prime}_{\rm Sy}$ is dominant at all heights for $B^{\prime}=B^{\prime}_{\rm eq}$ (see Fig. 3). The synchrotron emission is so intense that the characteristic emission time is lower than the acceleration time of the bubble, $t^{\prime}_{\rm acc}$. Therefore, electrons of energy $E^{\prime}_{\rm Sy}$ immediately radiate via synchrotron emission the energy acquired via particle acceleration triggered by the jet-bubble interaction. Electrons in the flow frame are accelerated at a rate $\xi qB^{\prime}c$, where $\xi$, a free parameter representing the acceleration efficiency, is fixed to 0.1. Such an efficiency is sufficient to produce observable 100 MeV synchrotron photons, although for $\xi$ well below 0.1 synchrotron emission would not reach gamma-ray energies (X-rays of luminosities similar to those of IC in Table 1 could still be produced).

Synchrotron 100 MeV photons are produced when a single bubble enters the jet and is accelerated, following a flare-like lightcurve (see Barkov et al., 2012). From Khangulyan et al. (2013), the apparent total energy emitted by a bubble in the synchrotron fast cooling regime is:

[TABLE]

where we take $\bar{F_{\rm e}}=0.2$ (Barkov et al., 2010), and $F_{\rm rad}$ is the efficiency with which the particle loses energy through synchrotron radiation. Doppler boosting is already accounted for in Eq. 15 for the energy radiated by one single bubble.

Unlike in the case of the IC emission, synchrotron $L_{\rm peak}$ is 1–3 orders of magnitude higher than the average luminosity $L_{\rm pop}$, with

[TABLE]

where we take $F_{\rm e,max}=0.4$ (Barkov et al., 2010), and $r_{\rm b}$ is the radius of the bubble once it has expanded and reached relativistic velocities (Khangulyan et al., 2013). The apparent duration of this intense emission can be roughly estimated as $t_{\rm peak}=E_{\rm b}/L_{\rm b}$, which shows that the emission is highly variable, like intense, short flares occurring a few times per century or millennia.

For the high $\xi$ adopted, the synchrotron 100 MeV emission could reach luminosities of $\sim 10^{44}$ erg s*-1* for the most powerful jets considered, even under low jet Lorentz factors. As the emission is intense but short, we would have duty cycles of $PR\cdot t_{\rm obs}\sim 10^{-4}-10^{-2}$, depending on $L_{\rm j}$ and $\Gamma_{\rm j}$.

5.2.3 Synchrotron emission at low magnetic fields

Synchrotron emission at $\sim 100$ MeV is highly dependent on the value of the magnetic field: high magnetic fields yield an intense emission radiated in a very short amount of time, yielding a low duty cycle.

In order to maximize the duty cycle, we take the lowest possible value of the magnetic field, $B=B_{\rm min}$, for which electrons with $E^{\prime}_{\rm Sy}$ still cool more efficiently through synchrotron emission than through non-radiative losses. Looking at Fig. 3, this would result in values of $t^{\prime}_{\rm Sy}$ just below $t^{\prime}_{\rm ad}$ and/or $t^{\prime}_{\rm esc}$. The obtained magnetic field values, $B_{\rm min}$, are listed in Table 1. For these values of $B$, luminosity is up to $\sim 10^{42}$ erg s*-1*, but the duty cycle can be increased up to $PR\cdot t_{\rm obs}>1$.

We note that in this case with $B=B_{\rm min}$, $t^{\prime}_{\rm sy}>t^{\prime}_{\rm acc}$, meaning that the synchrotron emission is computed as described in Sect. 5.2.1, that is, as is done for inverse Compton radiation.

6 Discussion

The dynamical evolution of a bubble of stellar material expelled within the jet is semi-quantitatively well-described by the analytical estimates used in Sect. 4, and only differs from simulations by a moderate numerical factor. Our numerical estimates on gamma-ray energy production are not significantly affected by these differences as long as the bubble Lorentz factor eventually reaches $\sim\Gamma_{\rm j}$. Only in the case of intense Synchrotron emission at high magnetic fields, could there be a small reduction of luminosity, of a factor of a few, if indeed the bubble radius were overestimated in the long run.

Results presented in Table 1 indicate that the emission generated by stars penetrating the jet can be relatively persistent at high energies, through either inverse Compton or through synchrotron emission in the case of low magnetic fields. With large magnetic fields, emission at 100 MeV could be dominated by bright and infrequent events on top of the persistent, lower IC radiation.

The steady state emission of the whole population is unlikely to be detectable. We note that the similar values of $L_{\rm st}$ for all the explored jet configurations are due to the fact that it does not strongly depend on either $\Gamma_{\rm j}$ or $L_{\rm j}$. The small differences are due to jet geometry, which influences the amount of stars within the jet (e.g. a jet of $L_{\rm j}=10^{43}$ erg s*-1* recollimates at low heights), or due to differences in $F_{\rm rad}$ (e.g. under high $B$ values, for $L_{\rm j}=10^{45}$ the impact of synchrotron losses on electrons at $E^{\prime}_{\rm IC}$ is noticeable).

Inverse Compton emission of the bubbles at $E^{\prime}_{\rm IC}$ seems difficult to detect for the explored jet configurations unless $\eta\rightarrow 1$, or for a very nearby source. This is because non-radiative losses, or even synchrotron losses, dominate the process.

Synchrotron emission at $E^{\prime}_{\rm Sy}$ with high magnetic fields can result in bright, detectable flares, although the high luminosity implies a short duration. If we consider lower magnetic fields, persistent emission can be achieved, and luminosities of the order of $10^{42}$ erg s*-1* for $\Gamma_{\rm j}=10$.

There are some factors that may easily increase the radiation in the scenario studied here. For instance, for significantly lower values of $\rho_{\rm ISM}$ at $z\sim$ kpc, say \gtrsim 0.1\leavevmode\nobreak\cm*-3*, the emission energetics would grow by a factor of several due to the associated larger bubble mass (see Sect. 3.2), limited now by $M_{\rm in}$ (see Fig. 1). Furthermore, a younger, more massive red giant population, or the sporadic presence of an asymptotic giant branch (AGB) star within the jet, could also enhance the expected emission.

6.1 Younger red giant populations and AGB stars

Our results show that the emission produced by wind bubbles as they penetrate the jet is generally dominated by the most evolved stars within our modelled population (with $\dot{M}\sim 1.4\times 10^{-8}$ M*⊙* yr*-1*). The mass of the bubbles scales with mass-loss rate as $\propto\dot{M}^{3/2}$, and the luminosity emitted by the bubble interacting with the jet is proportional to its mass (except for the case of synchrotron emission of 100 MeV photons at $B_{\rm eq}$). Therefore, our results are scalable with mass-loss rate.

If we had a population of red giants of $M_{\rm RG}=1.5$ M*⊙* with a total mass comparable to that of the Galaxy, it would imply an increase of luminosity of a factor of $\sim 20$ in the case of low-B synchrotron emission, and of a factor of $\sim 60$ in the case of Inverse Compton emission. In the case of high-B emission, emission depends on the final radius of the bubble, and it would increase by a factor of $\sim 6$. In all cases, considering these less abundant stars would lead to a decrease of event occurences of $\sim 0.4$.

We note that in the rare occasion an AGB star penetrates an AGN jet, its wind bubble can inject into the latter up to $\sim 10^{-2}$ M⊙. This could potentially lead to a long duration event with IC luminosity $\sim 10^{44}$ erg s*-1*.

6.2 Caveats of the radiation estimates

In this work we estimate the energy radiated by accelerated particles at a given energy, $E^{\prime}_{\rm IC,sy}$, where the maximum of emission is expected to be. In assuming that all the available energy that goes into particle acceleration is radiated or lost at the mentioned energy, we are overestimating the overall emission. Depending on the energy distribution of the particles, this simple approximation can differ from a more precise calculation. In the case of inverse Compton emission, we expect this overestimation to be of a factor of a few, and in the case of synchrotron emission, of up to one order of magnitude (Vieyro et al., 2017; Vieyro et al. in prep). This has to be taken into account when reading the luminosities in Table 1. The approximations adopted are valid for electron energy distributions $\propto E^{2}$, typical for astrophysical sources. Harder electron distributions would still lead to similar results to those obtained but for very extreme cases, while steeper distributions would imply an even higher overestimate of the gamma-ray luminosity.

The radiation estimates were focused on gamma rays and a broad band study would require more detailed modelling. Nevertheless, it is worth exploring in the future the radio and X-ray synchrotron emission in the jet-bubble interactions.

7 Summary

In this work we have studied the gamma-ray emission produced when the bubbles formed by red giant winds penetrate the jet of a blazar in an elliptical galaxy, and described their dynamical evolution both analytically and through one illustrative simulation. We have shown that the analytical approximation is reasonably accurate at the present stage in the relativistic regime. We have found that the gamma rays produced by the red giant wind bubbles interacting with jets may reach detectable levels if Lorentz factors are high, non-thermal particles are generated very efficiently, and particle acceleration takes place at very high rates (to produce synchrotron 100 MeV photons). This predicted emission could be higher under the presence of an important population of younger red giants, in the rare event an AGB star enters the jet, or for relatively low values of $\rho_{\rm ISM}$. Unless $B=B_{\rm eq}$, duty cycles are not far below one, and in some cases, a few dominant bubbles could contribute simultaneously to the gamma-ray emission. For $B=B_{\rm eq}$, short bright synchrotron 100 MeV flares may be detected, with year or sub-year scale duration. Most of the known blazars in the local universe are hosted by elliptical galaxies. Therefore, provided that electrons are efficiently accelerated in bubble-jet interactions, it is plausible that some persistent or transient gamma-ray emission detected from the nearest blazars could originate in events like those described here.

Acknowledgments

We want to thank the anonymous referee for constructive and useful comments. We acknowledge support by the Spanish Ministerio de Economía y Competitividad (MINECO/FEDER, UE) under grants AYA2013-47447-C3-1-P, AYA2016-76012-C3-1-P, with partial support by the European Regional Development Fund (ERDF/FEDER), MDM-2014-0369 of ICCUB (Unidad de Excelencia ‘María de Maeztu’), and the Catalan DEC grant 2017 SGR 643. N.T-A. acknowledges support from MINECO through FPU14/04887 grant.