© ESO 2018

1. Introduction

High-mass microquasars (HMMQ) are X-ray binaries that host a massive star and a compact object (CO) from which jets are produced. The stellar wind can strongly influence the jet propagation, both because the jet has to propagate surrounded by wind material, and because the wind lateral impact may significantly bend the jet away from the star. Several authors have used numerical and analytical methods to study, on the scales of the binary, the interaction of HMMQ jets with stellar winds and their radiative consequences (e.g., Romero et al. 2003; Romero & Orellana 2005; Perucho & Bosch-Ramon 2008, Perucho & Bosch-Ramon 2012; Owocki et al. 2009; Araudo et al. 2009; Perucho et al. 2010; Yoon & Heinz 2015; Bosch-Ramon & Barkov 2016; Yoon et al. 2016).

On scales larger than the binary system, orbital motion should also affect the dynamics of the jet, making it follow a helical trajectory (e.g., Bosch-Ramon 2013; Bosch-Ramon & Barkov 2016). There are a few well-established microquasars in which evidence of a helical jet has been observed; for example, SS 433, 1E 1740.7−2942 and Cygnus X-3 (Abell & Margon 1979; Mioduszewski et al. 2001; Stirling et al. 2002; Miller-Jones et al. 2004; Luque-Escamilla et al. 2015; see also Sell et al. 2010 for the possible case of Circinus X-1). In the case of SS 433, the jet helical geometry likely originates in the accretion disk (e.g., Begelman et al. 2006); in the case of 1E 1740.7−2942, the system likely hosts a low-mass star and the role of its wind could be minor; in the case of Cygnus X-3, the apparent jet helical shape has an unclear origin. Bosch-Ramon & Barkov (2016) proposed that Cygnus X-1 and Cygnus X-3 could be affected by orbital motion and the stellar wind, although the uncertainties on the jet and wind properties make quantitative predictions difficult. Unlike the case of Cygnus X-3, no clear evidence of a helical or bent jet has yet been found for Cygnus X-1 (Stirling et al. 2001).

In this work we study the implications of the impact of the stellar wind and orbital motion for the nonthermal emission of a HMMQ, from radio to gamma rays. Using an analytical prescription for the jet dynamics based on the results of Bosch-Ramon & Barkov (2016), the nonthermal radiation from a helical jet (spectra, light curves, and morphology) is computed numerically, considering a leptonic model with synchrotron and inverse Compton (IC) emission.

The paper is organized as follows: in Sect. 2, the interaction between the stellar wind and the jet is described. The technical aspects of the model are explained in Sect. 3. The results of the calculations are presented in Sect. 4 and a discussion is provided in Sect. 5. Unless stated otherwise, the convention Q_x = Q/10^x is used throughout the paper, with Q in cgs units.

2. Jet-wind interaction

The physical scenario considered in this work consists of a massive star with a strong stellar wind and an accreting CO from which two jets are launched in opposite directions perpendicular to the orbital plane. The binary has a relatively compact orbit, taken circular for simplicity, with a period of P = 5 days and a separation of d_O = 3 × 10¹² cm (M₁ + M₂ ≈ 43 M_⊙). The luminosity of the star is taken to be L_⋆ = 10³⁹ erg s⁻¹, and its temperature T_⋆ = 4 × 10⁴ K (R_⋆ ≈ 10 R_⊙). The stellar wind is assumed to be spherically symmetric, with a mass-loss rate and a velocity of ${\dot{M}}_{\mathrm{w}}={10}^{-6}$ M_⊙ yr⁻¹ and ν_w = 2 × 10⁸ cm s⁻¹, respectively (typical for O-type stars; e.g., Muijres et al. 2012). Such a system would have similar parameters to those of Cygnus X-1, although it would be significantly wider and with a weaker wind than Cygnus X-3 (but the adopted wind and jet momentum rate relation will be similar; see below).

The orbit is assumed to lie on the xy-plane, with the star at the center of coordinates. The CO would be on the x-axis, orbiting the star counter clockwise, and the jet (counter jet) would initially point along $\widehat{z}$ ($-\widehat{z}$). The axes are defined in a frame corotating with the CO around the companion star. The jets are assumed to have a conical shape, with an initial half-opening angle θ_j0 = 0.1 rad ≈5.7°. In this work, we do not need to adopt a specific value for θ_j0, although to make our predictions qualitatively valid, it should be smaller than the jet deflection angle imposed by the wind impact (see below). We take an initial jet Lorentz factor of γ_j = 2, although the final results are not sensitive to this parameter as long as the jet is strongly supersonic and mildly relativistic. The jet power is taken as L_j = 3 × 10³⁶ erg s⁻¹.

Figure 1 sketches the scenario under study. The jet goes through three stages, the first of which is not present in the figure: Initially, the jet moves perpendicularly to the orbital plane. Then, the wind impact produces an asymmetric recollimation shock on the jet, which gets inclined away from the star (stage 2). Further from the binary, the combined effect of the stellar wind and orbital motion leads to a (orbit-related) force that bends the already inclined jet, now against orbital motion (clockwise in Fig. 1). In this way, the jet becomes helical, with asymmetric shocks being expected because of the orbit-related force (stage 3). Here it is assumed that electrons are efficiently accelerated in those shocks, ascribed here to one specific point: the onset of the helical jet, which neglects the fact that acceleration may occur along the helical jet region. Relativistic protons may also be generated in those shocks, but their radiation efficiency would be much lower in general, and therefore they are not considered in our calculations (see e.g., Bosch-Ramon & Khangulyan 2009, for a discussion of the different electron and proton cooling timescales in HMMQ).

Fig. 1. Schematic representation of the onset of the helical jet. The position, with respect to the star, of the change in the jet initial (bent) direction is r0, where the jet bending starts to be dominated by the added effects of the stellar wind and orbital motion. The axes are defined in a frame corotating with the CO.

For the radiation calculations, a fraction η_NT = 10⁻² of L_j is injected in the form of nonthermal electrons right where the jet starts its helical path. This fraction is only constrained by η_NT <  1, although the predicted luminosities are proportional to its value and easy to translate if a different fraction is adopted. Nevertheless, adopting η_NT = 10⁻² does not require particularly efficient acceleration processes. The inclination of the system with respect to the line of sight of the observer, assumed to be at a reference distance of d = 3 kpc, is characterized by the angle i. The values for the different fixed parameters of our model are listed in Table 1.

2.1. Wind impact: jet bending

Interaction between the stellar wind and the jet can lead to bending of the latter in the x-axis direction (stage 2 above). This bending can be characterized by an angle ϕ with respect to the jet initial direction (i.e., the z-axis; see Fig. 1). Following Bosch-Ramon & Barkov (2016) we introduce a nondimensional parameter, χ_j, which corresponds to the ratio between the wind momentum rate intercepted by the jet, and the jet momentum rate:$$\begin{array}{c}\hfill {\chi }_{\mathrm{j}}\approx \frac{{\mathrm{\theta }}_{\mathrm{j}0}{\dot{P}}_{\mathrm{w}}}{4\pi {\dot{P}}_{\mathrm{j}}}=\frac{{\mathrm{\theta }}_{\mathrm{j}0}{\dot{M}}_{\mathrm{w}}{\nu }_{\mathrm{w}}({\gamma }_{\mathrm{j}}-1)c}{4\pi L_{\mathrm{j}}{\gamma }_{\mathrm{j}}{\beta }_{\mathrm{j}}}\approx 0.9\frac{{\mathrm{\theta }}_{\mathrm{j}0,-1}{\dot{M}}_{\mathrm{w},-6}{\nu }_{\mathrm{w},8.3}}{L_{\mathrm{j},36.5}{\beta }_{\mathrm{j}}}\frac{({\gamma }_{\mathrm{j}}-1)}{{\gamma }_{\mathrm{j}}},\end{array}$$(1)

where ${\dot{M}}_{\mathrm{w}}$ is in units of M_⊙ yr⁻¹ in the rightmost expression, ${\dot{P}}_{\mathrm{w}}={\dot{M}}_{\mathrm{w}}{\nu }_{\mathrm{w}}$ is the stellar wind momentum rate, and ${\dot{P}}_{\mathrm{j}}$ is the jet momentum rate, which can be expressed in terms of L_j, γ_j, and the initial jet velocity in c units β_j0:$$\begin{array}{c}\hfill {\dot{P}}_{\mathrm{j}}=\frac{L_{\mathrm{j}}{\gamma }_{\mathrm{j}}{\beta }_{\mathrm{j}0}}{c({\gamma }_{\mathrm{j}}-1)}·\end{array}$$(2)

An approximate computation of the bending angle ϕ, with a discrepancy of less than 10% with respect to the numerical case, yields the following expression (for more details, see the Appendix in Bosch-Ramon & Barkov 2016):$$\begin{array}{c}\hfill \mathrm{\varphi }=\frac{{\pi }^2{\chi }_{\mathrm{j}}}{2\pi {\chi }_{\mathrm{j}}+4{\chi }_{\mathrm{j}}^{1/2}+{\pi }^2}·\end{array}$$(3)

For bending angles ϕ <  θ_j0, the lateral impact of the stellar wind may produce a weak recollimation shock or a sound wave in the jet, and the large-scale jet evolution should not be strongly affected (although the role of instability growth may still be important; see, e.g., Perucho & Bosch-Ramon 2008; Perucho & Bosch-Ramon 2012; Perucho et al. 2010 under uniform and clumpy stellar winds). However, when ϕ is significantly larger than θ_j0, a strong asymmetric recollimation shock is expected, with the jet substantially deviating from its initial direction. In such a case, the jet is likely to turn into a helical-shaped structure at larger scales (stage 3, above). Consequently, bending due to stellar wind impact becomes important for the jet evolution when ϕ >  θ_j0.

2.2. Orbital effects: helical jet

The orbital motion of the system combined with a jet inclined by the stellar wind (stage 2) leads, through the orbit-related force, to a jet helical pattern. This effect becomes clearer for significant bending angles. Otherwise, for ϕ <  θ_j0, the conical jet expansion will dominate the geometry of the resulting structure. We note that a helical pattern would be expected from an inclined jet even under ballistic jet propagation, but the presence of the stellar wind makes jet propagation nonballistic. Instability growth should also play a role, but here we consider its effects only phenomenologically through a slower jet flow, caused by the wind mass entrainment and hence deceleration associated to instabilities.

One can estimate the distance x_turn that a nonballistic jet travels in the x-axis direction before acquiring a negative velocity y-component, i.e., bending against the orbital motion. With the additional assumptions that the half-opening angle of the jet remains constant after being bent (see below), and that x_turn is much larger than the orbital separation, we get (see Bosch-Ramon & Barkov 2016, for the derivation of this expression)$$\begin{array}{c}\hfill x_{\mathrm{turn}}=\frac{K}{\sqrt{2{\chi }_{\mathrm{j}}}}\frac{{\nu }_{\mathrm{w}}}{\mathrm{\omega }}\approx 1.3\times {10}^{13}\frac{K}{\sqrt{{\chi }_{\mathrm{j}}}}\frac{{\nu }_{\mathrm{w},8.3}}{{\mathrm{\omega }}_{-5}}\mathrm{cm},\end{array}$$(4)

where ω is the orbital angular velocity, with the normalization corresponding to a period of several days (e.g., ω = 2π/5 days ≈ 1.4 × 10⁻⁵ rad s⁻¹), and K is a constant of O(1), introduced to account for the specific details of the wind–jet interaction. Here we fix K = 1.

Knowing x_turn and ϕ, we can compute the height z at which the jet starts forming the helical structure: z_turn = x_turn/tanϕ. For the parameters given above, and following Eqs. (3) and (4), the location vector of the beginning of the helical structure is r₀ = (1.58, 0, 3.56)×10¹³ cm, with ϕ = 19.8°. We note that x_turn = 1.28 × 10¹³ cm ≫d_O, as required. The distance from r₀ to the star is D_⋆ = 3.88 × 10¹³ cm and to the CO is l₀ = 3.77 × 10¹³ cm. The jet curvature due to orbital motion for distances smaller than l₀ is neglected, that is, y_turn = 0.

The point r₀ (the -helical- jet base hereafter) is the position where nonthermal electrons are injected. The evolution of these electrons is followed from this point onwards along the helical jet assuming that no more significant electron acceleration happens further downstream in the jet¹. Electron acceleration may also happen closer to the binary system (see Sect. 5), but here we focus specifically on the nonthermal emission produced in the helical jet region.

From r₀ onwards along the jet, the stellar wind should have a very important dynamical role. As the jet becomes helical due to the presence of wind plus orbital motion, the former is also pushing the latter in the radial direction. Therefore, the radial velocity (ν_r) of the jet flow may be as slow as the stellar wind itself, or higher if the shocked stellar wind is effectively accelerated radially by the jet. In addition, effects of wind entrainment are likely to affect the jet flow itself, slowing the latter down to an effective velocity ν_j <  β_j0c. We adopt here ν_r and ν_j as phenomenological parameters that fulfill ν_w ≲ ν_r ≲ ν_j. They provide a simplified prescription for the helical jet geometry, defining the pitch of the helical jet trajectory. It is also assumed that this trajectory is confined to the conical surface, characterized, in the nonrotating frame, by the direction determined by the vector r₀ along the orbit in that frame (i.e., a ring traced by the points r₀, plus the radial direction from the star), which should be a relatively good approximation as long as l₀ ≫ d_O.

The half-opening angle of the helical jet flow (θ_j) need not be constant, as the nonballistic nature of the jet trajectory, with weak-shock heating, rarefaction waves, wind (re)confinement, and mass-loading, may lead to (re)collimation or (re)widening of the helical jet, depending on the details of the jet and wind interaction. For simplicity, at this stage we assume θ_j = θ_j0. In addition, we also take ν_j = 5 × 10⁹ cm s⁻¹ and ν_r = 2 × 10⁹ cm s⁻¹; this is somewhat arbitrary, but except for the radio morphology (see Sect. 5) their actual values do not qualitatively change the results. Hydrodynamical simulations are required to properly characterize ν_r, ν_j, θ_j, the role of instability growth, mass-loading, and so on. Nevertheless, we consider our dynamical model realistic on scales of a few D_⋆, as large-scale perturbations should grow up to a size ≳D_⋆ to disrupt the helical structure.

3. Model description

We use a semi-analytical code to compute the energy and spatial distribution of the accelerated electrons. As the orbital period is significantly longer than the crossing time of the helical jet length (l_max/ν_j), the nonthermal electrons are assumed to be in the steady state at each orbital phase.

3.1. Energy losses

Analytical expressions for the energy losses are used in this work. For synchrotron cooling, assuming, for simplicity, an isotropic magnetic field B in the flow frame, one obtains, in cgs units (Longair, 1981):$$\begin{array}{c}\hfill {\dot{E}}_{\mathrm{sync}}\approx -\frac{4}{3}c{\sigma }_{\mathrm{T}}{\mathrm{\omega }}_{\mathrm{mag}}{\gamma }^2\approx 1.6\times {10}^{-5}{B}_{-2}^2{E}^2{\mathrm{erg}\mathrm{s}}^{-1},\end{array}$$(5)

where σ_T is the Thomson cross-section, ω_mag = B²/8π, and γ is the Lorentz factor of an electron with energy E. The magnetic field at the jet base is parametrized through the ratio of magnetic pressure to stellar photon energy density, η_B, as$$\begin{array}{c}\hfill \frac{B^2}{8\pi }={\eta }_{\mathrm{B}}\frac{L_{\star }}{4\pi c{D}_{\star }^2}·\end{array}$$(6)

This way of normalizing B allows an easy comparison between the expected luminosity of synchrotron and IC emission; for IC in the Thomson regime, the two processes contribute similarly for electrons of the same energy if η_B ∼ 1. However, for electrons of energy ≳ m_ec²/3k_BT_⋆ (with k_B being the Boltzmann constant), which scatter photons in the IC Klein-Nishina regime, η_B ≳ 0.1 already leads to global strong synchrotron dominance, which is corroborated by our calculations. One can relate η_B with the ratio of magnetic pressure to jet energy density (${\overline{\eta }}_{\mathrm{B}}$) for the same B via: ${\overline{\eta }}_{\mathrm{B}}={\theta }_{\mathrm{j}}^2({\eta }_{\mathrm{B}}/4)(L_{\star }/L_{\mathrm{j}})({\nu }_{\mathrm{j}}/c)$ (≈0.14 η_B for the values adopted here).

Losses due to IC are computed as in Khangulyan et al. (2014) for an isotropic electron population embedded in the stellar photon field taken as a black body of temperature T_⋆. This expression is valid in both the Thomson and the Klein-Nishina regimes and has an accuracy of ∼ 5%:$$\begin{array}{cc}& {\dot{E}}_{\mathrm{IC}}\approx -\frac{3{\sigma }_{\mathrm{T}}\pi m_e^2{c}^2{R}_{\star }^2{k}_{\mathrm{B}}^2{T}_{\star }^2}{2R^2{h}^3}G(t)\hfill \\ \hfill & \approx -3.8\times {10}^{-3}\frac{R_{\star ,12}^2{T}_{\star ,4.6}^2}{R_{14}^2}\frac{Eln(1+5.3E)}{1+183E}{\mathrm{erg}\mathrm{s}}^{-1},\hfill \end{array}$$(7)

with R_⋆ being the radius of the star, R the distance from the star to the emitter, $t=4E{k}_{\text{B}}{T}_{\star }/m_e^2{c}^4$, and$$\begin{array}{c}\hfill G(t)=\frac{4.62tln(1+0.156t)}{1+5.62t}·\end{array}$$(8)

For a conical jet, adiabatic losses are$$\begin{array}{c}\hfill {\dot{E}}_{\mathrm{ad}}=-\frac{2}{3}\frac{{\nu }_{\perp }}{r_{\mathrm{j}}}E\approx -6.7\times {10}^{-6}\frac{{\nu }_{\perp ,8}}{r_{\mathrm{j},13}}E{\mathrm{erg}\mathrm{s}}^{-1},\end{array}$$(9)

where ν_⊥ ≈ ν_j θ_j is the jet expansion velocity and r_j is the helical jet radius.

The total energy-loss rate is$$\begin{array}{c}\hfill \dot{E}={\dot{E}}_{\mathrm{sync}}+{\dot{E}}_{\mathrm{IC}}+{\dot{E}}_{\mathrm{ad}}.\end{array}$$(10)

3.2. Electron distribution

The nonthermal electrons are injected in the jet base with an energy distribution that follows a power law of index -2, typical for efficient nonthermal sources (and for acceleration in nonrelativistic strong shocks; Drury 1983), plus a cutoff at high energies set by energy and escape losses:$$\begin{array}{c}\hfill Q(E)\propto E^{-2}exp(-\frac{E}{E_{\mathrm{cut}}})·\end{array}$$(11)

We note that a much steeper Q(E) would mean inefficient electron acceleration, whereas a much harder one would imply exceptionally efficient acceleration in which most of energy is in the highest-energy electrons.

The normalization of the electron injection is taken as a fraction of the jet power available for the acceleration of nonthermal electrons, L_NT = η_NTL_j:$$\begin{array}{c}\hfill {\int }_{E_{\min }}^{E_{\max }}EQ(E)\mathrm{d}E=L_{\mathrm{NT}}.\end{array}$$(12)

The value of E_min is fixed here to 1 MeV for simplicity, although a very high E_min, greater than 0.1 − 1 GeV for example, would affect the more compact radio emission. The value of E_max is obtained from electron acceleration constraints.

Electrons are accelerated gaining energy at a rate ${\dot{E}}_{\mathrm{acc}}={\eta }_{\mathrm{acc}}eBc$, where η_acc is the acceleration efficiency, and e the electron charge. We take η_acc = 0.01, which would correspond, for example, to acceleration by a high-velocity (nonrelativistic) shock under Bohm diffusion (Drury 1983). Low acceleration efficiencies would reduce the maximum energy of synchrotron photons, for which the photon energy is ε ∝ E², much more than that of IC photons, produced in the Klein-Nishina regime at the highest spectral end, for which ε ∼ E (for IC in Thomson: ε ∝ E²).

Equating ${\dot{E}}_{\mathrm{acc}}=|\dot{E}|$ yields the maximum energy that electrons can reach in the accelerator (at r₀): E_acc. However, electrons can diffuse away from the accelerator before they reach E_acc. Therefore, the diffusion timescale, $t_{\text{d}}=r_{\text{j0}}^2/2D_{\text{B}}$, is also to be compared with the acceleration timescale, $t_{\mathrm{acc}}=E_{\mathrm{d}}/{\dot{E}}_{\mathrm{acc}}$, where r_j0 is the characteristic jet base radius, and D_B = E_dc/3qB is the diffusion coefficient assuming it proceeds in the Bohm regime. This gives a maximum energy value $E_{\mathrm{d}}=qB{r}_{\mathrm{j}}\sqrt{3{\eta }_{\mathrm{acc}}/2}$; for a larger diffusion coefficient D (still ∝E), $E_{\mathrm{d}}\propto \sqrt{D_{\mathrm{B}}/D}$. The E_cut value is obtained as the smallest among E_acc and E_d, whereas E_max can be taken as several times E_cut.

To compute the electron energy distribution along the jet, the latter is divided in cylindrical segments with increasing radius. The segment lengths are determined (with some exception; see below) by ν_j multiplied by one fifth of the local shortest cooling time, which is derived from $t_{\mathrm{loss}}=E_{\max }/|\dot{E}(E_{\max })|$. In this way, the energy and spatial evolution of the fastest evolving electrons is reasonably well sampled. The magnetic field is assumed to be mostly perpendicular to the flow motion, typical for jets far from their origin, and therefore B ∝ r_j0/r_j under frozen conditions and constant ν_j. Each individual jet division is treated as a homogeneous emitter, which in the worse case is correct within a ∼3% error (spatial scales grow segment by segment by a factor 1 + l_max/N/l₀ ∼ 1.03, where N is the number of segments).

For η_B → 1, synchrotron cooling times become very short (t_loss ≪ t_adv, with t_adv being the advection timescale), and the large N makes calculations very long. In those extreme cases, we have adopted an approximation which consists in limiting N and then introducing a correction factor, × t_loss/t_adv, to the injected population at the first segment. This simplified approach can overestimate the high energy part of the electron distribution by a factor around two, but speeds up calculations by several orders of magnitude for extreme B cases. Nevertheless, we indicate that η_B → 1 is too high, as in the region of interest the jet power is likely largely dominated by kinetic energy, and is considered here for illustrative purposes only.

To find the energy evolution up to a given segment, one computes the initial energy E₀ that electrons of energy E had when they were injected at the jet base. This is done iteratively by making small time steps back in time (and in space) to sample the energy evolution of electrons backwards from a given segment:$$\begin{array}{c}\hfill E_{i-1}=E_i-\dot{E}(E_i)\mathrm{d}t.\end{array}$$(13)

Here, E_i − 1 is the energy that an electron had at a time t − dt, and E_i is the electron energy at time t. We note that E_i − 1 >  E_i due to the negative sign of $\dot{E}$. The linear approximation used in Eq. (13) is only valid as long as E_i − 1 − E_i ≪ E_i.

Once E₀ is known, the electron energy distribution at the jet base, N₀(E₀), of a segment k with length dl_k, is computed as$$\begin{array}{c}\hfill N_0(E_0)=Q(E_0)t_k^{\mathrm{adv}},\end{array}$$(14)

with $t_k^{\text{adv}}=d{l}_k/{\nu }_{\text{j}}$ (if small enough; otherwise $t_k^{\text{adv}}\to t_{\text{loss}}$). The electron energy distribution at a segment k is computed from the distribution when at the segment k–1:$$\begin{array}{c}\hfill N_k(E_k)=N_{k-1}(E_{k-1})\frac{{\dot{E}}_k(E_{k-1})}{{\dot{E}}_k(E_k)},\end{array}$$(15)

where E_k and E_k − 1 are energies related by Eq. (13) in segments k and k–1, respectively². By applying this procedure down to the jet base, one gets the general expression for the electron energy distribution at each jet segment:$$\begin{array}{c}\hfill N_k(E_k)=N_0(E_0)\prod _{i=k}^1\frac{{\dot{E}}_i(E_{i-1})}{{\dot{E}}_i(E_i)},\end{array}$$(16)

keeping in mind that E_k ≳ m_ec² and E₀ <  E_max.

3.3. Emission of radiation

The power per photon energy unit radiated by an electron of energy E via synchrotron emission is computed following (Pacholczyk 1970)$$\begin{array}{c}\hfill P_{\epsilon }(E)=\frac{\sqrt{3}{e}^{3}Bsin\mathrm{\theta }}{h{m}_{\mathrm{e}}{c}^2}F(x),\end{array}$$(17)

where θ is the pitch angle between the electron velocity and B vectors, and x = ε/ε_c, with ε_c being the characteristic synchrotron photon energy defined as$$\begin{array}{c}\hfill {\epsilon }_{\mathrm{c}}=\frac{3}{4\pi }\frac{ehBsin\mathrm{\theta }}{m_e^3{c}^5}{E}^2.\end{array}$$(18)

The function F(x) includes an integral of the K_5/3(ζ) Bessel function, but in this work we adopt an approximation valid within a few percent error in the range 0.1 <  x <  10 (e.g., Melrose 1980; Aharonian 2004):$$\begin{array}{c}\hfill F(x)=x{\int }_x^{\infty }{K}_{5/3}(\zeta )\mathrm{d}\zeta \approx 1.85x^{1/3}{e}^{-x}.\end{array}$$(19)

The geometry of B is not well known, so an isotropic distribution in the fluid frame is assumed, taking for simplicity $Bsin\theta \approx \sqrt{<B^2>}=B\sqrt{2/3}$. The synchrotron spectral energy distribution (SED) for an isotropic population of electrons with distribution N(E) is:$$\begin{array}{c}\hfill \epsilon L_{\epsilon }^{\mathrm{sync}}=\epsilon {\int }_0^{\infty }{P}_{\epsilon }(E)N(E)\mathrm{d}E\end{array}$$(20)

The IC emission is computed considering that relativistic electrons interact with a beam of photons with energy ε₀. The IC kernel for an angle θ with respect to the initial beam direction is (Aharonian & Atoyan 1981):$$\begin{array}{c}\hfill \frac{{\mathrm{d}}^{2}N(\mathrm{\theta },\epsilon )}{\mathrm{d}\epsilon \mathrm{d}\mathrm{\Omega }}=\frac{3{\sigma }_T{m}_e^2{c}^4}{16\pi {\epsilon }_0{E}^2}[1+\frac{z^2}{2(1-z)}-\frac{2z}{b_{\mathrm{\theta }}(1-z)}+\frac{2z^2}{b_{\mathrm{\theta }}^2{(1-z)}^2}],\end{array}$$(21)

where $b_{\theta }=2(1-\cos \theta ){\epsilon }_{0}E/m_e^2{c}^4$, z = ε/E, and ε can vary in the limits ε₀ <  ε <  Eb_θ/(1 + b_θ). The (apparent) IC SED at a given direction is obtained by convolving Eq. (21) with the electron energy distribution N(E) and the energy distribution density of stellar photons n(ε₀):$$\begin{array}{c}\hfill \epsilon L_{\epsilon }^{\mathrm{IC}}=4\pi c{\epsilon }^2{\int }_0^{\infty }\mathrm{d}E{\int }_0^{\infty }\mathrm{d}{\epsilon }_0\frac{{\mathrm{d}}^{2}N(\mathrm{\theta },\epsilon )}{\mathrm{d}\epsilon \mathrm{d}\mathrm{\Omega }}N(E)n({\epsilon }_0).\end{array}$$(22)

Target photons other than those of stellar origin, like those coming from an accretion disk, are unimportant for the jet regions studied in this work. Synchrotron self-Compton is also negligible even when one assumes that most of the nonthermal energy is released as synchrotron radiation.

3.4. Absorption of radiation

The main absorption mechanism considered in this work is electron-positron pair production from gamma rays interacting with stellar photons. The cross section for this process is (Gould & Schréder 1967)$$\begin{array}{c}\hfill {\sigma }_{\gamma \gamma }=\frac{3}{16}{\sigma }_{\mathrm{T}}(1-{\beta }^2)[(3-{\beta }^4)ln\left(\frac{1+\beta }{1-\beta }\right)-2\beta (2-{\beta }^2)],\end{array}$$(23)

where$$\begin{array}{c}\hfill \beta =\sqrt{1-\frac{2m_e^2{c}^4}{\epsilon {\epsilon }_0(1-cos{\mathrm{\theta }}_{\gamma \gamma })}},\end{array}$$(24)ε and ε₀ are the energies of the high-energy and stellar photons, respectively, and θ_γγ is the angle between their propagation directions.

The optical depth of gamma-ray absorption is (Gould & Schréder 1967)$$\begin{array}{c}\hfill {\tau }_{\gamma \gamma }(\epsilon )={\int }_0^l\mathrm{d}l[1-cos{\mathrm{\theta }}_{\gamma \gamma }(l)]{\int }_{{\epsilon }_{0,\min }}^{\infty }\mathrm{d}{\epsilon }_{0}n({\epsilon }_0,l){\sigma }_{\gamma \gamma }({\epsilon }_0,\epsilon ,{\mathrm{\theta }}_{\gamma \gamma }),\end{array}$$(25)

where l is the distance to the observer covered by the gamma ray, and ${\epsilon }_{0,\text{min}}=2m_e^2{c}^4/\epsilon (1-\cos {\theta }_{\gamma \gamma })$ is the minimum energy of the gamma rays absorbed by the target stellar photons.

For the cases studied here, synchrotron self-absorption is estimated in radio by computing the optically thick-to-optically thin transition frequency (e.g., Bosch-Ramon 2009), and is found to be negligible in general; it has therefore not been included in the calculations. Free-free radio absorption in the stellar wind is also negligible at the scales of the emitting region. Nonthermal ultraviolet (UV) photons are also strongly absorbed on their way to the observer, but we have not included this effect as UV photons from the star would be largely dominant in any case. In a very compact system with a very powerful stellar wind, such as, for example, Cygnus X-3, radio absorption could still be an issue at a distance D_⋆ from the binary, although we defer specific source studies for future work.

4. Results

The results are obtained using the parameter values given in Table 1, together with the sets of values presented in Table 2 for i, the orbital phase α, and η_B, whose variations have strong effects on the nonthermal emission. Considering the free parameters, the electron energy distribution is affected only by η_B, whereas SEDs, light curves, and the emitter geometry are also affected by α and i. Most SEDs are calculated for α = 0.25 as a typical example, which corresponds to one of the orbit nodes. Varying α, we also compute some SEDs and IC light curves to illustrate orbital variations of this component. A comparison of SEDs and light curves is done between the extended and the point-like (one-zone) emitter, to check which features arise when an extended (helical) jet is considered. Finally, maps at 5 GHz are computed to illustrate the expected radio morphologies from a helical jet.

4.1. Nonthermal electrons

The electron energy distribution for the jet³, ×E² to emphasize the energy content at different energy scales, is shown in Fig. 2 for η_B = 10⁻⁴ and 1, and up to three different jet lengths starting from r₀: D_⋆, 3D_⋆ and 25D_⋆ (whole jet). As seen in the figure, a significant fraction of electrons below ∼1 GeV survive beyond 3D_⋆, and for the lowest B, TeV electrons can also reach a larger distance from the binary (see also, e.g., Khangulyan et al. 2008). As expected, the energy distribution is steep for electron energies dominated by synchrotron and Thomson IC losses. When IC dominates, the SED becomes harder at ∼10 GeV because of the Thomson-to-Klein-Nishina transition of the electron distribution, getting steeper again only at the highest energies due to synchrotron losses (and the effect of E_cut).

Fig. 2. Electron energy distribution of the jet ×E2 up to D⋆ (dashed lines), 3D⋆ (dotted lines) and 25D⋆ (the whole jet; solid lines), for ηB = 10−4 and 1.

4.2. Spectral energy distribution

Figure 3 shows the SEDs of the synchrotron and the IC emission for the same cases studied in Fig. 2, showing now within each panel the jet and the counter-jet contributions separately. An inclination of i = 60° is adopted, and α = 0.25. As expected from the severe drop in radiation efficiencies, the contribution to the emission of the regions beyond 3D_⋆ is minor for all the cases studied, although radio emission is still non-negligible beyond that distance, as seen in Sect. 4.5. It is worth noting that most of the IC radiation comes from the counter-jet because of the larger IC θ-values. Doppler boosting is negligible for the adopted ν_j-value and is not included.

Fig. 3. Synchrotron and IC SEDs of the jet (blue) and the counter-jet (green), for i = 60° and α = 0.25, up to D⋆ (dashed lines), 3D⋆ (dotted lines), and 25D⋆ (the whole jet; solid lines), for ηB = 10−4 (right panel) and 1 (left panel).

The effect of the magnetic field on the synchrotron and IC SEDs is shown in Fig. 4 for i = 30° and α = 0.25, and for the jet and the counter-jet separately. As expected, the synchrotron to IC luminosity ratio depends strongly on η_B. Figure 4 also shows the effect of gamma-ray absorption (stronger for i = 30° than for 60°), which is quite modest given the relatively large distances from the star and the α-value considered in this plot.

Fig. 4. Synchrotron and IC SEDs of the jet (solid lines) and the counter-jet (dashed lines), for i = 30° and α = 0.25, and ηB = 10−4 (blue lines), 10−2 (green lines), and 1 (red lines). The unabsorbed IC emission is also shown (black dotted lines).

Figure 5 shows the effect of the orbital phase on the IC spectrum. With i = 30° and η_B = 10⁻², the SEDs for α = 0, 0.25 and 0.5 are plotted. For α = 0.5, which corresponds to the CO superior conjunction (i.e., the CO is behind the star), there is a drop of up to five orders of magnitude in the counter-jet emission above 10 GeV due to gamma-ray absorption, which makes the jet dominate the radiation output in this range, effectively smoothing the overall impact of gamma-ray absorption for this orbital phase.

Fig. 5. IC SEDs of the jet (solid lines) and the counter-jet (dashed lines), for i = 30°, ηB = 10−2, α = 0 (red lines), 0.25 (blue lines) and 0.5 (green lines). The unabsorbed IC emission is also shown (black dotted lines).

4.3. Light curve

Figure 6 shows the IC light curves for η_B = 10⁻². The peak around α = 0.5 and 0.1 − 100 GeV for i = 60° is explained by the larger IC θ-angles, and the small dip for i = 30° manifests some level of gamma-ray absorption, which starts at ∼10 GeV in a relatively hard IC SED (see Fig. 5). For energies > 100 GeV, absorption has a very strong impact on the fluxes by a factor of several, although the jet component smoothens somewhat the light curve around α = 0.5 (superior conjunction of the CO) because of reduced gamma-ray absorption due to a smaller θ_γγ-value (see Fig. 5). On the other hand, the counter-jet still leads to significant IC fluxes around α = 0 (inferior conjunction; see Fig. 5), unlike an emitter located very close to the CO, which would lead to stronger differences in flux.

Fig. 6. Light curves of IC emission (jet + counter-jet) at 0.1–100 GeV (left panel) and > 100 GeV (right panel), for i = 30° (red lines) and 60° (blue lines), and ηB = 10−2. The unabsorbed light curves are also shown (dotted lines). The CO is in inferior (superior) conjunction with the star when the orbital phase is 0 (0.5).

4.4. Extended versus one-zone emitter

Figure 7 shows a comparison between the SEDs of an extended (helical) jet and a one-zone emitter, for the jet and the counter-jet, for i = 30°, α = 0.25, and η_B = 10⁻². The one-zone emitter is assumed to be just one segment located at r₀ and of length D_⋆/3, which sets the electron escape time. Because of fast radiation losses, the one-zone approximation works well in the high-energy part of the synchrotron and IC SEDs, as expected. Qualitatively, the extended and the one-zone emitters look quite similar, although differences of up to a factor of 2 are seen, mainly but not only in the jet synchrotron and IC components at low energies, as low-energy electrons escaped from the one-zone region are still radiatively relevant. We note that inaccurate source knowledge and model simplifications imply systematic uncertainties likely larger than a factor of 2, which means that the SED alone cannot help to discriminate a helical jet model from that of a one-zone emitter.

Fig. 7. Comparison of the (absorbed) SEDs of an extended (helical) jet (solid lines) and a one-zone emitter (dashed lines) for the jet and the counter-jet for i = 30°, α = 0.25, and ηB = 10−2.

Figure 8 compares the light curves of an extended and a one-zone emitter in the cases presented in Fig. 6. The fluxes are larger (lower) for the extended emitter for 0.1 − 100 GeV (> 100 GeV), and are symmetric around α = 0.5. The different behavior between the two energy bands is caused by the fact that electrons emitting > 100 GeV photons via IC reach farther from the star than those emitting 0.1 − 100 GeV photons (see Sect. 4.1), both for the counter-jet and the jet, decreasing the IC flux as they meet a more diluted target photon field. We note that for higher jet powers, the helical jet would start further away (see Eqs. (1) and (4)), in which case the light curve would become asymmetric for the extended jet.

Fig. 8. As in Fig. 6 but comparing the extended jet (solid lines) and the one-zone case (dashed lines), and for the absorbed gamma-ray component only.

4.5. Radio emission

Radio sky maps at 5 GHz for η_B = 10⁻² and different orbital phases are shown in Figs. 9 and 10 for i = 30° and i = 60°, respectively. These maps are obtained following three steps: first, the helical jet structure is projected into the plane of the sky. Subsequently, the emission from each segment is convolved with a Gaussian with standard deviation σ = r_j/2 to approximately distribute the emission over the actual sky-projected area of each segment. Finally, this emission is convolved again with a Gaussian with FWHM = 1 mas in order to mimic the response of a radio, very long baseline interferometer (see, e.g., Walker et al. 1995 for a VLBA description). With the sensitivity of current instrumentation, typically of a few tens of μJy beam⁻¹, the extended emission could be resolved, but only from the parts of the jet that are relatively close to its base. The bending of the jet appears more dramatic the smaller the inclination due to projection effects. The total received flux for each map is 7.4 mJy. As the flux increases proportionally to η_NT, one could detect the jet up to larger distances for higher η_NT. To illustrate this, Fig. 11 shows a comparison between two sky maps with different nonthermal energy fractions, η_NT = 10⁻² and η_NT = 5 × 10⁻². We note that the radio fluxes also depend on η_B as $\propto {\eta }_{\mathrm{B}}^{3/4}$.

Fig. 9. Sky maps at 5 GHz for i = 30°, ηB = 10−2, and orbital phases of 0, 0.25, 0.5 and 0.75 from left to right. Contour levels correspond to fluxes of 0.032, 0.1, 0.32 and 1 mJy beam−1. The beam size is 1 mas.

Fig. 10. As in Fig. 9 but for i = 60°.

Fig. 11. Sky maps at 5 GHz for i = 60°, α = 0.25, ηNT = 10−2 (left) and 5 × 10−2 (right). Contour levels are (mJy beam): 0.032, 0.1, 0.32, 1, 3.2. The beam size is 1 mas.

5. Summary and discussion

As argued by Bosch-Ramon & Barkov (2016), a nonballistic helical jet region is the likely outcome of a HMMQ jet interacting with the stellar wind under the effect of orbital motion. In the present work, we have shown that such a region is predicted to produce significant fluxes from radio to gamma rays if L_NT ∼ η_NTL_j ≳ 10³⁴η_NT, −2 erg s⁻¹, for a source at a few kiloparsecs. These significant fluxes are explained by the high synchrotron and IC efficiencies at the helical jet onset location, r₀, for typical HMMQ L_⋆-values. Since realistic η_B-values are expected to be well below equipartition with the radiation field, the IC component is likely to dominate the nonthermal emission, peaking around 10 GeV. Specific gamma-ray light curve features are also predicted: a non-negligible impact of gamma-ray absorption, combined with angular effects for this process and IC. Peculiar changing radio morphologies, which can trace the jet helical structure, are expected as well.

We note that hydrodynamical instabilities and wind-jet mixing make our results mostly valid for the inner helical jet region, say up to a few D_⋆. Beyond that point, the jet helical geometry is likely to become blurred, eventually turning into a bipolar, relatively wide, and collimated supersonic outflow; a mixture of jet and wind material. Although the flow may still be mildly relativistic for very energetic ejections, in general the resulting “jet” is expected to be much faster than the stellar wind, but nonrelativistic.

Our results are based on a number of strong simplifications, in particular regarding the stability of the helical jet. Therefore, detailed numerical simulations of the jet-wind interaction on middle to large scales are necessary for more accurate predictions. Nevertheless, despite the simplifications adopted, the main features of the scenario are expected to be rather robust, at least at a semi-quantitative level. These are: non-negligible fluxes; specific light curves characterized by the value of r₀, the interplay of jet and counter-jet emission, and angular IC and gamma-ray absorption effects; and curved jet radio morphologies with characteristic orbital evolution. This is so because the first- or even zeroth-order level of the dynamical, radiative, and geometrical effects involved are taken into account in our calculations.

Two important system parameters, eccentricity and d_O, have been assumed constant in our calculations. Despite being difficult to ascertain the impact of eccentricity in a HMMQ without numerical calculations, simulations of pulsar-star wind colliding in eccentric binaries (Barkov & Bosch-Ramon 2016; Bosch-Ramon et al. 2017) suggest that the effect of eccentricity should be important only for rather eccentric binaries. In that case, the nature of the interaction structure may not be helical at all, being not only strongly asymmetric, and inclined towards the apastron side of the orbit, but also likely much more sensitive to disruption. The value of d_O, and thus P, may also change. This would leave the helical jet geometry unaffected: both the helical vertical step and x_turn are $\propto P\propto d_{\mathrm{O}}^{3/2}$. On the other hand, the luminosity would change as $\propto P^{-1}\propto d_{\mathrm{O}}^{-3/2}$ except for emission produced under the fast radiative cooling regime of electrons (in our setup the highest-energy end of the synchrotron and IC emission), which is independent of d_O. The helical jet radiation in a very wide system would be hardly detectable unless the flow is slow, which allows electrons to radiate more energy. On the other hand, in a very compact system, radio emission would be hard to resolve for a slow helical jet flow, as the source would be too small.

The emission on the scales of the binary, produced for instance by jet internal shocks or by the first wind-induced recollimation shock in the jet (stage 2), could be important. In which case: (i) assuming that the same nonthermal luminosity is injected in both regions, the radio emission from the binary scales would be strongly absorbed via synchrotron self-absorption and wind free-free absorption; and (ii) the X- and gamma-ray emission would be at a similar level if produced in the fast radiation cooling regime, and at a higher level otherwise because of the higher B-value and IC-target density. Gamma-ray absorption could however strongly attenuate the gamma-ray luminosity > 100 GeV for most of the orbit. Therefore, the helical jet could contribute significantly to the overall nonthermal radiation of HMMQ, and particularly in radio and > 100 GeV, even if electrons are also accelerated closer to the jet base. This was the main motivation for this work: to perform a first exploration of specific radiation features of the helical jet region, so that it could be disentangled from the other emitting sites. For the HMMQ Cygnus X-1 and Cygnus X-3, their uncertain wind and jet parameters make it difficult to make concrete predictions; for example, a non-ballistic jet region may not form at all (Yoon et al. 2016; Bosch-Ramon & Barkov 2016 and references therein), and an individual detailed source analysis is out of the scope of this work. In spite of this, our findings are not in contradiction with the radio morphology and the gamma-ray light curves of these sources (e.g., Stirling et al. 2001; Mioduszewski et al. 2001; Miller-Jones et al. 2004; Zanin et al. 2016; Zdziarski et al. 2018). Future detailed modeling together with deep, high-resolution radio observations and the high-energy light curves could provide hints of helical jet emission.

The wind-jet interaction in HMMQ resembles, to a significant extent, the wind-wind interaction in high-mass binaries hosting young pulsars. The latter sources are probably extended mostly on the orbital plane and along the orbit semi-major axis, whereas the bipolar helical-jet structure is focused mostly in a direction perpendicular to the orbit plane (Bosch-Ramon et al. 2015, 2017; Bosch-Ramon & Barkov 2016). Nevertheless, the interactions of the stellar wind with a relativistic pulsar wind or jet, on binary, middle, and large scales, share many qualitative properties, and this may mask a fundamentally different engine accretion vs. a pulsar wind when observing sources of yet unknown CO. These similarities are: extended radio structures with orbital evolution found at mas scales; strong radio absorption on the binary scales; similar hydrodynamics, and therefore impact of instabilities, adiabatic losses, and Doppler boosting effects⁴, for the shocked flows outside the binary region; and gamma-ray light curves strongly affected by IC and gamma-ray absorption angular effects, in a likely extended/multi-zone emLitter. Radio observations, high-quality multiwavelength data, numerical simulations, and realistic radiation calculations (e.g., Aharonian et al. 2006; Perucho & Bosch-Ramon 2008, 2012; Albert et al. 2009; Fermi LAT Collaboration 2009; Perucho et al. 2010; Moldón et al. 2011; Moldón et al. 2011; Moldón et al. 2012; Bosch-Ramon et al. 2015, 2017; Dubus et al. 2015; Zanin et al. 2016; Yoon & Heinz 2015; Yoon et al. 2016; de la Cita et al. 2017; Wu et al. 2018; Barkov & Bosch-Ramon 2018) are therefore needed to characterize radiation scenarios in high-mass binaries with unknown CO, such that empirical observables can help to disentangle the CO nature.

Acknowledgments

We would like to thank the referee for a constructive and useful report that helped to improve the manuscript. We acknowledge support by the Spanish Ministerio de Economía y Competitividad (MINECO/FEDER, UE) under grants AYA2013-47447-C3-1-P and 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. EM acknowledges support from MINECO through grant BES-2016-076342.

References

All Tables

All Figures

	Fig. 1. Schematic representation of the onset of the helical jet. The position, with respect to the star, of the change in the jet initial (bent) direction is r0, where the jet bending starts to be dominated by the added effects of the stellar wind and orbital motion. The axes are defined in a frame corotating with the CO.
In the text

	Fig. 2. Electron energy distribution of the jet ×E2 up to D⋆ (dashed lines), 3D⋆ (dotted lines) and 25D⋆ (the whole jet; solid lines), for ηB = 10−4 and 1.
In the text

	Fig. 3. Synchrotron and IC SEDs of the jet (blue) and the counter-jet (green), for i = 60° and α = 0.25, up to D⋆ (dashed lines), 3D⋆ (dotted lines), and 25D⋆ (the whole jet; solid lines), for ηB = 10−4 (right panel) and 1 (left panel).
In the text

	Fig. 4. Synchrotron and IC SEDs of the jet (solid lines) and the counter-jet (dashed lines), for i = 30° and α = 0.25, and ηB = 10−4 (blue lines), 10−2 (green lines), and 1 (red lines). The unabsorbed IC emission is also shown (black dotted lines).
In the text

	Fig. 5. IC SEDs of the jet (solid lines) and the counter-jet (dashed lines), for i = 30°, ηB = 10−2, α = 0 (red lines), 0.25 (blue lines) and 0.5 (green lines). The unabsorbed IC emission is also shown (black dotted lines).
In the text

	Fig. 6. Light curves of IC emission (jet + counter-jet) at 0.1–100 GeV (left panel) and > 100 GeV (right panel), for i = 30° (red lines) and 60° (blue lines), and ηB = 10−2. The unabsorbed light curves are also shown (dotted lines). The CO is in inferior (superior) conjunction with the star when the orbital phase is 0 (0.5).
In the text

	Fig. 7. Comparison of the (absorbed) SEDs of an extended (helical) jet (solid lines) and a one-zone emitter (dashed lines) for the jet and the counter-jet for i = 30°, α = 0.25, and ηB = 10−2.
In the text

	Fig. 8. As in Fig. 6 but comparing the extended jet (solid lines) and the one-zone case (dashed lines), and for the absorbed gamma-ray component only.
In the text

	Fig. 9. Sky maps at 5 GHz for i = 30°, ηB = 10−2, and orbital phases of 0, 0.25, 0.5 and 0.75 from left to right. Contour levels correspond to fluxes of 0.032, 0.1, 0.32 and 1 mJy beam−1. The beam size is 1 mas.
In the text

	Fig. 10. As in Fig. 9 but for i = 60°.
In the text

	Fig. 11. Sky maps at 5 GHz for i = 60°, α = 0.25, ηNT = 10−2 (left) and 5 × 10−2 (right). Contour levels are (mJy beam): 0.032, 0.1, 0.32, 1, 3.2. The beam size is 1 mas.
In the text