# Modeling Particle Transport in Astrophysical Outflows and Simulations of Associated Emissions from Hadronic Microquasar Jets

Authors: D. A. Papadopoulos; O. T. Kosmas; S. Ganatsios
Journal: Advances in High Energy Physics
 (2022)

---

## Abstract

In this work, after improving the formulation of the model on particle transport within astrophysical plasma outflows and constructing the appropriate algorithms, we test the reliability and effectiveness of our method through numerical simulations on well-studied galactic microquasars as the SS 433 and the Cyg X-1 systems. Then, we concentrate on predictions of the associated emissions, focusing on detectable high-energy neutrinos andγ-rays originated from the extragalactic M33 X-7 system, which is an X-ray binary discovered in 2006, located in the neighboring galaxy Messier 33, and has not yet been modeled in detail. The particle and radiation energy distributions, produced from magnetized hadronic astrophysical jets in the context of our method, are assumed to originate from decay and scattering processes taking place among the secondary particles created when hot (relativistic) protons of the jet scatter on thermal (cold) ones (p-p interaction mechanism inside the jet). These distributions are computed by solving the system of coupled integrodifferential transport equations of multiparticle processes (reactions chain) following the inelastic proton-proton (p-p) collisions. For the detection of such high-energy neutrinos as well as multiwavelength (radio, X-ray, and gamma-ray) emissions, extremely sensitive space telescopes and other γ-ray and neutrino detection instruments are in operation or have been designed like the CTA, IceCube, ANTARES, KM3NeT, and IceCube-Gen-2.


---

## Body

## 1. Introduction

During the last few decades, collimated astrophysical outflows have been observed to emerge from a wide variety of galactic and extragalactic compact structures Mirabel_Rodriguez_1999, Romero_Torres_et_al_2003, Reynoso08, Reynoso09, RomeroReview. Among such objects, the stellar scale class of microquasar (MQ) and X-ray binary systems (XRBs) possesses prominent positions [1–5]. These two-body cosmic structures consist of a collapsed stellar remnant, a stellar-mass black hole or a neutron star (compact object), and a companion (donor) main-sequence star in coupled orbit around their center of mass.Due to the strong gravitational field pertaining around the compact object, mass from the companion star is accreted onto the equatorial region of the black hole forming an accretion disc. In the cases when the black hole is rotating rapidly and the accretion disc is geometrically rather thick and hot, ejection of two powerful oppositely directed mostly relativistic mass outflows (jets) occurs perpendicular to the accretion disk [6]. Such systems constitute excellent “laboratories” for the investigation of astrophysical outflows (jets). Currently, for the detection of multiwavelength (radio, X-ray, and gamma-ray) as well as high-energy neutrino emissions, sensitive detection instruments are in operation or have been designed like the IceCube, ANTARES, KM3NeT, CTA, and IceCube-Gen2 [7–13].From a theoretical and phenomenological viewpoint, due to the presence of rather strong magnetic fields, the MQ jets are treated as magnetohydrodynamical flows emanating from the vicinity of the compact object (usually a stellar-mass black hole) [14–16], and this is also assumed to be the case in our present work. From the observed characteristics of MQs, researchers have concluded that they share a lot of similarities in their physical properties with the class of Active Galactic Nuclei (AGN), which are mostly located at the central region of galaxies. The latter, however, is enormously larger in scale compared to microquasars, and in addition, the evolution of AGNs is appreciably slower than that of MQs that makes the observation of many phenomena significantly difficult [17–22].In this article, we focus on the magnetized astrophysical outflows (jets) that are characterized by hadronic content (p,π±, light nuclei, etc.) in their jets. We consider them as fluid flows emanating from a central source at the jet’s origin [4, 5, 23]. We, furthermore, assume that the compact object is a spinning black hole with mass up to few tens (30-50) of the Sun’s mass [1, 5, 24].Initially, we make an effort to improve the model employed to perform numerical simulations for particle and radiation emissions from hadronic MQ jets [25, 26]. We mention that some authors employ leptohadronic models which consider that primary electrons (and/or positrons) carry significant portion of the kinetic power inside the jets (see Section 2) [6, 24, 27]. In this paper, we adopt that the main mechanism producing the high-energy band of the electromagnetic spectrum (GeV-TeV γ-rays) and high-energy neutrinos is the proton-proton (p-p) collision taking place within the hadronic jets, i.e., the inelastic scattering of hot (nonthermal) protons on thermal (cold) ones [24].The scattering, diffusion, decay, etc., of the secondary particles (π±, K±, μ±, e±, etc.) produced afterwards, from a mathematical modeling point of view, are governed by a system of coupled integrodifferential transport equations. One of the purposes of this work is to attempt, for a first time, to formulate in a compact way this differential system of transport equations (satisfied by the primary protons and the secondary multiparticles, multispecies) and derive advantageous algorithms [28] to perform the required simulations for persistent emissions (transient episodes or flares are ignored [6]) emanating from the hadronic component of the MQ jets. These calculations are based on the particle distributions extracted from the solutions of the time-independent (steady-state approximation) system of integrodifferential transport equations mentioned above [29–32].Moreover, we extend the calculations of References [5, 14–16] so as to include contributions to neutrino emissivity (intensity) originating from the secondary muons (μ±) produced from the decay of charged pions (π±) which in our previous works had been ignored due to the long time consuming required for such simulations. Obviously, the mathematical problem of such a study becomes also more complicated (e.g., the number of coupled equations in the above-mentioned differential system increases rapidly in the case of considering neutrinos coming out of kaon (K±) decays) Lipari07.After fixing the model parameters and testing the derived algorithms on the reproducibility of some known properties of the well-studied SS 433 galactic microquasar [14–16], we perform detailed simulations for the galactic Cyg X-1 system as well as for the extragalactic M33 X-7 MQ [33]. The latter system has not been studied up to now from a jet emission viewpoint because this is a rather recently discovered (in 2006) X-ray binary system located in the neighboring galaxy of our Milky Way Galaxy, known as Messier 33 [33].The rest of the paper is organized as follows. In Section2, we describe briefly the main characteristics of hadronic jets in galactic microquasars (SS433 and Cyg X-1) as well as the extragalactic M33 X-7 system. Then, in Section 3, we present the formulation of our improved method related to the differential system of transport equations (integrodifferential system of coupled equations) and derive the appropriate algorithms. In Section 4, we present and discuss our results referred to high-energy neutrino and γ-ray energy spectra emitted from the extragalactic M33 X-7 MQ. Finally, in Section 5, we summarize the main conclusions extracted from the present investigation.


## 2. Brief Description of the Hadronic MQ Systems

In this section, we summarize briefly some basic properties of the hadronic MQ systems (like those mentioned before) and the characteristics of hadronic models that describe reliably microquasar jet emissions in the high-energy band of the electromagnetic spectrum (GeV-TeV gamma-rays) as well as high-energy neutrino emission [1, 4–6]. In the last few decades, from the investigation of a great number of MQs and X-ray emitting binary systems, stellar-mass black holes have been observed, with masses in the region of our interest determined mainly from the dynamics and structure properties of their companion stars [17, 18, 20].The galactic binary system SS 433, located 5.5 kpc from the Earth in Aquila constellation, displays two mildly relativistic jets (with bulk velocityυb≈0.26c) that are oppositely directed and precess in cones [34]. This system consists of a compact object (black hole) and an A-type companion (donor) star in coupled orbit with a period P~13.1 days [21]. For the masses of the component stars of this system, in this work, we adopted the observations of the INTEGRAL [22], i.e., MBH=9Mʘ and Mdon=30Mʘ for the black hole and the companion star, respectively. According to recent estimations for SS433, it suggests a low limit of the compact object of 5-9Mʘ placing it as a stellar-mass black hole [19].The second galactic system addressed in this work, Cygnus X-1, is an X-ray source in the constellation Cygnus which is located 1.86 kpc from the Earth [3]. The 19.2Mʘ O-type companion star [35] is in coupled orbit with the 14.8Mʘ compact object (stellar-mass black hole) [35], with orbital period P~5.6 days [36].As mentioned before, one of our main goals in this work is to calculate high-energyγ-ray and neutrino emissions from the binary system M33 X-7 [37, 38] which is located in the nearby galaxy Messier 33, the only known up to now black hole that is in an eclipsing binary system, with orbital period P~3.45 days [33]. The distance from Earth of this system is between ~840 kpc [39] and ~960 kpc [40]. In addition, the black hole mass is MBH=15.65±1.45Mʘ which orbits its Mdon=70.0±6.9Mʘ O-type companion star.In Table1, we tabulate in addition some other important parameters employed in this work for the systems under investigation, SS 433, Cygnus X-1, and M33 X-7.Table 1
Parameters of the model for M33 X-7, Cygnus X-1, and SS 433.
ParameterSymbolSS 433Cygnus X-1M33 X-7Black hole massMBH9.0Mʘ14.8Mʘ15.65MʘDistance from Earthd5.5 kpc1.86 kpc840-960 kpcDonor star massMdon30Mʘ19.2Mʘ70MʘDonor star type—A-typeO-typeO-typeOrbital periodP13.1 days5.6 days3.45 daysJet’s kinetic powerLk1039ergs−11038ergs−11038ergs−1Jet’s launching pointz01.3×109 cm108 cm1.9×108 cmBulk velocity of jet particlesυb0.26c0.6c0.8cJet’s bulk Lorentz factorΓb1.041.251.66Jet’s half-opening angleξ0.6∘1.5∘7∘Jet’s viewing angleθ78.05∘27.1∘74.6∘


### 2.1. Dynamics and Energetic Evolution of Hadronic MQ Jets

In the model considered in this work for the description of MQ jets, an accretion disk in the equatorial region of the compact object is present, and a fraction of the accreted material is expelled in two oppositely directed jets [4, 5, 41]. We assume an approximately conical mass outflow (jet) with a half-opening angle ξ (for the M33 X-7 MQ, we adopt the value ξ=7∘) and a radius given by rz=ztanξ (the coordinate z-axis coincides with the cone axis, jet’s direction). The injection point (plane) of the jet is at a distance z0 from the compact object. When z=z0, the radius of the jet is given by r0=z0tanξ (see Table 1).The initial jet radius isr0=rz0≈5Rsch, where Rsch=2GMBH/c2 (Rsch is the known Schwarzschild radius, the critical black hole radius for which the escape speed of particles becomes equal to the speed of light c), with G denoting Newton’s gravitational constant. Then, we find, for example, that the injection point is at z0=r0/tanξ≃1.9×108 cm, for the M33 X-7. For the sake of comparison, we mention that for Cygnus X-1, the injection point is at z0≃108 cm, while for SS 433, it is z0≃1.3×109 cm.It should be also noted that, for all systems studied in this work, the extension of the jet is determined through a maximum value ofz, denoted by zmax, which is assumed to be equal to zmax=5z0 (this point may be assumed to be at the boundaries of the jet with the ambient region).In discussing the energetic evolution and dynamics of MQ jets, the kinetic energy density of the jet,ρkz, is related to its kinetic luminosity, Lk, through the expression [1]
(1)ρkz=Lkπυbrz2,where υb is the bulk velocity of the jet particles (mostly protons). Within the context of the jet-accretion coupling hypothesis, only around 10% of the Eddington luminosity goes into the jet [42]. Here, we adopt Lk=1038ergs−1, for a black hole of mass ~15.7Mʘ, i.e., for the M33 X-7 system.Furthermore, assuming equipartition between the magnetic energy and the kinetic energy in the jet, it implies thatρmag=ρk [43, 44], and hence, the magnetic field that collimates the astrophysical jet plasma is given by [27]
(2)Bz≡8πρmagz=8πρkz.In general, for the kinetic power in the jet, many authors consider that a fraction is carried by primary protons and the rest is carried by electrons, i.e.,L=Lp+Le. Then, the relation between the proton and electron power is determined through a parameter α in such a way that Lp=αLe [6, 24, 27].In this article, however, we adopt the case ofα=100 which means that we consider a proton-dominated jet. Other quantities needed for the purposes of this paper are listed in Table 1.


### 2.2. p-p Collision Mechanism inside MQ Jets

The collision of relativistic protons with the cold ones inside the jet (p-p collision mechanism) produces high-energy charged particles (pionsπ±, kaons K±, muons μ±, etc.) and neutral particles (π0, K0, K~0, η particles, etc.). Important primary reactions of this type are
(3)p+p⟶p+p+aπ0+bπ++π−,p+p⟶p+n+π++aπ0+bπ++π−,p+p⟶n+n+2π++aπ0+bπ++π−,where a and b denote the pion multiplicities [1], and similarly for kaon production [45].Charged pionsπ± (and kaons K±), afterwards, decay to charged leptons (muons μ±, and electrons or positrons, e±) as well as neutrinos and antineutrinos as
(4)π+⟶μ++νμ,π+⟶e++νe,π−⟶μ−+ν¯μ,π−⟶e−+ν¯e.Furthermore, muons also decay giving neutrinos and electrons (or positrons) as(5)μ+⟶e++νe+ν¯μ,μ−⟶e−+ν¯e+νμ.We mention that neutral pions (π0), η-particles, etc., decay producing γ-rays according to the reactions
(6)π0⟶γ+γ,π0⟶γ+e−+e+,(7)η⟶γ+γ,η⟶π0+e−+e+.In our present study, we consider neutrinos generated from both types of charged-pion decays (for neutrinos coming from kaon decays, the reader is referred, e.g., to Reference [45] and also from both types of charged muon decays [25, 26]. Moreover, we consider γ-rays produced from the π0 decays of Equation (6). Thus, high-energy neutrinos are inevitably accompanied by pionic gamma-rays, a phenomenon well known as multimessenger emission from microquasar jets.


### 2.3. Accelerating and Cooling Rates of the Jet Processes

In general, the primary charged particles (p and e−) gain energy during moving within the magnetic field B (Fermi acceleration). The acceleration rate of the (initially cold) protons to an energy E=Ep, known as shock acceleration, is defined as tacc−1=E−1dE/dt and is given by the relation [1]
(8)tacc−1≡1EdEdt≈ηceBEp.η=0.1 is the acceleration efficiency which means that only 10% of the cold/thermal protons may acquire relativistic energies. This is equivalent to existence of an efficient accelerator at the base of the jets where shocks are rather relativistic [46].In the p-p collision mechanism, the cross section and the rate of inelastic p-p scattering (relativistic protons scatter on cold ones),σppinel and tpp−1, respectively, play a crucial role (the inelasticity coefficient is Kpp≈1/2). The corresponding cross section, σppinel, as a function of the fast proton energy Ep, reads [47, 48]
(9)σppinelEp=0.25L2+1.88L+34.31−EthEp42×10−27cm2,where L=lnEp/1000GeV and Eth denotes the minimum energy of (thermal) protons, Epmin, which is equal to Epmin≡Eth=1.22 GeV, see e.g., Appendix of Reference [5].The corresponding rate,tpp−1, of the inelastic p-p scattering is given in Table 2 in terms of the cross section σppinel and the number density of cold particles (jet protons), npz, at a distance z from the black hole which may be written as [1]
(10)npz=1−qrelΓbmpc2ρkz.Table 2
Accelerating, cooling, and decay rates for particles moving inside hadronic MQ jets.
Rate parameterRate symbolBasic rate definitionProton accelerating ratestacc−1ηceB/EpProton-proton collision ratetpp−1npzσppinelEpKppSynchrotron radiation ratetsyn−14/3me/m3γσTB2/mec8πAdiabatic expansion ratetad−12/3υb/zPion-proton collisiontπp−10.5cnpzσπpinelEpProton knock out (escape) ratetp−1=tesc,p−1tesc,p−1=c/zmax−zPion knock out ratetπ−1tesc,π−1+tdec,π−1Muon knock out ratetμ−1tesc,μ−1+tdec,μ−1qrel=0.1 denotes the portion of the relativistic protons.In Table2, in addition to the p-p collision rate and accelerating rates of protons, we tabulate the most significant cooling rates (synchrotron, adiabatic, etc.) for protons, pions, and muons as well.In more detail, the well-known synchrotron radiation is emitted from charged particles of energyE=γmc2, with γ being the particle’s Lorentz factor and m its mass, moving inside the magnetized plasma. Thus, the rate of the synchrotron radiation, tsyn−1, is very important in MQ systems with strong magnetic fields (see Table 2).Furthermore, due to the adiabatic expansion of the jet plasma, all comoving particles inside it lose energy with a ratetad−1 [44]. This rate, known as the adiabatic cooling rate, depends on the particle velocities υbtanξ.In addition,tπp−1 gives the rate of pion-proton (π−p) inelastic scattering inside the jet. The corresponding σπpinel cross section is related to the inelastic p-p scattering cross section with the relationship σπpE=2/3σppinelE, where the factor 2/3 comes out of the fact that protons are made of three valence quarks, while the pions of only two quarks [49].It is worth mentioning that some processes taking place inside the jet lead to particles’ knock out, either by escaping from the jet (and/or from the energy range of interest) or by their decay processes [50]. Thus, the corresponding rates are (i) the escaping rate, tesc−1, given by
(11)tesc−1≈czmax−z.zmax−z0 represents the length of the acceleration zone, and (ii) the decay rate of the particle in question, tdec−1, which is known from direct measurements of particle’s life time.The rates of the knock out (or catastrophic) processes inside the jet for pions,tπ−1, and muons, tμ−1, include contributions from the decay process and escaping process as
(12)tπ−1E,z=tesc,π−1z+tdec,π−1E,(13)tμ−1E,z=tesc,μ−1z+tdec,μ−1E,where tdec,π−1=2.6×10−8γπ−1s−1, for pions, and tdec,μ−1=2.2×10−6γμ−1s−1, for muons. We note that Equation (12) holds also for protons due to the fact that tdec,p−1=0 (or equivalently the proton life time is infinite), so as tp−1E,z=tesc,p−1z, see Table 2.The rates discussed above enter the set of basic coupled transport equations (system of kinematic equations) that describe the particle distributions as we will discuss in Section3.


#### 2.3.1. Energy Variation of the Various Rates of Jet Particles

In Figures1–3, the variation of the aforementioned main rates, versus the particle energy E, for protons, pions, and muons, respectively, and for the SS 433 and M33 X-7 MQ systems, is illustrated. In other words, in these figures, the region of dominance of the accelerating, cooling, etc., rates throughout the energy range of interest, Eth≡Emin≤E≤Emax, is demonstrated.Figure 1
Cooling rates for protons in the jets of M33 X-7 (a) and SS 433 (b) at the base of the jetsz0. The plots of cooling rates show the synchrotron emission (solid lines), the adiabatic cooling (dotted lines), and the inelastic p-p collision rate (dashed lines).
(a)(b)Figure 2
Cooling rates for pions in the jets of M33 X-7 (a) and SS 433 (b) at the base of the jets.
(a)(b)Figure 3
Cooling rates for muons in the jets of M33 X-7 (a) and SS 433 (b) at the base of the jets.
(a)(b)More specifically, Figure1 indicates the cooling rates at the base of the jet for protons, Figure 2 shows those for pions and Figure 3 shows the corresponding ones for muons, in the extragalactic binary system M33 X-7 (top of these figures). For comparison, the corresponding results for the galactic binary system SS 433 (bottom of these figures) are also shown.The plots in Figure2 for pions, show, respectively, the variation versus the energy E of the rate for synchrotron emission (solid lines), for the pion-proton inelastic collision (dashed lines), for the adiabatic cooling (dotted lines), and for the decay of pions (dot-dashed lines). They refer to the M33 X-7 (top) and the SS 433 (bottom) microquasar systems. Also, Figure 3 illustrates the cooling rates of muons for the synchrotron emission (solid lines), for the adiabatic expansion (dotted lines), and for the decay of muons (dot-dashed lines). They are similar to those for pions (Figure 2) except for the corresponding muon-proton and muon-pion scattering channels which are ignored in this work.From Figures2 and 3, it becomes obvious that the particle synchrotron losses dominate the high-energy region. In the case of protons (Figure 1), however, due to the large proton mass, synchrotron losses are not dominant up to very high energies (~107 GeV) for the MQ systems studied. On the other hand, for pions and muons, the decay losses dominate for lower energies, due to their short life time (see Table 2).The main difference between these two MQ systems (SS 433 and M33 X-7) concerns the synchrotron cooling rates. In general, according to Equations (1) and (2), for systems with wide half-opening angle (e.g., M33 X-7, with ξ=7∘), the magnetic energy density is lower than that of systems with small ξ (e.g., the SS 433 with ξ=0.6∘). Hence, the magnetic field for M33 X-7 is lower which becomes obvious by comparing Equation (2) for the two systems. The latter conclusion justifies the lower synchrotron loss rate for the M33 X-7 system.


## 2.1. Dynamics and Energetic Evolution of Hadronic MQ Jets

In the model considered in this work for the description of MQ jets, an accretion disk in the equatorial region of the compact object is present, and a fraction of the accreted material is expelled in two oppositely directed jets [4, 5, 41]. We assume an approximately conical mass outflow (jet) with a half-opening angle ξ (for the M33 X-7 MQ, we adopt the value ξ=7∘) and a radius given by rz=ztanξ (the coordinate z-axis coincides with the cone axis, jet’s direction). The injection point (plane) of the jet is at a distance z0 from the compact object. When z=z0, the radius of the jet is given by r0=z0tanξ (see Table 1).The initial jet radius isr0=rz0≈5Rsch, where Rsch=2GMBH/c2 (Rsch is the known Schwarzschild radius, the critical black hole radius for which the escape speed of particles becomes equal to the speed of light c), with G denoting Newton’s gravitational constant. Then, we find, for example, that the injection point is at z0=r0/tanξ≃1.9×108 cm, for the M33 X-7. For the sake of comparison, we mention that for Cygnus X-1, the injection point is at z0≃108 cm, while for SS 433, it is z0≃1.3×109 cm.It should be also noted that, for all systems studied in this work, the extension of the jet is determined through a maximum value ofz, denoted by zmax, which is assumed to be equal to zmax=5z0 (this point may be assumed to be at the boundaries of the jet with the ambient region).In discussing the energetic evolution and dynamics of MQ jets, the kinetic energy density of the jet,ρkz, is related to its kinetic luminosity, Lk, through the expression [1]
(1)ρkz=Lkπυbrz2,where υb is the bulk velocity of the jet particles (mostly protons). Within the context of the jet-accretion coupling hypothesis, only around 10% of the Eddington luminosity goes into the jet [42]. Here, we adopt Lk=1038ergs−1, for a black hole of mass ~15.7Mʘ, i.e., for the M33 X-7 system.Furthermore, assuming equipartition between the magnetic energy and the kinetic energy in the jet, it implies thatρmag=ρk [43, 44], and hence, the magnetic field that collimates the astrophysical jet plasma is given by [27]
(2)Bz≡8πρmagz=8πρkz.In general, for the kinetic power in the jet, many authors consider that a fraction is carried by primary protons and the rest is carried by electrons, i.e.,L=Lp+Le. Then, the relation between the proton and electron power is determined through a parameter α in such a way that Lp=αLe [6, 24, 27].In this article, however, we adopt the case ofα=100 which means that we consider a proton-dominated jet. Other quantities needed for the purposes of this paper are listed in Table 1.


## 2.2. p-p Collision Mechanism inside MQ Jets

The collision of relativistic protons with the cold ones inside the jet (p-p collision mechanism) produces high-energy charged particles (pionsπ±, kaons K±, muons μ±, etc.) and neutral particles (π0, K0, K~0, η particles, etc.). Important primary reactions of this type are
(3)p+p⟶p+p+aπ0+bπ++π−,p+p⟶p+n+π++aπ0+bπ++π−,p+p⟶n+n+2π++aπ0+bπ++π−,where a and b denote the pion multiplicities [1], and similarly for kaon production [45].Charged pionsπ± (and kaons K±), afterwards, decay to charged leptons (muons μ±, and electrons or positrons, e±) as well as neutrinos and antineutrinos as
(4)π+⟶μ++νμ,π+⟶e++νe,π−⟶μ−+ν¯μ,π−⟶e−+ν¯e.Furthermore, muons also decay giving neutrinos and electrons (or positrons) as(5)μ+⟶e++νe+ν¯μ,μ−⟶e−+ν¯e+νμ.We mention that neutral pions (π0), η-particles, etc., decay producing γ-rays according to the reactions
(6)π0⟶γ+γ,π0⟶γ+e−+e+,(7)η⟶γ+γ,η⟶π0+e−+e+.In our present study, we consider neutrinos generated from both types of charged-pion decays (for neutrinos coming from kaon decays, the reader is referred, e.g., to Reference [45] and also from both types of charged muon decays [25, 26]. Moreover, we consider γ-rays produced from the π0 decays of Equation (6). Thus, high-energy neutrinos are inevitably accompanied by pionic gamma-rays, a phenomenon well known as multimessenger emission from microquasar jets.


## 2.3. Accelerating and Cooling Rates of the Jet Processes

In general, the primary charged particles (p and e−) gain energy during moving within the magnetic field B (Fermi acceleration). The acceleration rate of the (initially cold) protons to an energy E=Ep, known as shock acceleration, is defined as tacc−1=E−1dE/dt and is given by the relation [1]
(8)tacc−1≡1EdEdt≈ηceBEp.η=0.1 is the acceleration efficiency which means that only 10% of the cold/thermal protons may acquire relativistic energies. This is equivalent to existence of an efficient accelerator at the base of the jets where shocks are rather relativistic [46].In the p-p collision mechanism, the cross section and the rate of inelastic p-p scattering (relativistic protons scatter on cold ones),σppinel and tpp−1, respectively, play a crucial role (the inelasticity coefficient is Kpp≈1/2). The corresponding cross section, σppinel, as a function of the fast proton energy Ep, reads [47, 48]
(9)σppinelEp=0.25L2+1.88L+34.31−EthEp42×10−27cm2,where L=lnEp/1000GeV and Eth denotes the minimum energy of (thermal) protons, Epmin, which is equal to Epmin≡Eth=1.22 GeV, see e.g., Appendix of Reference [5].The corresponding rate,tpp−1, of the inelastic p-p scattering is given in Table 2 in terms of the cross section σppinel and the number density of cold particles (jet protons), npz, at a distance z from the black hole which may be written as [1]
(10)npz=1−qrelΓbmpc2ρkz.Table 2
Accelerating, cooling, and decay rates for particles moving inside hadronic MQ jets.
Rate parameterRate symbolBasic rate definitionProton accelerating ratestacc−1ηceB/EpProton-proton collision ratetpp−1npzσppinelEpKppSynchrotron radiation ratetsyn−14/3me/m3γσTB2/mec8πAdiabatic expansion ratetad−12/3υb/zPion-proton collisiontπp−10.5cnpzσπpinelEpProton knock out (escape) ratetp−1=tesc,p−1tesc,p−1=c/zmax−zPion knock out ratetπ−1tesc,π−1+tdec,π−1Muon knock out ratetμ−1tesc,μ−1+tdec,μ−1qrel=0.1 denotes the portion of the relativistic protons.In Table2, in addition to the p-p collision rate and accelerating rates of protons, we tabulate the most significant cooling rates (synchrotron, adiabatic, etc.) for protons, pions, and muons as well.In more detail, the well-known synchrotron radiation is emitted from charged particles of energyE=γmc2, with γ being the particle’s Lorentz factor and m its mass, moving inside the magnetized plasma. Thus, the rate of the synchrotron radiation, tsyn−1, is very important in MQ systems with strong magnetic fields (see Table 2).Furthermore, due to the adiabatic expansion of the jet plasma, all comoving particles inside it lose energy with a ratetad−1 [44]. This rate, known as the adiabatic cooling rate, depends on the particle velocities υbtanξ.In addition,tπp−1 gives the rate of pion-proton (π−p) inelastic scattering inside the jet. The corresponding σπpinel cross section is related to the inelastic p-p scattering cross section with the relationship σπpE=2/3σppinelE, where the factor 2/3 comes out of the fact that protons are made of three valence quarks, while the pions of only two quarks [49].It is worth mentioning that some processes taking place inside the jet lead to particles’ knock out, either by escaping from the jet (and/or from the energy range of interest) or by their decay processes [50]. Thus, the corresponding rates are (i) the escaping rate, tesc−1, given by
(11)tesc−1≈czmax−z.zmax−z0 represents the length of the acceleration zone, and (ii) the decay rate of the particle in question, tdec−1, which is known from direct measurements of particle’s life time.The rates of the knock out (or catastrophic) processes inside the jet for pions,tπ−1, and muons, tμ−1, include contributions from the decay process and escaping process as
(12)tπ−1E,z=tesc,π−1z+tdec,π−1E,(13)tμ−1E,z=tesc,μ−1z+tdec,μ−1E,where tdec,π−1=2.6×10−8γπ−1s−1, for pions, and tdec,μ−1=2.2×10−6γμ−1s−1, for muons. We note that Equation (12) holds also for protons due to the fact that tdec,p−1=0 (or equivalently the proton life time is infinite), so as tp−1E,z=tesc,p−1z, see Table 2.The rates discussed above enter the set of basic coupled transport equations (system of kinematic equations) that describe the particle distributions as we will discuss in Section3.


### 2.3.1. Energy Variation of the Various Rates of Jet Particles

In Figures1–3, the variation of the aforementioned main rates, versus the particle energy E, for protons, pions, and muons, respectively, and for the SS 433 and M33 X-7 MQ systems, is illustrated. In other words, in these figures, the region of dominance of the accelerating, cooling, etc., rates throughout the energy range of interest, Eth≡Emin≤E≤Emax, is demonstrated.Figure 1
Cooling rates for protons in the jets of M33 X-7 (a) and SS 433 (b) at the base of the jetsz0. The plots of cooling rates show the synchrotron emission (solid lines), the adiabatic cooling (dotted lines), and the inelastic p-p collision rate (dashed lines).
(a)(b)Figure 2
Cooling rates for pions in the jets of M33 X-7 (a) and SS 433 (b) at the base of the jets.
(a)(b)Figure 3
Cooling rates for muons in the jets of M33 X-7 (a) and SS 433 (b) at the base of the jets.
(a)(b)More specifically, Figure1 indicates the cooling rates at the base of the jet for protons, Figure 2 shows those for pions and Figure 3 shows the corresponding ones for muons, in the extragalactic binary system M33 X-7 (top of these figures). For comparison, the corresponding results for the galactic binary system SS 433 (bottom of these figures) are also shown.The plots in Figure2 for pions, show, respectively, the variation versus the energy E of the rate for synchrotron emission (solid lines), for the pion-proton inelastic collision (dashed lines), for the adiabatic cooling (dotted lines), and for the decay of pions (dot-dashed lines). They refer to the M33 X-7 (top) and the SS 433 (bottom) microquasar systems. Also, Figure 3 illustrates the cooling rates of muons for the synchrotron emission (solid lines), for the adiabatic expansion (dotted lines), and for the decay of muons (dot-dashed lines). They are similar to those for pions (Figure 2) except for the corresponding muon-proton and muon-pion scattering channels which are ignored in this work.From Figures2 and 3, it becomes obvious that the particle synchrotron losses dominate the high-energy region. In the case of protons (Figure 1), however, due to the large proton mass, synchrotron losses are not dominant up to very high energies (~107 GeV) for the MQ systems studied. On the other hand, for pions and muons, the decay losses dominate for lower energies, due to their short life time (see Table 2).The main difference between these two MQ systems (SS 433 and M33 X-7) concerns the synchrotron cooling rates. In general, according to Equations (1) and (2), for systems with wide half-opening angle (e.g., M33 X-7, with ξ=7∘), the magnetic energy density is lower than that of systems with small ξ (e.g., the SS 433 with ξ=0.6∘). Hence, the magnetic field for M33 X-7 is lower which becomes obvious by comparing Equation (2) for the two systems. The latter conclusion justifies the lower synchrotron loss rate for the M33 X-7 system.


## 2.3.1. Energy Variation of the Various Rates of Jet Particles

In Figures1–3, the variation of the aforementioned main rates, versus the particle energy E, for protons, pions, and muons, respectively, and for the SS 433 and M33 X-7 MQ systems, is illustrated. In other words, in these figures, the region of dominance of the accelerating, cooling, etc., rates throughout the energy range of interest, Eth≡Emin≤E≤Emax, is demonstrated.Figure 1
Cooling rates for protons in the jets of M33 X-7 (a) and SS 433 (b) at the base of the jetsz0. The plots of cooling rates show the synchrotron emission (solid lines), the adiabatic cooling (dotted lines), and the inelastic p-p collision rate (dashed lines).
(a)(b)Figure 2
Cooling rates for pions in the jets of M33 X-7 (a) and SS 433 (b) at the base of the jets.
(a)(b)Figure 3
Cooling rates for muons in the jets of M33 X-7 (a) and SS 433 (b) at the base of the jets.
(a)(b)More specifically, Figure1 indicates the cooling rates at the base of the jet for protons, Figure 2 shows those for pions and Figure 3 shows the corresponding ones for muons, in the extragalactic binary system M33 X-7 (top of these figures). For comparison, the corresponding results for the galactic binary system SS 433 (bottom of these figures) are also shown.The plots in Figure2 for pions, show, respectively, the variation versus the energy E of the rate for synchrotron emission (solid lines), for the pion-proton inelastic collision (dashed lines), for the adiabatic cooling (dotted lines), and for the decay of pions (dot-dashed lines). They refer to the M33 X-7 (top) and the SS 433 (bottom) microquasar systems. Also, Figure 3 illustrates the cooling rates of muons for the synchrotron emission (solid lines), for the adiabatic expansion (dotted lines), and for the decay of muons (dot-dashed lines). They are similar to those for pions (Figure 2) except for the corresponding muon-proton and muon-pion scattering channels which are ignored in this work.From Figures2 and 3, it becomes obvious that the particle synchrotron losses dominate the high-energy region. In the case of protons (Figure 1), however, due to the large proton mass, synchrotron losses are not dominant up to very high energies (~107 GeV) for the MQ systems studied. On the other hand, for pions and muons, the decay losses dominate for lower energies, due to their short life time (see Table 2).The main difference between these two MQ systems (SS 433 and M33 X-7) concerns the synchrotron cooling rates. In general, according to Equations (1) and (2), for systems with wide half-opening angle (e.g., M33 X-7, with ξ=7∘), the magnetic energy density is lower than that of systems with small ξ (e.g., the SS 433 with ξ=0.6∘). Hence, the magnetic field for M33 X-7 is lower which becomes obvious by comparing Equation (2) for the two systems. The latter conclusion justifies the lower synchrotron loss rate for the M33 X-7 system.


## 3. Formulation of the Key-Role Transport Equations inside Jet Plasmas

In this section, we make a first step towards improving the mathematical formulation describing the creation, scattering, decay, diffusion, and emission of particles (including neutrinos) and electromagnetic radiation in astrophysical outflows (jets). The basic key-role tool of this formalism is the general transport equation which satisfies each of the particles considered [6].We should note that such a formalism may be reliably applicable for multispecies (multi-particle) emissions from MQs and X-ray binaries, as well as from central regions of galaxies like the AGN systems, in which the compact object is a massive or supermassive black hole, see, e.g., [1, 50]. The geometry of the latter sources is assumed to be rather spherical, while in the first class of systems that are studied in this work, the sources are assumed conical and the particles (neutrinos) or radiation are emitted to the forward or backward direction.From a mathematical point of view, the transport equation which satisfies any kind of particles moving into the (dark) jets of MQs is an integrodifferential equation, as is explained the following.


### 3.1. The Transport Equation for Particles Moving inside MQ Jets

The general transport equation describes the concentration (distribution) of particle ofj-kind, NjE,r,t, where j=p,e±,π±,μ±, etc., as a function of the time t, the particle’s energy E, and the position r inside the (conical) jet. In essence, this is a phenomenological macroscopic equation of particle (or radiation) transport describing astrophysical outflows, and it is written as
(14)∂Nj∂t−∇·Dj∇Nj+∂bjNj∂E−12∂2djNj∂E2=QjE,r,t−pjNj+∑k∫PjkE′,ENjE,r,tdE,j=p,π±,μ±⋯,where the parameters Dj, bj, dj, pj, and Pjk may depend upon the space and time coordinates and also on the energy E. The latter equation coincides with the continuity equation for particles of type-j, with j=1,2,3 and 1,2,3≡p,π,μ.The termQjE,r,t in the r.h.s. of Equation (14) is equal to the intensity of the source producing the particles-j, which is also known as the injection function of particles-j. This means that QjE,r,tdEd3rdt represents the number of particles kind-j provided by the sources in a volume element d3r, in the energy range between E and E+dE during the time dt. In the case when the j-type particles are products of a chain reaction (as it holds in our present work assuming the p-p reaction chain described in Section 2.2), the function QjE,r,t couples the j-reaction with its parent reaction, i.e., the equation of particles kind-(j−1).Further, the term proportional topj exists in cases when catastrophic (or knock out) processes take place that cause catastrophic energy losses. Then, this gives the probability per unit time for the losses in question to occur. Thus, the coefficient pj is written as
(15)pj=1Tj≡tj−1,with Tj is the mean life time of the particles of kind-j. Also, pjNj denotes the number of particles “knocked out” per unit time.The coefficientsbj in the third term (r.h.s.) are equal to the mean energy increment of the particle-j per unit time, i.e.,
(16)bj=dEdt.For the special case whenbj=bjE, this coefficient is related to energy losses of various cooling processes. Also, the coefficients PjkE, with k>j, in the latter summation of the transport equation, are related to the existence of fragmentation of the primary (or secondary) particles which in this work are assumed as nonexisting.The parameterDj is known as the diffusion coefficient which, in general, is a function of the coordinates r and time t if the concentration of the jet plasma and its macroscopic motion are inhomogeneous in the volume of the jet. By assuming that Dj=0, we may reliably describe the regular motion of particles j-kind along the lines of force of the magnetic field B.For simplicity, in this work, we make the realistic assumptions that the above coefficients depend only on the particle’s energyE, i.e., we assume that DjE, bjE, and pjE. Then, the system of general transport equations (14) reads
(17)∂∂tNpE,z,tNπE,z,tNμE,z,t−∇2DpNpE,z,tDπNπE,z,tDμNμE,z,t+∂∂EbpNpE,z,tbπNπE,z,tbμNμE,z,t−12∂2∂E2dpNpE,z,tdπNπE,z,tdμNμE,z,t+tesc−1NpE,z,ttπ−1NπE,z,ttμ−1NμE,z,t=QpE,z,tQπE,z,tQμE,z,t.Note that, in the latter system of coupled differential equations, we write down only three particles (protons, pions, and muons), without considering the complete reaction family tree, since we have not distinguished particles with different charge asπ± and μ±, and we have not considered the possibility of left-right symmetry, μL,R, of muons as we have done in Section 4. This means that, by considering all these particles, the system of Equation (17) will have seven lines.Usually, we ignore the term involving second derivative with respect toE means that the acceleration term is negligible, an approximation adopted even by recent leptohadronic steady-state models [31].The solution of the general transport equation for particles of one kind (when the last term of this equation can be omitted) is described in detail in Reference [16]. In the latter work, however, neutrinos produced only from π± have been considered. In the present paper, we proceed further and include neutrino emissions also from the muon decays, so we need to solve the system of Equation (17) by including the third line too.The complete form of the system of Equation (17) is treated with the method of Reference [28] in order to find the exact solutions. In Reference [26], Equation (17) is treated semianalytically. In the present work, we restrict ourselves to simplified forms resulting by neglecting various phenomena (processes) taking place inside the astrophysical outflows (jet plasma), even though some of the assumed omissions may be considered as rather crud approximations.In the following, we discuss the cases of Equation (17) satisfying the conditions of the one-zone approximation [16, 51], i.e., the cases when the particle distributions are independent of time (steady-state approximation).


#### 3.1.1. Transport Equation Assuming Absence of Energy Losses

At first, in the calculations of our present work, we start by writing the simplest solutions of the system of transport equations obtained by assuming that all energy losses are absent, i.e.,bj=0. Then, by denoting Np,0, Nπ,0, and Nμ,0 the corresponding energy distributions for protons, pions, and muons, respectively, the system of transport equations takes the trivial form
(18)tp−1ENp,0E,ztπ−1ENπ,0E,ztμ−1ENμ,0E,z=QpE,zQπE,zQμE,z.The calculated distributionsNj,0, with j=p, π±, and μ± of the latter equations are discussed in the next section.


#### 3.1.2. Steady-State Transport Equations with Particle Losses when Moving inside the Jet Plasma

Under the conditions of the steady-state approximation, and assuming that the various energy losses are absent, the system of transport equations (Equation (17)) takes the form
(19)∂∂EbpNpE,zbπNπE,zbμNμE,z+tesc−1NpE,ztπ−1NπE,ztμ−1NμE,z=QpE,zQπE,zQμE,z.In order to calculate the time-independent neutrino and gamma-ray emissivities, we need, first, to calculate the distributions of protons, pions, and muons,Np, Nπ, and Nμ, respectively, from Equation (19). For these computations, we used a code written in the C programming language, mainly following the assumptions of References [6, 24, 27].We note that, in the l.h.s. of Equation (19), the rates of the knock out processes inside the jet, tj−1, for the particles j=p,π,μ, are provided by Equation (12) and the values given below that equation. Moreover, the coefficients bj entering Equation (19) are obtained from expressions of the form bjE=−Etloss,j−1E by inserting the corresponding energy loss rates tloss,j−1 for each particle as
(20)bpE=−Etloss,p−1E=−Etsyn−1+tad−1+tpp−1,bπE=−Etloss,π−1E=−Etsyn−1+tad−1+tπp−1,bμE=−Etloss,μ−1E=−Etsyn−1+tad−1.It is worth noting that, in the latter equations, some additional energy loss terms may appear in cases when the rates of some other processes, ignored in our present work (like, e.g., the inverse Compton scattering and the pion-muon scattering), may be considered important for the description of neutrino andγ-ray emissions from other cosmic structures (AGNs, etc.).


### 3.2. Brief Description of the Model Assumptions

Before embarking to concrete calculations, we find it helpful for the reader to briefly summarize the simplifications made and the main approximations assumed in our method. At first, in purely hadronic models, only the protons from the primary charged particles (p,e+,e−) inside the jet are considered to contribute significantly in the emission of high-energy γ-rays and neutrinos, see References [4, 5] and discussion after Equation (2).As mentioned in Introduction, in the present calculations, we consider that the (GeV-TeV energy range of our interest) neutrinos andγ-rays are due to the p-p collision mechanism taking place inside the microquasar jets. The maximum nonthermal (primary) proton energy employed in our numerical integrations in the following (Section 3) is taken equal to Epmax≈106−107 GeV. We should stress that the p-p mechanism is assumed to dominate over others like, e.g., scattering of jet protons with the protons of the donor-star wind or the proton-photon scattering inside the jet (case of γ-ray production) [1].We, furthermore, note that in this work, we ignore neutrinos coming from the decay of the secondary kaons (K±) and γ-rays produced from η particles, neutral kaons, etc. As discussed in Section 2.2, the contributions originated from charged pion decays (π±) and the subsequent charged muon decays (μ±), for neutrino production, and those coming from neutral pion decays (π0), for γ-ray production, are assumed to be the main chain-reaction channels.


## 3.1. The Transport Equation for Particles Moving inside MQ Jets

The general transport equation describes the concentration (distribution) of particle ofj-kind, NjE,r,t, where j=p,e±,π±,μ±, etc., as a function of the time t, the particle’s energy E, and the position r inside the (conical) jet. In essence, this is a phenomenological macroscopic equation of particle (or radiation) transport describing astrophysical outflows, and it is written as
(14)∂Nj∂t−∇·Dj∇Nj+∂bjNj∂E−12∂2djNj∂E2=QjE,r,t−pjNj+∑k∫PjkE′,ENjE,r,tdE,j=p,π±,μ±⋯,where the parameters Dj, bj, dj, pj, and Pjk may depend upon the space and time coordinates and also on the energy E. The latter equation coincides with the continuity equation for particles of type-j, with j=1,2,3 and 1,2,3≡p,π,μ.The termQjE,r,t in the r.h.s. of Equation (14) is equal to the intensity of the source producing the particles-j, which is also known as the injection function of particles-j. This means that QjE,r,tdEd3rdt represents the number of particles kind-j provided by the sources in a volume element d3r, in the energy range between E and E+dE during the time dt. In the case when the j-type particles are products of a chain reaction (as it holds in our present work assuming the p-p reaction chain described in Section 2.2), the function QjE,r,t couples the j-reaction with its parent reaction, i.e., the equation of particles kind-(j−1).Further, the term proportional topj exists in cases when catastrophic (or knock out) processes take place that cause catastrophic energy losses. Then, this gives the probability per unit time for the losses in question to occur. Thus, the coefficient pj is written as
(15)pj=1Tj≡tj−1,with Tj is the mean life time of the particles of kind-j. Also, pjNj denotes the number of particles “knocked out” per unit time.The coefficientsbj in the third term (r.h.s.) are equal to the mean energy increment of the particle-j per unit time, i.e.,
(16)bj=dEdt.For the special case whenbj=bjE, this coefficient is related to energy losses of various cooling processes. Also, the coefficients PjkE, with k>j, in the latter summation of the transport equation, are related to the existence of fragmentation of the primary (or secondary) particles which in this work are assumed as nonexisting.The parameterDj is known as the diffusion coefficient which, in general, is a function of the coordinates r and time t if the concentration of the jet plasma and its macroscopic motion are inhomogeneous in the volume of the jet. By assuming that Dj=0, we may reliably describe the regular motion of particles j-kind along the lines of force of the magnetic field B.For simplicity, in this work, we make the realistic assumptions that the above coefficients depend only on the particle’s energyE, i.e., we assume that DjE, bjE, and pjE. Then, the system of general transport equations (14) reads
(17)∂∂tNpE,z,tNπE,z,tNμE,z,t−∇2DpNpE,z,tDπNπE,z,tDμNμE,z,t+∂∂EbpNpE,z,tbπNπE,z,tbμNμE,z,t−12∂2∂E2dpNpE,z,tdπNπE,z,tdμNμE,z,t+tesc−1NpE,z,ttπ−1NπE,z,ttμ−1NμE,z,t=QpE,z,tQπE,z,tQμE,z,t.Note that, in the latter system of coupled differential equations, we write down only three particles (protons, pions, and muons), without considering the complete reaction family tree, since we have not distinguished particles with different charge asπ± and μ±, and we have not considered the possibility of left-right symmetry, μL,R, of muons as we have done in Section 4. This means that, by considering all these particles, the system of Equation (17) will have seven lines.Usually, we ignore the term involving second derivative with respect toE means that the acceleration term is negligible, an approximation adopted even by recent leptohadronic steady-state models [31].The solution of the general transport equation for particles of one kind (when the last term of this equation can be omitted) is described in detail in Reference [16]. In the latter work, however, neutrinos produced only from π± have been considered. In the present paper, we proceed further and include neutrino emissions also from the muon decays, so we need to solve the system of Equation (17) by including the third line too.The complete form of the system of Equation (17) is treated with the method of Reference [28] in order to find the exact solutions. In Reference [26], Equation (17) is treated semianalytically. In the present work, we restrict ourselves to simplified forms resulting by neglecting various phenomena (processes) taking place inside the astrophysical outflows (jet plasma), even though some of the assumed omissions may be considered as rather crud approximations.In the following, we discuss the cases of Equation (17) satisfying the conditions of the one-zone approximation [16, 51], i.e., the cases when the particle distributions are independent of time (steady-state approximation).


### 3.1.1. Transport Equation Assuming Absence of Energy Losses

At first, in the calculations of our present work, we start by writing the simplest solutions of the system of transport equations obtained by assuming that all energy losses are absent, i.e.,bj=0. Then, by denoting Np,0, Nπ,0, and Nμ,0 the corresponding energy distributions for protons, pions, and muons, respectively, the system of transport equations takes the trivial form
(18)tp−1ENp,0E,ztπ−1ENπ,0E,ztμ−1ENμ,0E,z=QpE,zQπE,zQμE,z.The calculated distributionsNj,0, with j=p, π±, and μ± of the latter equations are discussed in the next section.


### 3.1.2. Steady-State Transport Equations with Particle Losses when Moving inside the Jet Plasma

Under the conditions of the steady-state approximation, and assuming that the various energy losses are absent, the system of transport equations (Equation (17)) takes the form
(19)∂∂EbpNpE,zbπNπE,zbμNμE,z+tesc−1NpE,ztπ−1NπE,ztμ−1NμE,z=QpE,zQπE,zQμE,z.In order to calculate the time-independent neutrino and gamma-ray emissivities, we need, first, to calculate the distributions of protons, pions, and muons,Np, Nπ, and Nμ, respectively, from Equation (19). For these computations, we used a code written in the C programming language, mainly following the assumptions of References [6, 24, 27].We note that, in the l.h.s. of Equation (19), the rates of the knock out processes inside the jet, tj−1, for the particles j=p,π,μ, are provided by Equation (12) and the values given below that equation. Moreover, the coefficients bj entering Equation (19) are obtained from expressions of the form bjE=−Etloss,j−1E by inserting the corresponding energy loss rates tloss,j−1 for each particle as
(20)bpE=−Etloss,p−1E=−Etsyn−1+tad−1+tpp−1,bπE=−Etloss,π−1E=−Etsyn−1+tad−1+tπp−1,bμE=−Etloss,μ−1E=−Etsyn−1+tad−1.It is worth noting that, in the latter equations, some additional energy loss terms may appear in cases when the rates of some other processes, ignored in our present work (like, e.g., the inverse Compton scattering and the pion-muon scattering), may be considered important for the description of neutrino andγ-ray emissions from other cosmic structures (AGNs, etc.).


## 3.1.1. Transport Equation Assuming Absence of Energy Losses

At first, in the calculations of our present work, we start by writing the simplest solutions of the system of transport equations obtained by assuming that all energy losses are absent, i.e.,bj=0. Then, by denoting Np,0, Nπ,0, and Nμ,0 the corresponding energy distributions for protons, pions, and muons, respectively, the system of transport equations takes the trivial form
(18)tp−1ENp,0E,ztπ−1ENπ,0E,ztμ−1ENμ,0E,z=QpE,zQπE,zQμE,z.The calculated distributionsNj,0, with j=p, π±, and μ± of the latter equations are discussed in the next section.


## 3.1.2. Steady-State Transport Equations with Particle Losses when Moving inside the Jet Plasma

Under the conditions of the steady-state approximation, and assuming that the various energy losses are absent, the system of transport equations (Equation (17)) takes the form
(19)∂∂EbpNpE,zbπNπE,zbμNμE,z+tesc−1NpE,ztπ−1NπE,ztμ−1NμE,z=QpE,zQπE,zQμE,z.In order to calculate the time-independent neutrino and gamma-ray emissivities, we need, first, to calculate the distributions of protons, pions, and muons,Np, Nπ, and Nμ, respectively, from Equation (19). For these computations, we used a code written in the C programming language, mainly following the assumptions of References [6, 24, 27].We note that, in the l.h.s. of Equation (19), the rates of the knock out processes inside the jet, tj−1, for the particles j=p,π,μ, are provided by Equation (12) and the values given below that equation. Moreover, the coefficients bj entering Equation (19) are obtained from expressions of the form bjE=−Etloss,j−1E by inserting the corresponding energy loss rates tloss,j−1 for each particle as
(20)bpE=−Etloss,p−1E=−Etsyn−1+tad−1+tpp−1,bπE=−Etloss,π−1E=−Etsyn−1+tad−1+tπp−1,bμE=−Etloss,μ−1E=−Etsyn−1+tad−1.It is worth noting that, in the latter equations, some additional energy loss terms may appear in cases when the rates of some other processes, ignored in our present work (like, e.g., the inverse Compton scattering and the pion-muon scattering), may be considered important for the description of neutrino andγ-ray emissions from other cosmic structures (AGNs, etc.).


## 3.2. Brief Description of the Model Assumptions

Before embarking to concrete calculations, we find it helpful for the reader to briefly summarize the simplifications made and the main approximations assumed in our method. At first, in purely hadronic models, only the protons from the primary charged particles (p,e+,e−) inside the jet are considered to contribute significantly in the emission of high-energy γ-rays and neutrinos, see References [4, 5] and discussion after Equation (2).As mentioned in Introduction, in the present calculations, we consider that the (GeV-TeV energy range of our interest) neutrinos andγ-rays are due to the p-p collision mechanism taking place inside the microquasar jets. The maximum nonthermal (primary) proton energy employed in our numerical integrations in the following (Section 3) is taken equal to Epmax≈106−107 GeV. We should stress that the p-p mechanism is assumed to dominate over others like, e.g., scattering of jet protons with the protons of the donor-star wind or the proton-photon scattering inside the jet (case of γ-ray production) [1].We, furthermore, note that in this work, we ignore neutrinos coming from the decay of the secondary kaons (K±) and γ-rays produced from η particles, neutral kaons, etc. As discussed in Section 2.2, the contributions originated from charged pion decays (π±) and the subsequent charged muon decays (μ±), for neutrino production, and those coming from neutral pion decays (π0), for γ-ray production, are assumed to be the main chain-reaction channels.


## 4. Injection Functions and Energy Distributions of Particles inside the Jet

In this section, the source functionsQjE,r,t entering the r.h.s of the system of coupled integrodifferential equations of transport type are discussed. Towards this aim, initially in the system of Equation (17), we insert phenomenological expressions obtained from Reference [45]. A mathematical semianalytic way of solving this system is proposed in Reference [52]. Also, following Reference [28], a method is derived for the numerical solution of the differential system of transport equations.


### 4.1. Proton Injection Function and Energy Distribution

#### 4.1.1. Relativistic Proton Injection Function

In the assumed p-p mechanism, the density of fast (relativistic) protons injected is more important, while the corresponding density of slow (thermal) protons is dynamically significant.A usual injection function for the relativistic protons, coming out of the acceleration mechanism, has been taken to be a power-law with exponent equal to two, i.e., a function of the energy of the form [53]
(21)QpE,z=Q0z0z31E2,where Q0 is a normalization constant obtained by specifying the power in the relativistic protons [24], see the appendix. The above form is valid for the jet’s frame of reference (the corresponding expression in the observer’s system is shown in the appendix).It should be noted that Equation (21) enters (through integration) the determination of the injection functions of the secondary particles (π±, μ±, see the following) and this justifies the term systems of coupled integrodifferential equations used in Section 3. The latter injection functions are also dependent on the rates and cross sections of the reactions that preceded in the chain processes following the inelastic p-p collision.


#### 4.1.2. Relativistic Proton Energy Distribution in the Jet

Under the assumptions discussed before, the system of transport equations (19) can be easily solved and the obtained proton distribution NpE,z is written as [2, 24]
(22)NpE,z=∫EpEpmaxbpE−1QpE′,ze−τpE,E′dE′,where
(23)τpE,E′=∫EE′tesc,p−1dE′′bpE′′.The quantityQpE,z corresponds to the relativistic injection function of protons at the observer’s frame (see the appendix). The minimum energy of protons is Epmin=Eth=1.22 GeV, while the maximum energy is assumed to be Epmax=106−107 GeV.Figure4 shows the proton energy distribution NpE at the base of the jet for M33 X-7 (solid line) and SS 433 (dashed line). As can be seen, the number of protons appears reduced (even by two orders of magnitude) for the wide half-opening angle system of M33 X-7 (ξ=7∘) as compared to the narrow of SS 433 (ξ=0.6∘). This is because the jet radius rz is bigger and hence the magnetic field B is smaller in the case of M33 X-7 than that of SS 433.Figure 4
Proton energy distribution for M33 X-7 and SS 433 microquasar systems.


### 4.2. Pion Injection Functions and Energy Distributions inside the Jet

#### 4.2.1. Pion Injection Functions

For pions produced through the inelastic p-p scattering, as injection function, we adopt the one given by [47] as
(24)QπE,z=nzc∫ε1σppinelExNpEx,zFπx,Exdxx,where ε=E/Epmax, x=E/Ep, and Fπ denotes the distribution of pions produced per p−p collision (see the appendix).


#### 4.2.2. Pion Energy Distributions

The steady-state energy distributionNπE,z for pions (created through the scattering of hot protons off the cold ones) results through a similar way to that of Equations (22) and (23). With the replacement p⟶π and tesc−1⟶tπ−1E,z, respectively, on the latter equations, the corresponding solution for pion energy distribution is, then, written as
(25)NπE,z=∫EEmaxbπE−1QπE′,ze−τπdE′,where the rate tπ−1 includes contributions from the decay and escape rates (see Equation (12)].


### 4.3. Muon Injection Functions and Energy Distributions inside the Jet

#### 4.3.1. Muon Injection Functions

As mentioned before, in this work in addition to neutrinos coming from the decay of charged pionsπ±, we evaluate neutrino emissivities originating from the decay of the secondary charged muons μ±. For the latter emissivity, we follow Reference [45], in order to take into account properly the muon energy loss. This means that it is necessary to consider the production of both left-handed and right-handed muons, μL− and μR+, separately, because μL− and μR+ have different decay spectra.Thus, the injection functions of the left-handed and right-handed muons are [45]
(26)QμL−,μR+Eμ,z=∫EμEmaxtdec,π−1EπNπEπ,zΘx−rπ×rπ1−xEπx1−rπ2dEπ,QμR−,μL+Eμ,z=∫EμEmaxtdec,π−1EπNπEπ,zΘx−rπ×x−rπEπx1−rπ2dEπ,with x=Eμ/Eπ and rπ=mμ/mπ2.


#### 4.3.2. Muon Energy Distributions

The corresponding muon energy distributions, coming out of the transport equation, result from Equation (25) via the replacement tesc−1⟶tμ−1E,z, for muons, and are then written as
(27)NμiE,z=∫EEmaxbμE−1QμiE′,ze−τμdE′,where the rate tμ−1 for muons is given in Equation (12).In Figure5, the pion (a) and muon (b) energy distributions for M33 X-7 are illustrated. In both cases of particles, the distribution NjE,z=z0, with j=π or μ, which considers the most important energy losses, is compared with that which neglects energy losses, Nj,0E,z=z0. The solid lines correspond to distributions obtained by considering energy losses, and the dashed lines correspond to those where the energy losses obtained have been neglected.Figure 5
Pion (a) and muon (b) distributions for M33 X-7. In both cases of particles-j, we compare the distribution which considers the important energy losses NjE,z=z0 and with that which neglects energy losses Nj,0E,z=z0 for j=π,μ. Both types of particle distributions refer to the base of the jet, i.e., at z=z0.
(a)(b)


## 4.1. Proton Injection Function and Energy Distribution

### 4.1.1. Relativistic Proton Injection Function

In the assumed p-p mechanism, the density of fast (relativistic) protons injected is more important, while the corresponding density of slow (thermal) protons is dynamically significant.A usual injection function for the relativistic protons, coming out of the acceleration mechanism, has been taken to be a power-law with exponent equal to two, i.e., a function of the energy of the form [53]
(21)QpE,z=Q0z0z31E2,where Q0 is a normalization constant obtained by specifying the power in the relativistic protons [24], see the appendix. The above form is valid for the jet’s frame of reference (the corresponding expression in the observer’s system is shown in the appendix).It should be noted that Equation (21) enters (through integration) the determination of the injection functions of the secondary particles (π±, μ±, see the following) and this justifies the term systems of coupled integrodifferential equations used in Section 3. The latter injection functions are also dependent on the rates and cross sections of the reactions that preceded in the chain processes following the inelastic p-p collision.


### 4.1.2. Relativistic Proton Energy Distribution in the Jet

Under the assumptions discussed before, the system of transport equations (19) can be easily solved and the obtained proton distribution NpE,z is written as [2, 24]
(22)NpE,z=∫EpEpmaxbpE−1QpE′,ze−τpE,E′dE′,where
(23)τpE,E′=∫EE′tesc,p−1dE′′bpE′′.The quantityQpE,z corresponds to the relativistic injection function of protons at the observer’s frame (see the appendix). The minimum energy of protons is Epmin=Eth=1.22 GeV, while the maximum energy is assumed to be Epmax=106−107 GeV.Figure4 shows the proton energy distribution NpE at the base of the jet for M33 X-7 (solid line) and SS 433 (dashed line). As can be seen, the number of protons appears reduced (even by two orders of magnitude) for the wide half-opening angle system of M33 X-7 (ξ=7∘) as compared to the narrow of SS 433 (ξ=0.6∘). This is because the jet radius rz is bigger and hence the magnetic field B is smaller in the case of M33 X-7 than that of SS 433.Figure 4
Proton energy distribution for M33 X-7 and SS 433 microquasar systems.


## 4.1.1. Relativistic Proton Injection Function

In the assumed p-p mechanism, the density of fast (relativistic) protons injected is more important, while the corresponding density of slow (thermal) protons is dynamically significant.A usual injection function for the relativistic protons, coming out of the acceleration mechanism, has been taken to be a power-law with exponent equal to two, i.e., a function of the energy of the form [53]
(21)QpE,z=Q0z0z31E2,where Q0 is a normalization constant obtained by specifying the power in the relativistic protons [24], see the appendix. The above form is valid for the jet’s frame of reference (the corresponding expression in the observer’s system is shown in the appendix).It should be noted that Equation (21) enters (through integration) the determination of the injection functions of the secondary particles (π±, μ±, see the following) and this justifies the term systems of coupled integrodifferential equations used in Section 3. The latter injection functions are also dependent on the rates and cross sections of the reactions that preceded in the chain processes following the inelastic p-p collision.


## 4.1.2. Relativistic Proton Energy Distribution in the Jet

Under the assumptions discussed before, the system of transport equations (19) can be easily solved and the obtained proton distribution NpE,z is written as [2, 24]
(22)NpE,z=∫EpEpmaxbpE−1QpE′,ze−τpE,E′dE′,where
(23)τpE,E′=∫EE′tesc,p−1dE′′bpE′′.The quantityQpE,z corresponds to the relativistic injection function of protons at the observer’s frame (see the appendix). The minimum energy of protons is Epmin=Eth=1.22 GeV, while the maximum energy is assumed to be Epmax=106−107 GeV.Figure4 shows the proton energy distribution NpE at the base of the jet for M33 X-7 (solid line) and SS 433 (dashed line). As can be seen, the number of protons appears reduced (even by two orders of magnitude) for the wide half-opening angle system of M33 X-7 (ξ=7∘) as compared to the narrow of SS 433 (ξ=0.6∘). This is because the jet radius rz is bigger and hence the magnetic field B is smaller in the case of M33 X-7 than that of SS 433.Figure 4
Proton energy distribution for M33 X-7 and SS 433 microquasar systems.


## 4.2. Pion Injection Functions and Energy Distributions inside the Jet

### 4.2.1. Pion Injection Functions

For pions produced through the inelastic p-p scattering, as injection function, we adopt the one given by [47] as
(24)QπE,z=nzc∫ε1σppinelExNpEx,zFπx,Exdxx,where ε=E/Epmax, x=E/Ep, and Fπ denotes the distribution of pions produced per p−p collision (see the appendix).


### 4.2.2. Pion Energy Distributions

The steady-state energy distributionNπE,z for pions (created through the scattering of hot protons off the cold ones) results through a similar way to that of Equations (22) and (23). With the replacement p⟶π and tesc−1⟶tπ−1E,z, respectively, on the latter equations, the corresponding solution for pion energy distribution is, then, written as
(25)NπE,z=∫EEmaxbπE−1QπE′,ze−τπdE′,where the rate tπ−1 includes contributions from the decay and escape rates (see Equation (12)].


## 4.2.1. Pion Injection Functions

For pions produced through the inelastic p-p scattering, as injection function, we adopt the one given by [47] as
(24)QπE,z=nzc∫ε1σppinelExNpEx,zFπx,Exdxx,where ε=E/Epmax, x=E/Ep, and Fπ denotes the distribution of pions produced per p−p collision (see the appendix).


## 4.2.2. Pion Energy Distributions

The steady-state energy distributionNπE,z for pions (created through the scattering of hot protons off the cold ones) results through a similar way to that of Equations (22) and (23). With the replacement p⟶π and tesc−1⟶tπ−1E,z, respectively, on the latter equations, the corresponding solution for pion energy distribution is, then, written as
(25)NπE,z=∫EEmaxbπE−1QπE′,ze−τπdE′,where the rate tπ−1 includes contributions from the decay and escape rates (see Equation (12)].


## 4.3. Muon Injection Functions and Energy Distributions inside the Jet

### 4.3.1. Muon Injection Functions

As mentioned before, in this work in addition to neutrinos coming from the decay of charged pionsπ±, we evaluate neutrino emissivities originating from the decay of the secondary charged muons μ±. For the latter emissivity, we follow Reference [45], in order to take into account properly the muon energy loss. This means that it is necessary to consider the production of both left-handed and right-handed muons, μL− and μR+, separately, because μL− and μR+ have different decay spectra.Thus, the injection functions of the left-handed and right-handed muons are [45]
(26)QμL−,μR+Eμ,z=∫EμEmaxtdec,π−1EπNπEπ,zΘx−rπ×rπ1−xEπx1−rπ2dEπ,QμR−,μL+Eμ,z=∫EμEmaxtdec,π−1EπNπEπ,zΘx−rπ×x−rπEπx1−rπ2dEπ,with x=Eμ/Eπ and rπ=mμ/mπ2.


### 4.3.2. Muon Energy Distributions

The corresponding muon energy distributions, coming out of the transport equation, result from Equation (25) via the replacement tesc−1⟶tμ−1E,z, for muons, and are then written as
(27)NμiE,z=∫EEmaxbμE−1QμiE′,ze−τμdE′,where the rate tμ−1 for muons is given in Equation (12).In Figure5, the pion (a) and muon (b) energy distributions for M33 X-7 are illustrated. In both cases of particles, the distribution NjE,z=z0, with j=π or μ, which considers the most important energy losses, is compared with that which neglects energy losses, Nj,0E,z=z0. The solid lines correspond to distributions obtained by considering energy losses, and the dashed lines correspond to those where the energy losses obtained have been neglected.Figure 5
Pion (a) and muon (b) distributions for M33 X-7. In both cases of particles-j, we compare the distribution which considers the important energy losses NjE,z=z0 and with that which neglects energy losses Nj,0E,z=z0 for j=π,μ. Both types of particle distributions refer to the base of the jet, i.e., at z=z0.
(a)(b)


## 4.3.1. Muon Injection Functions

As mentioned before, in this work in addition to neutrinos coming from the decay of charged pionsπ±, we evaluate neutrino emissivities originating from the decay of the secondary charged muons μ±. For the latter emissivity, we follow Reference [45], in order to take into account properly the muon energy loss. This means that it is necessary to consider the production of both left-handed and right-handed muons, μL− and μR+, separately, because μL− and μR+ have different decay spectra.Thus, the injection functions of the left-handed and right-handed muons are [45]
(26)QμL−,μR+Eμ,z=∫EμEmaxtdec,π−1EπNπEπ,zΘx−rπ×rπ1−xEπx1−rπ2dEπ,QμR−,μL+Eμ,z=∫EμEmaxtdec,π−1EπNπEπ,zΘx−rπ×x−rπEπx1−rπ2dEπ,with x=Eμ/Eπ and rπ=mμ/mπ2.


## 4.3.2. Muon Energy Distributions

The corresponding muon energy distributions, coming out of the transport equation, result from Equation (25) via the replacement tesc−1⟶tμ−1E,z, for muons, and are then written as
(27)NμiE,z=∫EEmaxbμE−1QμiE′,ze−τμdE′,where the rate tμ−1 for muons is given in Equation (12).In Figure5, the pion (a) and muon (b) energy distributions for M33 X-7 are illustrated. In both cases of particles, the distribution NjE,z=z0, with j=π or μ, which considers the most important energy losses, is compared with that which neglects energy losses, Nj,0E,z=z0. The solid lines correspond to distributions obtained by considering energy losses, and the dashed lines correspond to those where the energy losses obtained have been neglected.Figure 5
Pion (a) and muon (b) distributions for M33 X-7. In both cases of particles-j, we compare the distribution which considers the important energy losses NjE,z=z0 and with that which neglects energy losses Nj,0E,z=z0 for j=π,μ. Both types of particle distributions refer to the base of the jet, i.e., at z=z0.
(a)(b)


## 5. Results for Neutrino andγ-Ray Emissivities

In this section, we present and discuss simulated emissivities of high-energy neutrinos andγ-rays produced in the extragalactic microquasar M33 X-7 system which is located in the Messier 33 galaxy, at a distance ~840−960 kpc from the Earth [33]. The derived algorithms for this purpose have been tested on the well-studied galactic microquasars SS 433 and Cyg X-1 system as well.Because in our previous calculations on neutrino production from MQs [14–16] we neglected (due to the complexity and long time consuming) emissivity (intensity) of neutrinos originated from the secondary muons (μ±) produced from the charged pion (π±) decays, in this section, we present contributions also originating from the μ± channel (see Equation (5)).By using the concentrations for protonsNp, pions Nπ, and muons Nμ, the neutrino and γ-ray intensities are subsequently calculated as described in the following.


### 5.1. Neutrino Emission from the p-p Reaction Family Tree

After the above discussion, the total emissivityQνE,z produced from a MQ is the sum of contributions from the two sources: (i) the first comes from the direct π± decay (prompt neutrino production), and (ii) the second comes from the μ± decay (delayed neutrino production). Thus,
(28)QνE,z=Qπ⟶νE,z+Qμ⟶νE,z.For pion decays, the injection function is given by(29)Qπ⟶νE,z=∫EEmaxtdec,π−1EπNπEπ,zΘ1−rπ−xEπ1−rπdEπ,with x=E/Eπ, while for the four types of muon decays, the injection function reads
(30)Qμ⟶νE,z=∑i=14∫EEmaxtdec,μ−1EμNμiEμ,z×53−3x2+43x3+3x2−13−83x3hidEμEμ,with x=E/Eμ.In the last expression, the symbolsμi,i=1,2,3,4 correspond to μ1,2=μL−,+, μ3,4=μR−,+ [45], while h1,2=−h3,4=−1. Then, the calculation of each of the latter two integrals of Equations (29) and (30) provides separately the partial emissivity of neutrinos for the prompted and the delayed neutrino source, respectively. Needless to note that Earth and space telescope are not able to discriminate the two sources as it happens with laboratory neutrino sources.Subsequently, we easily obtain the neutrino intensity (in unitsGeV−1s−1) by the spatial integration
(31)IνE=∫VQνE,zd3r=πtanξ2∫z0zmaxQνE,zz2dz.Figure6(a) shows the neutrino intensity produced at the base of the jet from direct decays of secondary pions and muons coming from p-p collisions in the jets of M33 X-7 and SS 433, respectively. As can be seen, the number of produced neutrinos reduces significantly for energies E>103 GeV, following the behavior of pion and muon distribution. We can also see that the neutrinos produced in M33 X-7 are less than those produced in SS 433, due to the wider-half opening angle ξ as we explained before.Figure 6
Neutrino (a) and gamma-ray (b) intensities for M33 X-7 and SS 433.
(a)(b)It should be noted that, to perform the above calculations, we have updated and improved our codes to reduce the time consuming to a reasonable level, so as the integral involved in Equation (30) to be easily obtained, since in going the step from Equation (29) to Equation (30), the time consuming increases rapidly.


### 5.2. Gamma-Ray Emission from the p-p Reaction Chain

As we have discussed in Section2, the p-p collisions inside the jets, produce secondary neutral particles (pions, eta particle, etc.) that decay to give γ-rays. The dominant of these channels goes through the reaction
(32)π0⟶γ+γ.ForEγ≥100GeV, we consider the γ-ray emissivity at a distance z along the jets (in units GeV−1s−1) to be given by [24] as
(33)Qγ=c∫Eγ/Epmax1σppinelEγxNpEγx,zFγx,Eγxdxx,where Fγ is the spectrum of the produced γ-rays [47, 48] with energy x=Eγ/Ep for a primary proton energy Ep. The function Fγ is given in the appendix.Subsequently, the corresponding spectral intensity ofγ-rays can be obtained from the following spatial integration over the jet’s volume V as
(34)IγEγ=∫VQγEγ,zd3r=πtanξ2∫z0zmaxQγEγ,zz2dz.Figure6(b) shows the γ-ray intensity for energies Eγ>100 GeV. As can be seen, the produced γ-ray intensity reduces rather steadily following the behavior of the proton distribution in both cases of MQ systems. Again for the M33 X-7, it is lower mainly due to the wider half-opening angle (i.e., weaker magnetic field) but also due to other parameters.


### 5.3. Gamma-Ray and Neutrino Fluxes Arriving at the Earth

In the final stage of this work, we calculate neutrino andγ-ray fluxes for the studied systems and compare them with existing observational data or sensitivities of the KM3NeT, IceCude, of the ANTARES detectors (for the neutrino fluxes), and of the CTA, MAGIC, and HESS telescopes (for the γ-ray fluxes). Towards this aim, we define the differential flux of neutrinos and γ-rays arriving at the Earth as [27]
(35)dΦdE=14πd2IE,where d denotes the distance of the studied microquasar from Earth (see Table 1). Then, the differential fluxes weighted by the energy squared, E2dΦ/dE, can readily be obtained.Differential neutrino fluxes (left) and gamma-ray fluxes (right) weighted over the energy squaredE2dΦ/dE for the M33 X-7 (solid lines) and SS 433 (dashed lines). They are compared with the corresponding sensitivities of the neutrino detectors (left panel) ANTARES [54], IceCube [55], and KM3NeT [56] and the respective gamma-ray sensitivities (right) from MAGIC and HESS telescopes as well as the foreseen sensitivity of the next-generation Cherenkov Telescope Array (CTA) [57].In Figure7, we illustrate the predictions of our model for the quantity E2dΦ/dE produced at the base of the jets in the binary system M33 X-7 (solid lines), and the well-studied system SS 433 (dashed lines). In Figure 7(a), the comparison of these predictions for the neutrino fluxes with the corresponding sensitivities of the detectors ANTARES [54], IceCube [55], and KM3NeT [56] is shown. Also, in Figure 7(b), the comparison of these predictions for the γ-ray fluxes with the corresponding sensitivities of the MAGIC and HESS telescopes, as well as the foreseen sensitivity of the next-generation Cherenkov Telescope Array (CTA) [57], is presented. In this subfigure, the observations are compared with our calculations obtained as described above.Figure 7
Differential neutrino fluxes (a) and gamma-ray fluxes (b) weighted over the energy squaredE2dΦ/dE for the M33 X-7 (solid lines) and SS 433 (dashed lines). They are compared with the corresponding sensitivities of the neutrino detectors (a) ANTARES [54], IceCube [55], and KM3NeT [56] and the respective gamma-ray sensitivities (b) from MAGIC and HESS telescopes as well as the foreseen sensitivity of the next-generation Cherenkov Telescope Array (CTA) [57].
(a)(b)As can be seen, our predictions are in good agreement with the available observations for SS 433 in theγ-ray energy range 102<E<104−105 GeV. However, due to the much larger distance from Earth of the M33 X-1 system, relevant observations from present-day telescopes are not feasible. The corresponding γ-ray flux arriving at the Earth according to our calculations is about two orders of magnitude weaker than the sensitivity of the operating telescopes.


## 5.1. Neutrino Emission from the p-p Reaction Family Tree

After the above discussion, the total emissivityQνE,z produced from a MQ is the sum of contributions from the two sources: (i) the first comes from the direct π± decay (prompt neutrino production), and (ii) the second comes from the μ± decay (delayed neutrino production). Thus,
(28)QνE,z=Qπ⟶νE,z+Qμ⟶νE,z.For pion decays, the injection function is given by(29)Qπ⟶νE,z=∫EEmaxtdec,π−1EπNπEπ,zΘ1−rπ−xEπ1−rπdEπ,with x=E/Eπ, while for the four types of muon decays, the injection function reads
(30)Qμ⟶νE,z=∑i=14∫EEmaxtdec,μ−1EμNμiEμ,z×53−3x2+43x3+3x2−13−83x3hidEμEμ,with x=E/Eμ.In the last expression, the symbolsμi,i=1,2,3,4 correspond to μ1,2=μL−,+, μ3,4=μR−,+ [45], while h1,2=−h3,4=−1. Then, the calculation of each of the latter two integrals of Equations (29) and (30) provides separately the partial emissivity of neutrinos for the prompted and the delayed neutrino source, respectively. Needless to note that Earth and space telescope are not able to discriminate the two sources as it happens with laboratory neutrino sources.Subsequently, we easily obtain the neutrino intensity (in unitsGeV−1s−1) by the spatial integration
(31)IνE=∫VQνE,zd3r=πtanξ2∫z0zmaxQνE,zz2dz.Figure6(a) shows the neutrino intensity produced at the base of the jet from direct decays of secondary pions and muons coming from p-p collisions in the jets of M33 X-7 and SS 433, respectively. As can be seen, the number of produced neutrinos reduces significantly for energies E>103 GeV, following the behavior of pion and muon distribution. We can also see that the neutrinos produced in M33 X-7 are less than those produced in SS 433, due to the wider-half opening angle ξ as we explained before.Figure 6
Neutrino (a) and gamma-ray (b) intensities for M33 X-7 and SS 433.
(a)(b)It should be noted that, to perform the above calculations, we have updated and improved our codes to reduce the time consuming to a reasonable level, so as the integral involved in Equation (30) to be easily obtained, since in going the step from Equation (29) to Equation (30), the time consuming increases rapidly.


## 5.2. Gamma-Ray Emission from the p-p Reaction Chain

As we have discussed in Section2, the p-p collisions inside the jets, produce secondary neutral particles (pions, eta particle, etc.) that decay to give γ-rays. The dominant of these channels goes through the reaction
(32)π0⟶γ+γ.ForEγ≥100GeV, we consider the γ-ray emissivity at a distance z along the jets (in units GeV−1s−1) to be given by [24] as
(33)Qγ=c∫Eγ/Epmax1σppinelEγxNpEγx,zFγx,Eγxdxx,where Fγ is the spectrum of the produced γ-rays [47, 48] with energy x=Eγ/Ep for a primary proton energy Ep. The function Fγ is given in the appendix.Subsequently, the corresponding spectral intensity ofγ-rays can be obtained from the following spatial integration over the jet’s volume V as
(34)IγEγ=∫VQγEγ,zd3r=πtanξ2∫z0zmaxQγEγ,zz2dz.Figure6(b) shows the γ-ray intensity for energies Eγ>100 GeV. As can be seen, the produced γ-ray intensity reduces rather steadily following the behavior of the proton distribution in both cases of MQ systems. Again for the M33 X-7, it is lower mainly due to the wider half-opening angle (i.e., weaker magnetic field) but also due to other parameters.


## 5.3. Gamma-Ray and Neutrino Fluxes Arriving at the Earth

In the final stage of this work, we calculate neutrino andγ-ray fluxes for the studied systems and compare them with existing observational data or sensitivities of the KM3NeT, IceCude, of the ANTARES detectors (for the neutrino fluxes), and of the CTA, MAGIC, and HESS telescopes (for the γ-ray fluxes). Towards this aim, we define the differential flux of neutrinos and γ-rays arriving at the Earth as [27]
(35)dΦdE=14πd2IE,where d denotes the distance of the studied microquasar from Earth (see Table 1). Then, the differential fluxes weighted by the energy squared, E2dΦ/dE, can readily be obtained.Differential neutrino fluxes (left) and gamma-ray fluxes (right) weighted over the energy squaredE2dΦ/dE for the M33 X-7 (solid lines) and SS 433 (dashed lines). They are compared with the corresponding sensitivities of the neutrino detectors (left panel) ANTARES [54], IceCube [55], and KM3NeT [56] and the respective gamma-ray sensitivities (right) from MAGIC and HESS telescopes as well as the foreseen sensitivity of the next-generation Cherenkov Telescope Array (CTA) [57].In Figure7, we illustrate the predictions of our model for the quantity E2dΦ/dE produced at the base of the jets in the binary system M33 X-7 (solid lines), and the well-studied system SS 433 (dashed lines). In Figure 7(a), the comparison of these predictions for the neutrino fluxes with the corresponding sensitivities of the detectors ANTARES [54], IceCube [55], and KM3NeT [56] is shown. Also, in Figure 7(b), the comparison of these predictions for the γ-ray fluxes with the corresponding sensitivities of the MAGIC and HESS telescopes, as well as the foreseen sensitivity of the next-generation Cherenkov Telescope Array (CTA) [57], is presented. In this subfigure, the observations are compared with our calculations obtained as described above.Figure 7
Differential neutrino fluxes (a) and gamma-ray fluxes (b) weighted over the energy squaredE2dΦ/dE for the M33 X-7 (solid lines) and SS 433 (dashed lines). They are compared with the corresponding sensitivities of the neutrino detectors (a) ANTARES [54], IceCube [55], and KM3NeT [56] and the respective gamma-ray sensitivities (b) from MAGIC and HESS telescopes as well as the foreseen sensitivity of the next-generation Cherenkov Telescope Array (CTA) [57].
(a)(b)As can be seen, our predictions are in good agreement with the available observations for SS 433 in theγ-ray energy range 102<E<104−105 GeV. However, due to the much larger distance from Earth of the M33 X-1 system, relevant observations from present-day telescopes are not feasible. The corresponding γ-ray flux arriving at the Earth according to our calculations is about two orders of magnitude weaker than the sensitivity of the operating telescopes.


## 6. Summary and Conclusion

In this work, we address neutrino andγ-ray emissions from microquasars and X-ray binary stars (XRBs) that consist of a stellar-mass black hole (compact object) and a main-sequence donor star. We assume that these emissions originate from decay and scattering processes of the secondary particles produced through the p-p scattering mechanism, i.e., the inelastic collision of relativistic protons of the jet with the thermal ones.Such high-energy neutrinos andγ-rays are detectable by operating terrestrial and space telescopes. Among the operating detectors are the under ice IceCube (at the South Pole), the ANTARES and KM3Net (under the Mediterranean sea), and many others.In the near future, there are more accurate measurements of gamma-ray and neutrino fluxes by the coming observatories, such as the CTA (Acharya et al. 2018), which is expected to shed more light on the nature and the emission sources. From the perspective of neutrino detection, the addition of more years of data with continuous operation of IceCube will improve the sensitivity of the search for galactic sources of cosmic neutrinos. Furthermore, the next-generation detection instrument, IceCube-Gen2, a substantial expansion of IceCube, will be 10 times larger. This next-generation neutrino observatory with five times the effective area of IceCube is expected to improve the neutrino source search sensitivity by the same order [58]. With higher neutrino statistics, identifying galactic sources will become more promising.On the other hand, theoretically, by modeling the solution of the system of coupled transport equations, we were able to perform detailed calculations for various processes taking place inside the jets of galactic and extragalactic system M33 X-7 assuming hadronic content in their jets. These cooling rates enter the proton, pion, and muon energy distributions through which one obtains neutrino andγ-ray intensities.Our predictions for SS 433 in theγ-ray energy range 102<E<104−105 GeV agree well with the observations. On the other hand, due to the much larger distance from Earth of the M33 X-1 MQ, relevant observations from present-day telescopes are not feasible since the corresponding fluxes arriving at the Earth are about two orders of magnitude weaker than the available sensitivities.



---

*Source: 8146675-2022-02-08.xml*