Modeling particle transport in astrophysical outflows and simulations of associated emissions
aa r X i v : . [ a s t r o - ph . H E ] J a n MNRAS , 1–11 (2020) Preprint 12 January 2021 Compiled using MNRAS L A TEX style file v3.0
Modeling particle transport in astrophysical outflows and simulations ofassociated emissions
D. A. Papadopoulos, , ★ O. T. Kosmas † and S. Ganatsios ‡ Department of Electrical and Computer Engineering, University of Western Macedonia, Kozani, Greece Modelling and Simulation Center, MACE, University of Manchester, Sackville Street, Manchester, UK Theoretical Physics Section, University of Ioannina, GR-45110 Ioannina, Greece
Last updated 2020 June 10; in original form 2013 September 5
ABSTRACT
In this work, after making an attempt to improve the formulation of the model on particle transport within astrophysical plasmaoutflows and constructing the appropriate algorithms, we test the reliability and effectiveness of our method through numericalsimulations on well-studied Galactic microquasars as the SS 433 and the Cyg X-1 systems. Then, we concentrate on predictionsof the associated emissions, focusing on detectable high energy neutrinos and 𝛾 -rays originated from the extra-galactic M33 X-7system, which is a recently discovered X-ray binary located in the neighboring galaxy Messier 33 and has not yet been modeled indetail. The particle and radiation energy distributions, produced from magnetized 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 arecomputed by solving the system of coupled integro-differential transport equations of multi-particle processes (reactions chain)following the inelastic proton-proton (p-p) collisions. For the detection of such high energy neutrinos as well as multi-wavelength(radio, X-ray and gamma-ray) emissions, extremely sensitive detection instruments are in operation or have been designed likethe CTA, IceCube, ANTARES, KM3NeT, IceCube-Gen-2, and other space telescopes. Key words: radiation mechanisms: general – neutrinos – binaries: general – stars: winds, outflows – ISM: jets and outflows –gamma-rays: general
During the last few decades, collimated astrophysical outflows havebeen observed to emerge from a wide variety of Galactic and ex-tragalactic compact structures (Spencer 1984; Mirabel & Rodríguez1999; Romero et al. 2003; Reynoso, Romero & Christiansen 2008;Reynoso & Romero 2009; Romero et al. 2016). Among such ob-jects, the stellar scale class of microquasar (MQ) and X-Ray binarysystems (XRBs) posses prominent positions (Romero et al. 2016;Smponias & Kosmas 2011, 2014; Romero & Vila 2008; Reid et al.2011). These two-body cosmic structures consist of a collapsed stel-lar remnant, a stellar mass black hole or a neutron star (compactobject), and a companion (donor) main sequence star in coupledorbit around their center of mass.Due to the strong gravitational field pertaining around the compactobject, mass from the companion star is accreted onto the equato-rial region of the black hole forming an accretion disc. In the caseswhen the black hole is rotating rapidly and the accretion disc is ge-ometrically rather thick and hot, ejection of two powerful oppositelydirected mostly relativistic mass outflows (jets) occur perpendicularto the accretion disk (Vieyro & Romero 2017). Such systems con-stitute excellent "laboratories" for the investigation of astrophysical ★ E-mail: [email protected] (DAP) † E-mail: [email protected] (OTK) ‡ E-mail: [email protected] (SG) outflows (jets). Currently, for the detection of multi-wavelength (ra-dio, X-ray and gamma-ray) as well as high energy neutrino emissions,appreciably sensitive detection tools are in operation or have beendesigned like the IceCube, ANTARES, KM3NeT, CTA, IceCube-Gen2, etc. (Saito et al. 2009; Actis et al. 2011; Aartsen et al. 2015,2016, 2018; Adrian-Martinez et al. 2016; Albert et al. 2020).From a theoretical and phenomenological view point, due tothe presence of rather strong magnetic fields, the MQ jets aretreated as magneto-hydrodymamical flows emanating from the vicin-ity of the compact object (usually a stellar mass black hole)(Smponias & Kosmas 2015, 2017; Kosmas & Smponias 2018), andthis is also assumed to be the case in our present work. From theobserved characteristics of MQs, researchers have concluded thatthey share a lot of similarities in their physical properties withthe class of Active Galactic Nuclei (AGN), which are mostly lo-cated at the central region of galaxies. The latter, however, are enor-mously larger in scale compared to microquasars and, in addition,the evolution of AGNs is appreciably slower than that of MQs thatmakes the observation of many phenomena significantly difficult(Remillard & McClintock 2006; Charles & Coe 2006; Orosz 2003;Fabrika 2004; Cherepashchuk et al. 2005).In this article, we focus on the magnetized astrophysical outflows(jets) that are characterized by hadronic content (p, 𝜋 ± , light nu-clei, etc.) in their jets. We consider them as fluid flows emanatingfrom a central source at the jet’s origin (Mirabel & Rodríguez 1999; © Smponias & Kosmas 2011, 2014). We, furthermore, assume that thecompact object is a spinning black hole with mass up to few tens(30-50) of the Sun’s mass (Smponias & Kosmas 2014; Romero et al.2016; Reynoso, Romero & Christiansen 2008).Initially, we make an effort to improve the model em-ployed to perform numerical simulations for particle and ra-diation emissions from hadronic MQ jets (Papadopoulos 2020;Papadopoulos, Papavasileiou & Kosmas 2020). We adopt that themain mechanism producing high energy neutrinos and 𝛾 -rays, is theproton-proton (p-p) collision taking place within the hadronic jets,i.e. the inelastic scattering of hot (non-thermal) protons on thermal(cold) ones (Reynoso, Romero & Christiansen 2008).The scattering, diffusion, decay, etc., of the secondary particles( 𝜋 ± , 𝐾 ± , 𝜇 ± , 𝑒 ± , etc.) produced afterwards, from a mathematicalmodeling point of view, are governed by a system of coupled integro-differential transport equations. One of the purposes of this workis to attempt, for a first time, to formulate in a compact way thisdifferential system of transport equations (satisfied by the primaryprotons and the secondary multi-particles, multi-species) and deriveadvantageous algorithms to perform the required simulations for theaforementioned emissivities (Tsoulos, Kosmas & Stavrou 2017).Moreover, we are intended to extend the calculations of Ref.(Smponias & Kosmas 2014, 2015, 2017; Kosmas & Smponias 2018)so as to include contributions to neutrino emissivity (intensity) orig-inating from the secondary muons ( 𝜇 ± ) produced from the decay ofcharged pions ( 𝜋 ± ) which in our previous works had been ignored dueto the long time consuming required for such simulations. Obviously,the mathematical problem of such a study becomes also more compli-cated (e.g. the number of coupled equations in the above mentioneddifferential system increases rapidly in the case of considering neu-trinos coming out of kaon ( 𝐾 ± ) decays (Lipari, Lusignoli & Meloni2009).After fixing the model parameters and testing the derived algo-rithms on the reproducibility of some known properties of the wellstudied SS 433 Galactic micro-quasar (Smponias & Kosmas 2015,2017; Kosmas & Smponias 2018), we perform detailed simulationsfor the Galactic Cyg X-1 system as well as for the extragalactic M33X-7 MQ (Pietsch et al. 2006). The latter system has not been studiedup to now from a jet emission view point because this is a ratherrecently discovered X-ray binary system located in the neighboringgalaxy of our Milky Way Galaxy, known as Messier 33 (Pietsch et al.2006).The rest of the paper is organized as follows. In Sect. 2, we de-scribe briefly the main characteristics of hadronic jets in Galacticmicro-quasars (SS433 and Cyg X-1) as well as the extragalacticM33 X-7 system. Then (Sect. 3), we present the formulation of ourimproved method related to the differential system of transport equa-tions (integro-differential system of coupled equations) and derivethe appropriate algorithms. In Sect. 4 we present and discuss ourresults referred to high energy neutrino and 𝛾 -ray energy-spectraemitted from the extragalactic M33 X-7 MQ. Finally (Sect. 5), wesummarize the main conclusions extracted from the present investi-gation. In this section we summarize briefly some basic properties of thehadronic MQ systems (like those mentioned before) and the char-acteristics of hadronic models that describe reliably microquasarjet emissions (Romero et al. 2016; Smponias & Kosmas 2011, 2014; Vieyro & Romero 2017). The last few decades, from the investigationof a great number of MQs and X-ray emitting binary systems, stellar-mass black holes have been observed, with masses in the regionof our interest determined mainly from the dynamics and structureproperties of their companion stars (Remillard & McClintock 2006;Charles & Coe 2006; Orosz 2003).The Galactic binary system SS433, located 5.5 kpc from the Earthin Aquila constellation, displays two mildly relativistic jets (withbulk velocity 𝜐 𝑏 ≈ . 𝑃 ∼ . 𝑀 𝐵𝐻 = 𝑀 ⊙ and 𝑀 𝑑𝑜𝑛 = 𝑀 ⊙ for the black hole and the com-panion star, respectively.The second Galactic system addressed in this work, Cygnus X-1, isan X-ray source in the constellation Cygnus which is located 1.86 kpcfrom the Earth (Reid et al. 2011). The 19 . 𝑀 ⊙ O-type companionstar (Orosz et al. 2011) is in coupled orbit with the 14 . 𝑀 ⊙ com-pact object (stellar mass black hole) (Orosz et al. 2011), with orbitalperiod 𝑃 ∼ . 𝛾 -ray and neutrino emissions from the recentlydiscovered binary system M33 X-7 (Long et al. 1981) which is lo-cated in the nearby galaxy Messier 33, the only known up to nowblack hole that is in an eclipsing binary system, with orbital period 𝑃 ∼ .
45 days (Pietsch et al. 2006). The distance from Earth of thissystem is between ∼
840 kpc (Orosz et al. 2008) and ∼
960 kpc(Bonanos, Stanek et al. 2006). In addition, the black hole mass is 𝑀 𝐵𝐻 = ( . ± . ) 𝑀 ⊙ which orbits its 𝑀 𝑑𝑜𝑛 = ( . ± . ) 𝑀 ⊙ O-type companion star.In Table 1, we tabulate in addition some other important parametersemployed in this work for the systems under investigation, SS 433,Cygnus X-1 and M33 X-7.
In the model considered in this work for the description of MQjets, an accretion disk in the equatorial region of the compactobject is present, and a fraction of the accreted material is ex-pelled in two oppositely directed jets (Falcke & Biermann 1995;Smponias & Kosmas 2011, 2014). We assume an approximately con-ical mass outflow (jet) with a half-opening angle 𝜉 (for the M33 X-7MQ, we adopt the value 𝜉 = ◦ ) and a radius given by 𝑟 ( 𝑧 ) = 𝑧 tan 𝜉 (the coordinate z-axis coincides with the cone axis, jet’s direction).The injection point (plane) of the jet is at a distance 𝑧 from the com-pact object. When 𝑧 = 𝑧 the radius of the jet is given by 𝑟 = 𝑧 tan 𝜉 (see Table I).The initial jet radius is 𝑟 = 𝑟 ( 𝑧 ) ≈ 𝑅 𝑠𝑐ℎ , where 𝑅 𝑠𝑐ℎ = 𝐺 𝑀 𝐵𝐻 / 𝑐 ( 𝑅 𝑠𝑐ℎ is the known Schwarzschild radius, the criticalblack hole radius for which the escape speed of particles becomesequal to the speed of light 𝑐 ), with 𝐺 denoting the Newton’s gravi-tational constant. Then, we find, for example, that the injection pointis at 𝑧 = 𝑟 / tan 𝜉 ≃ . × cm, for the M33 X-7. For the sake ofcomparison, we mention that for Cygnus X-1, the injection point isat 𝑧 ≃ cm, while for SS 433 it is 𝑧 ≃ . × 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 𝑧 𝑚𝑎𝑥 , which is assumed to be equal to 𝑧 𝑚𝑎𝑥 = 𝑧 MNRAS , 1–11 (2020)
Table 1.
Parameters of the model for M33 X-7, Cygnus X-1 and SS433Parameter Symbol SS433 Cygnus X-1 M33 X-7Black Hole mass 𝑀 𝐵𝐻 . 𝑀 ⊙ . 𝑀 ⊙ . 𝑀 ⊙ Distance from Earth 𝑑 𝑀 𝑑𝑜𝑛 𝑀 ⊙ . 𝑀 ⊙ 𝑀 ⊙ Donor star type - A-type O-type O-typeOrbital Period 𝑃 𝐿 𝑘 𝑒𝑟 𝑔 𝑠 − 𝑒𝑟 𝑔 𝑠 − 𝑒𝑟 𝑔 𝑠 − Jet’s launching point 𝑧 . × cm 10 cm 1 . × cmBulk velocity of jet particles 𝜐 𝑏 Γ 𝑏 𝜉 . ◦ . ◦ ◦ Jet’s viewing angle 𝜃 . ◦ . ◦ . ◦ (this point may be assumed to be at the boundaries of the jet with theambient region).In discussing the energetic evolution and dynamics of MQ jets,the kinetic energy density of the jet, 𝜌 𝑘 ( 𝑧 ) , is related to its kineticluminosity, 𝐿 𝑘 , through the expression (Romero et al. 2016) 𝜌 𝑘 ( 𝑧 ) = 𝐿 𝑘 𝜋𝜐 𝑏 [ 𝑟 ( 𝑧 )] , (1)where 𝜐 𝑏 is the bulk velocity of the jet particles (mostly pro-tons). Within the context of the jet-accretion coupling hypoth-esis, only around 10% of the Eddington luminosity goes intothe jet (Körding, Fender & Migliari 2006). Here we adopt 𝐿 𝑘 = 𝑒𝑟𝑔 𝑠 − , for a black hole of mass ∼ . 𝑀 ⊙ , i.e. for the M33X-7 system.Furthermore, assuming equipartition between the magnetic en-ergy and the kinetic energy in the jet, it implies that 𝜌 𝑚𝑎𝑔 = 𝜌 𝑘 (Bosch-Ramon, Romero & Paredes 2006), and hence the mag-netic field that colimates the astrophysical jet-plasma is given by(Reynoso & Romero 2009) 𝐵 ( 𝑧 ) ≡ q 𝜋𝜌 𝑚𝑎𝑔 ( 𝑧 ) = p 𝜋𝜌 𝑘 ( 𝑧 ) . (2)In general, for the kinetic power in the jet, many authors con-sider that a fraction is carried by primary protons and the rest byelectrons, i.e., 𝐿 = 𝐿 𝑝 + 𝐿 𝑒 . Then, the relation between the pro-ton and electron power is determined through a parameter 𝛼 insuch a way that 𝐿 𝑝 = 𝛼𝐿 𝑒 (Reynoso, Romero & Christiansen 2008;Reynoso & Romero 2009; Vieyro & Romero 2017).In this article, however, we adopt the case of 𝛼 =
100 which meansthat we consider proton-dominated jet. Other quantities needed forthe purposes of this paper are listed in Table I.
The collision of relativistic protons with the cold ones inside the jet(p-p collision mechanism), produces high energy charged particles(pions 𝜋 ± , kaons 𝐾 ± , muons 𝜇 ± , etc.) and neutral particles ( 𝜋 , 𝐾 ,˜ 𝐾 , 𝜂 particles, etc.). Important primary reactions of this type are 𝑝 + 𝑝 → 𝑝 + 𝑝 + 𝑎𝜋 + 𝑏 ( 𝜋 + + 𝜋 − ) ,𝑝 + 𝑝 → 𝑝 + 𝑛 + 𝜋 + + 𝑎𝜋 + 𝑏 ( 𝜋 + + 𝜋 − ) ,𝑝 + 𝑝 → 𝑛 + 𝑛 + 𝜋 + + 𝑎𝜋 + 𝑏 ( 𝜋 + + 𝜋 − ) , where 𝑎 and 𝑏 denote the pion multiplicities (Romero et al. 2016),and similarly for kaon production (Lipari, Lusignoli & Meloni2009). Charged pions 𝜋 ± (and kaons 𝐾 ± ), afterwards, decay to chargedleptons (muons 𝜇 ± , and electrons or positrons, 𝑒 ± ) as well as neutri-nos and anti-neutrinos as 𝜋 + → 𝜇 + + 𝜈 𝜇 , 𝜋 + → 𝑒 + + 𝜈 𝑒 𝜋 − → 𝜇 − + ¯ 𝜈 𝜇 , 𝜋 − → 𝑒 − + ¯ 𝜈 𝑒 . (3)Furthermore, muons also decay giving neutrinos and electrons (orpositrons) as 𝜇 + → 𝑒 + + 𝜈 𝑒 + ¯ 𝜈 𝜇 , 𝜇 − → 𝑒 − + ¯ 𝜈 𝑒 + 𝜈 𝜇 . (4)We mention that, neutral pions ( 𝜋 ), 𝜂 -particles, etc. decay pro-ducing 𝛾 -rays according to the reactions 𝜋 → 𝛾 + 𝛾 , 𝜋 → 𝛾 + 𝑒 − + 𝑒 + (5) 𝜂 → 𝛾 + 𝛾 , 𝜂 → 𝜋 + 𝑒 − + 𝑒 + . (6)In our present study, we consider neutrinos generated from both typesof charged-pion decays (for neutrinos coming from kaon decays thereader is referred, e.g. to Ref. (Lipari, Lusignoli & Meloni 2009)and also from both types of charged muon decays (Papadopoulos2020; Papadopoulos, Papavasileiou & Kosmas 2020). Moreover, weconsider 𝛾 -rays produced from the 𝜋 -decays of Eq. (5). Thus, high-energy neutrinos are inevitably accompanied by pionic gamma-rays,a phenomenon well known as multi-messenger emission from mi-croquasar jets. In general, the primary charged particles (p and 𝑒 − ) gain energyduring moving within the magnetic field 𝐵 (Fermi acceleration). Theacceleration rate of the (initially cold) protons to an energy 𝐸 = 𝐸 𝑝 ,known as shock acceleration, is defined as 𝑡 − 𝑎𝑐𝑐 = 𝐸 − 𝑑𝐸 / 𝑑𝑡 and isgiven by the relation 𝑡 − 𝑎𝑐𝑐 ≡ 𝐸 𝑑𝐸𝑑𝑡 ≈ 𝜂 𝑐𝑒𝐵𝐸 𝑝 (7)( 𝜂 = . 𝜎 𝑖𝑛𝑒𝑙𝑝𝑝 and 𝑡 − 𝑝𝑝 respectively, play crucial role (the in-elasticity coefficient is 𝐾 𝑝𝑝 ≈ / 𝜎 𝑖𝑛𝑒𝑙𝑝𝑝 , as a function MNRAS000
The collision of relativistic protons with the cold ones inside the jet(p-p collision mechanism), produces high energy charged particles(pions 𝜋 ± , kaons 𝐾 ± , muons 𝜇 ± , etc.) and neutral particles ( 𝜋 , 𝐾 ,˜ 𝐾 , 𝜂 particles, etc.). Important primary reactions of this type are 𝑝 + 𝑝 → 𝑝 + 𝑝 + 𝑎𝜋 + 𝑏 ( 𝜋 + + 𝜋 − ) ,𝑝 + 𝑝 → 𝑝 + 𝑛 + 𝜋 + + 𝑎𝜋 + 𝑏 ( 𝜋 + + 𝜋 − ) ,𝑝 + 𝑝 → 𝑛 + 𝑛 + 𝜋 + + 𝑎𝜋 + 𝑏 ( 𝜋 + + 𝜋 − ) , where 𝑎 and 𝑏 denote the pion multiplicities (Romero et al. 2016),and similarly for kaon production (Lipari, Lusignoli & Meloni2009). Charged pions 𝜋 ± (and kaons 𝐾 ± ), afterwards, decay to chargedleptons (muons 𝜇 ± , and electrons or positrons, 𝑒 ± ) as well as neutri-nos and anti-neutrinos as 𝜋 + → 𝜇 + + 𝜈 𝜇 , 𝜋 + → 𝑒 + + 𝜈 𝑒 𝜋 − → 𝜇 − + ¯ 𝜈 𝜇 , 𝜋 − → 𝑒 − + ¯ 𝜈 𝑒 . (3)Furthermore, muons also decay giving neutrinos and electrons (orpositrons) as 𝜇 + → 𝑒 + + 𝜈 𝑒 + ¯ 𝜈 𝜇 , 𝜇 − → 𝑒 − + ¯ 𝜈 𝑒 + 𝜈 𝜇 . (4)We mention that, neutral pions ( 𝜋 ), 𝜂 -particles, etc. decay pro-ducing 𝛾 -rays according to the reactions 𝜋 → 𝛾 + 𝛾 , 𝜋 → 𝛾 + 𝑒 − + 𝑒 + (5) 𝜂 → 𝛾 + 𝛾 , 𝜂 → 𝜋 + 𝑒 − + 𝑒 + . (6)In our present study, we consider neutrinos generated from both typesof charged-pion decays (for neutrinos coming from kaon decays thereader is referred, e.g. to Ref. (Lipari, Lusignoli & Meloni 2009)and also from both types of charged muon decays (Papadopoulos2020; Papadopoulos, Papavasileiou & Kosmas 2020). Moreover, weconsider 𝛾 -rays produced from the 𝜋 -decays of Eq. (5). Thus, high-energy neutrinos are inevitably accompanied by pionic gamma-rays,a phenomenon well known as multi-messenger emission from mi-croquasar jets. In general, the primary charged particles (p and 𝑒 − ) gain energyduring moving within the magnetic field 𝐵 (Fermi acceleration). Theacceleration rate of the (initially cold) protons to an energy 𝐸 = 𝐸 𝑝 ,known as shock acceleration, is defined as 𝑡 − 𝑎𝑐𝑐 = 𝐸 − 𝑑𝐸 / 𝑑𝑡 and isgiven by the relation 𝑡 − 𝑎𝑐𝑐 ≡ 𝐸 𝑑𝐸𝑑𝑡 ≈ 𝜂 𝑐𝑒𝐵𝐸 𝑝 (7)( 𝜂 = . 𝜎 𝑖𝑛𝑒𝑙𝑝𝑝 and 𝑡 − 𝑝𝑝 respectively, play crucial role (the in-elasticity coefficient is 𝐾 𝑝𝑝 ≈ / 𝜎 𝑖𝑛𝑒𝑙𝑝𝑝 , as a function MNRAS000 , 1–11 (2020) of the fast proton energy 𝐸 𝑝 , reads (Kelner, Aharonian & Bugayov2006, 2009) 𝜎 𝑖𝑛𝑒𝑙𝑝𝑝 ( 𝐸 𝑝 ) = ( . 𝐿 + . 𝐿 + . ) (cid:20) − (cid:16) 𝐸 𝑡ℎ 𝐸 𝑝 (cid:17) (cid:21) × − 𝑐𝑚 (8)where 𝐿 = 𝑙𝑛 ( 𝐸 𝑝 𝐺𝑒𝑉 ) , and 𝐸 𝑡ℎ denotes the minimum energyof (thermal) protons, 𝐸 𝑚𝑖𝑛𝑝 , which is equal to 𝐸 𝑚𝑖𝑛𝑝 ≡ 𝐸 𝑡ℎ = . 𝑡 − 𝑝𝑝 , of the inelastic p-p scattering is givenin Table 2 in terms of the cross section 𝜎 𝑖𝑛𝑒𝑙𝑝𝑝 and the number densityof cold particles (jet protons), 𝑛 𝑝 ( 𝑧 ) , at a distance 𝑧 from the blackhole which may be written as () 𝑛 𝑝 ( 𝑧 ) = ( − 𝑞 𝑟𝑒𝑙 ) Γ 𝑏 𝑚 𝑝 𝑐 𝜌 𝑘 ( 𝑧 ) , (9)( 𝑞 𝑟𝑒𝑙 =0.1 denotes the portion of the relativistic protons).In Table 2, in addition to the p-p collision rate and acceleratingrates of protons, we tabulate the most significant cooling rates (syn-chrotron,adiabatic, etc.) for protons, pions and muons as well.In more detail, the well known synchrotron radiation is emittedfrom charged particles of energy 𝐸 = 𝛾𝑚𝑐 , with 𝛾 being the par-ticle’s Lorentz factor and 𝑚 its mass, moving inside the magnetizedplasma. Thus, the rate of the synchrotron radiation, 𝑡 − 𝑠𝑦𝑛 , is veryimportant in MQ systems with strong magnetic fields (see Table 2).Furthermore, due to the adiabatic expansion of the jet plasma,all co-moving particles inside it loose energy with a rate 𝑡 − 𝑎𝑑 (Bosch-Ramon, Romero & Paredes 2006). This rate, known as adia-batic cooling rate, depends on the particle velocities 𝜐 𝑏 tan 𝜉 .In addition, 𝑡 − 𝜋 𝑝 gives the rate of pion-proton ( 𝜋 − 𝑝 ) inelasticscattering inside the jet. The corresponding 𝜎 𝑖𝑛𝑒𝑙𝜋 𝑝 cross section is re-lated to the inelastic p-p scattering cross section with the relationship 𝜎 𝜋 𝑝 ( 𝐸 ) = ( / ) 𝜎 𝑖𝑛𝑒𝑙𝑝𝑝 ( 𝐸 ) , where the factor 2 / 𝑡 − 𝑒𝑠𝑐 , given by 𝑡 − 𝑒𝑠𝑐 ≈ 𝑐𝑧 𝑚𝑎𝑥 − 𝑧 , (10)( 𝑧 𝑚𝑎𝑥 − 𝑧 represents the length of the acceleration zone), and (ii)the decay rate of the particle in question, 𝑡 − 𝑑𝑒𝑐 , which is known fromdirect measurements of particle’s life time.The rates of the knock out (or catastrophic) processes inside the jetfor pions, 𝑡 − 𝜋 , and muons, 𝑡 − 𝜇 , include contributions from the decayprocess and escaping process as 𝑡 − 𝜋 ( 𝐸, 𝑧 ) = 𝑡 − 𝑒𝑠𝑐, 𝜋 ( 𝑧 )+ 𝑡 − 𝑑𝑒𝑐, 𝜋 ( 𝐸 ) , 𝑡 − 𝜇 ( 𝐸, 𝑧 ) = 𝑡 − 𝑒𝑠𝑐,𝜇 ( 𝑧 )+ 𝑡 − 𝑑𝑒𝑐,𝜇 ( 𝐸 ) , (11)where 𝑡 − 𝑑𝑒𝑐, 𝜋 = [( . × − ) 𝛾 𝜋 ] − 𝑠 − , for pions, and 𝑡 − 𝑑𝑒𝑐,𝜇 = [( . × − ) 𝛾 𝜇 ] − 𝑠 − , for muons. We note that Eq. (11) holdsalso for protons due to the fact that 𝑡 − 𝑑𝑒𝑐,𝑝 = -6 -4 -2 t - [ s - ] E [GeV]
Proton cooling rates
Synchrotron emissionp-p collisionAdiabatic cooling -4 -2 t - [ s - ] E [GeV]
Proton cooling rates
Synchrotron emissionp-p collisionAdiabatic cooling
Figure 1.
Cooling rates for protons in the jets of M33 X-7 (top) and SS433(bottom) at the base of the jets 𝑧 . The plots of cooling rates show thesynchrotron emission (solid lines), the adiabatic cooling (dotted lines), andthe inelastic p-p collision rate (dashed lines). In Figs. 1, 2 and 3, the variation of the aforementioned main rates,versus the particle energy E, for protons, pions and muons, respec-tively, and for the SS 433 and M33 X-7 MQ systems, are illustrated.In other words, in these figures the region of dominance of the accel-erating, cooling, etc. rates throughout the energy range of interest, 𝐸 𝑡ℎ ≡ 𝐸 𝑚𝑖𝑛 𝐸 𝐸 𝑚𝑎𝑥 , are demonstrated.More specifically, Fig. 1 indicate the cooling rates at the base of thejet for protons, Fig. 2 those for pions and Fig. 3 the correspondingones for muons, in the extragalactic binary system M33 X-7 (topof these figures). For comparison, the corresponding results for theGalactic binary system SS 433 (bottom of these figures) are alsoshown.The plots in Fig. 2 for pions, show respectively the variation versusthe energy 𝐸 of the rate for synchrotron emission (solid lines), forthe pion-proton inelastic collision (dashed lines), for the adiabaticcooling (dotted lines), and for the decay of pions (dot-dashed lines).They refer to the M33 X-7 (top) and the SS 433 (bottom) micro-quasar systems. Also, Fig. 3 illustrates the cooling rates of muonsfor the synchrotron emission (solid lines), the adiabatic expansion(dotted lines), and for the decay of muons (dot-dashed lines). Theyare similar to those for pions (Fig. 2) except the corresponding muon-proton and muon-pion scattering channels which are ignored in thiswork.From Figs. 2 and 3 it becomes obvious that, the particle syn-chrotron losses dominate the high energy region. In the case of pro-tons (Fig. 1), however, due to the large proton mass, synchrotronlosses are not dominant up to very high energies ( ∼ GeV) for theMQ systems studied. On the other hand, for pions and muons, thedecay losses dominate for lower energies, due to their short life time(see Table 2).The main difference between these two MQ systems (SS 433 andM33 X-7) concerns the synchrotron cooling rates. In general, accord-ing to Eqs. (1) and (2), for systems with wide half-opening angle (e.g.M33 X-7, with 𝜉 = ◦ ), the magnetic energy density is lower thanthat of systems with small 𝜉 (e.g. the SS 433 with 𝜉 = . ◦ ). Hence,the magnetic field for M33 X-7 is lower which becomes obvious bycomparing Eq. (2) for the two systems. The latter conclusion justifiesthe lower synchrotron loss rate for the M33 X-7 system. MNRAS , 1–11 (2020)
Table 2.
Accelerating, cooling and decay rates for particles moving inside hadronic MQ jetsRate Parameter Rate Symbol Basic rate definitionProton accelerating rates 𝑡 − 𝑎𝑐𝑐 𝜂 𝑐𝑒𝐵𝐸 𝑝 Proton-proton collision rate 𝑡 − 𝑝𝑝 𝑛 𝑝 ( 𝑧 ) 𝜎 𝑖𝑛𝑒𝑙𝑝𝑝 ( 𝐸 𝑝 ) 𝐾 𝑝𝑝 Synchrotron radiation rate 𝑡 − 𝑠𝑦𝑛 (cid:0) 𝑚 𝑒 𝑚 (cid:1) 𝛾𝜎 𝑇 𝐵 𝑚 𝑒 𝑐 𝜋 Adiabatic expansion rate 𝑡 − 𝑎𝑑 𝜐 𝑏 𝑧 Pion-proton collision 𝑡 − 𝜋𝑝 . 𝑐 𝑛 𝑝 ( 𝑧 ) 𝜎 𝑖𝑛𝑒𝑙𝜋𝑝 ( 𝐸 𝑝 ) Proton escape rate 𝑡 − 𝑒𝑠𝑐 𝑐 /( 𝑧 𝑚𝑎𝑥 − 𝑧 ) Pion decay rate 𝑡 − 𝜋 𝑡 − 𝑒𝑠𝑐,𝜋 + 𝑡 − 𝑑𝑒𝑐,𝜋 Muon decay rate 𝑡 − 𝜇 𝑡 − 𝑒𝑠𝑐,𝜇 + 𝑡 − 𝑑𝑒𝑐,𝜇 -6 -4 -2 t - [ s - ] E [GeV]
Pion cooling rates
Synchrotron emissionpion-proton collisionAdiabatic coolingDecay -4 -2 t - [ s - ] E [GeV]
Pion cooling rates
Synchrotron emissionpion-proton collisionAdiabatic coolingDecay
Figure 2.
Cooling rates for pions in the jets of M33 X-7 (top) and SS433(bottom) at the base of the jets -2 t - [ s - ] E [GeV]
Muon cooling rates
Synchrotron emissionAdiabatic coolingDecay -2 t - [ s - ] E [GeV]
Muon cooling rates
Synchrotron emissionAdiabatic coolingDecay
Figure 3.
Cooling rates for muons in the jets of M33 X-7 (top) and SS433(bottom) at the base of the jets
In this section, we make a first step towards improving the mathemat-ical formulation describing the creation, scattering, decay, diffusionand emission of particles (including neutrinos) and electromagneticradiation in astrophysical outflows (jets). The basic key-role tool ofthis formalism is the general transport equation which satisfies eachof the particles considered (Vieyro & Romero 2017).We should note that, such a formalism may be reliably applica-ble for multi-species (multi-particle) emissions from MQs and X-raybinaries, as well as from central regions of galaxies like the AGN sys-tems, in which the compact object is a massive or supermassive blackhole, see e.g. (Romero et al. 2016; Ginzburg & Syrovatskii 1964).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, thesources are assumed conical and the particles (neutrinos) or radiationare emitted to the forward or backward direction.From a mathematical point of view, the transport equation whichsatisfies any kind of particles moving into the (dark) jets of MQs, isan integro-differential equation, as is explained below.
The general transport equation describes the concentration (distri-bution) of particle of j-kind, 𝑁 𝑗 ( 𝐸, r , 𝑡 ) , where 𝑗 = 𝑝, 𝑒 ± , 𝜋 ± , 𝜇 ± ,etc., as a function of the time 𝑡 , the particle’s energy 𝐸 , and the posi-tion r inside the (conical) jet. In essence, this is a phenomenologicalmacroscopic equation of particle (or radiation) transport describingastrophysical outflows, and it is written as 𝜕𝑁 𝑗 𝜕𝑡 − ∇ · (cid:0) 𝐷 𝑗 ∇ 𝑁 𝑗 (cid:1) + 𝜕 (cid:0) 𝑏 𝑗 𝑁 𝑗 (cid:1) 𝜕𝐸 − 𝜕 (cid:0) 𝑑 𝑗 𝑁 𝑗 (cid:1) 𝜕𝐸 = 𝑄 𝑗 ( 𝐸, r , 𝑡 ) − 𝑝 𝑗 𝑁 𝑗 + Õ 𝑘 ∫ 𝑃 𝑘𝑗 ( 𝐸 ′ , 𝐸 ) 𝑁 𝑗 ( 𝐸, r , 𝑡 ) 𝑑𝐸 , 𝑗 = 𝑝, 𝜋 ± , 𝜇 ± ... (12)where the parameters 𝐷 𝑗 , 𝑏 𝑗 , 𝑑 𝑗 , 𝑝 𝑗 and 𝑃 𝑘𝑗 may depend upon thespace and time coordinates and also on the energy 𝐸 . The latterequation coincides with the continuity equation for particles of type- 𝑗 , with j=1,2,3 and ( , , ) ≡ (p, 𝜋 , 𝜇 )]. The term 𝑄 𝑗 ( 𝐸, r , 𝑡 ) in the r.h.s. of Eq. (12) is equal to the intensityof the source producing the particles- 𝑗 , which is also known as theinjection function of particles-j. This means that 𝑄 𝑗 ( 𝐸, r , 𝑡 ) 𝑑𝐸 𝑑 r 𝑑𝑡 represents the number of particles kind- 𝑗 provided by the sources ina volume element 𝑑 r , in the energy range between 𝐸 and 𝐸 + 𝑑𝐸 MNRAS , 1–11 (2020) during the time 𝑑𝑡 . In the case when the 𝑗 -type particles are productsof a chain reaction (as it holds in our present work assuming thep-p reaction chain described in Sect. 2.2), the function 𝑄 𝑗 ( 𝐸, r , 𝑡 ) couples the 𝑗 -reaction with its parent reaction, i.e. the equation ofparticles kind-(j-1).Further, the term proportional to 𝑝 𝑗 exists in cases when catas-trophic (or knock out) processes take place that cause catastrophicenergy losses. Then, this gives the probability per unit time for thelosses in question to occur. Thus, the coefficient 𝑝 𝑗 is written as 𝑝 𝑗 = 𝑇 𝑗 ≡ 𝑡 − 𝑗 , (13)with 𝑇 𝑗 is the mean life time of the particles of kind- 𝑗 . Also 𝑝 𝑗 𝑁 𝑗 denotes the number of particles "knocked out" per unit time.The coefficients 𝑏 𝑗 in the third term (r.h.s.) are equal to the meanenergy increment of the particle- 𝑗 per unit time i.e. 𝑏 𝑗 = 𝑑𝐸𝑑𝑡 . (14) For the special case when 𝑏 𝑗 = 𝑏 𝑗 ( 𝐸 ) , this coefficient is relatedto energy losses of various cooling processes. Also, the coefficients 𝑃 𝑘𝑗 ( 𝐸 ) , with 𝑘 > 𝑗 , in the latter summation of the transport equation,are related with the existence of fragmentation of the primary (orsecondary) particles which in this work are assumed as non-existing.The parameter 𝐷 𝑗 is known as the diffusion coefficient which, ingeneral, is a functions of the coordinates r and time 𝑡 if the concentra-tion of the jet plasma and its macroscopic motion are inhomogeneousin the volume of the jet. By assuming that 𝐷 𝑗 = B .For simplicity, in this work we make the realistic assumptions that,the above coefficients depend only on the particle’s energy 𝐸 , i.e. weassume that 𝐷 𝑗 ( 𝐸 ) , 𝑏 𝑗 ( 𝐸 ) and 𝑝 𝑗 ( 𝐸 ) . Then, the system of generaltransport equations (12) reads 𝜕𝜕𝑡 © « 𝑁 𝑝 ( 𝐸, 𝑧, 𝑡 ) 𝑁 𝜋 ( 𝐸, 𝑧, 𝑡 ) 𝑁 𝜇 ( 𝐸, 𝑧, 𝑡 ) ª®¬ −∇ © « 𝐷 𝑝 ( 𝐸 ) 𝑁 𝑝 ( 𝐸, 𝑧, 𝑡 ) 𝐷 𝜋 ( 𝐸 ) 𝑁 𝜋 ( 𝐸, 𝑧, 𝑡 ) 𝐷 𝜇 ( 𝐸 ) 𝑁 𝜇 ( 𝐸, 𝑧, 𝑡 ) ª®¬ + 𝜕𝜕𝐸 © « 𝑏 𝑝 ( 𝐸 ) 𝑁 𝑝 ( 𝐸, 𝑧, 𝑡 ) 𝑏 𝜋 ( 𝐸 ) 𝑁 𝜋 ( 𝐸, 𝑧, 𝑡 ) 𝑏 𝜇 ( 𝐸 ) 𝑁 𝜇 ( 𝐸, 𝑧, 𝑡 ) ª®¬ − 𝜕 𝜕𝐸 © « 𝑑 𝑝 ( 𝐸 ) 𝑁 𝑝 ( 𝐸, 𝑧, 𝑡 ) 𝑑 𝜋 ( 𝐸 ) 𝑁 𝜋 ( 𝐸, 𝑧, 𝑡 ) 𝑑 𝜇 ( 𝐸 ) 𝑁 𝜇 ( 𝐸, 𝑧, 𝑡 ) ª®¬ + © « 𝑡 − 𝑒𝑠𝑐 𝑁 𝑝 ( 𝐸, 𝑧, 𝑡 ) 𝑡 − 𝜋 𝑁 𝜋 ( 𝐸, 𝑧, 𝑡 ) 𝑡 − 𝜇 𝑁 𝜇 ( 𝐸, 𝑧, 𝑡 ) ª®®¬ = © « 𝑄 𝑝 ( 𝐸, 𝑧, 𝑡 ) 𝑄 𝜋 ( 𝐸, 𝑧, 𝑡 ) 𝑄 𝜇 ( 𝐸, 𝑧, 𝑡 ) ª®¬ (15)Note that, in the latter system of coupled differential equations wewrite down only three particles (protons, pions and muons), with-out considering the complete reaction family tree, since we haven’tdistinguished particles with different charge as 𝜋 ± and 𝜇 ± , and wehaven’t considered the possibility of left-right symmetry, 𝜇 𝐿,𝑅 ofmuons as we have done below in Sect. 4. This means that, by con-sidering all these particles, the system of Eq. (15) will have sevenlines.Usually, we ignore the term involving second derivative with re-spect to 𝐸 means that the acceleration term is negligible.The solution of the general transport equation for particles of onekind (when the last term of this equation can be omitted) is describedin detail in Ref. (Kosmas & Smponias 2018). In the latter work,however, neutrinos produced only from 𝜋 ± have been considered. Inthe present paper we proceed further and include neutrino emissionsalso from the muon decays, so we need to solve the system of Eq.(15) by including the third line too.The complete form of the system of Eq. (15) is treated with themethod of Ref. (Tsoulos, Kosmas & Stavrou 2017) in order to findthe exact solutions. In Ref. (Papadopoulos, Papavasileiou & Kosmas2020), equation (15) is treated semi-analytically. In the present work,we restrict ourselves to simplified forms resulting by neglecting var-ious phenomena (processes) taking place inside the astrophysicaloutflows (jet plasma), even though some of the assumed omissionsmay be considered as rather crud approximations.Below we discuss the cases of Eq. (15) satisfying the con-ditions of the one-zone approximation (Khangulyan et al. 2007;Kosmas & Smponias 2018), i.e. the cases when the particle distribu-tions are independent of time (steady state approximation). At first, in the calculations of our present work, we start by writingthe simplest solutions of the system of transport equations obtainedby assuming that all energy losses are absent i.e. 𝑏 𝑗 =
0. Then, by de-noting 𝑁 𝑝, , 𝑁 𝜋, and 𝑁 𝜇, the corresponding energy distributionsfor protons, pions and muons, respectively, the system of transport equations takes the trivial form © « 𝑡 − 𝑝 ( 𝐸 ) 𝑁 𝑝, ( 𝐸, 𝑧 ) 𝑡 − 𝜋 ( 𝐸 ) 𝑁 𝜋, ( 𝐸, 𝑧 ) 𝑡 − 𝜇 ( 𝐸 ) 𝑁 𝜇, ( 𝐸, 𝑧 ) ª®®¬ = © « 𝑄 𝑝 ( 𝐸, 𝑧 ) 𝑄 𝜋 ( 𝐸, 𝑧 ) 𝑄 𝜇 ( 𝐸, 𝑧 ) ª®¬ (16)The calculated distributions 𝑁 𝑗, , with j=p, 𝜋 ± , 𝜇 ± of the latterequations are discussed in the next Section. Under the conditions of the steady state approximation, and assumingthat the various energy losses are absent, the system of transportequations Eq. (15) takes the form 𝜕𝜕𝐸 © « 𝑏 𝑝 ( 𝐸 ) 𝑁 𝑝 ( 𝐸, 𝑧 ) 𝑏 𝜋 ( 𝐸 ) 𝑁 𝜋 ( 𝐸, 𝑧 ) 𝑏 𝜇 ( 𝐸 ) 𝑁 𝜇 ( 𝐸, 𝑧 ) ª®¬ + © « 𝑡 − 𝑒𝑠𝑐 𝑁 𝑝 ( 𝐸, 𝑧 ) 𝑡 − 𝜋 𝑁 𝜋 ( 𝐸, 𝑧 ) 𝑡 − 𝜇 𝑁 𝜇 ( 𝐸, 𝑧 ) ª®®¬ = © « 𝑄 𝑝 ( 𝐸, 𝑧 ) 𝑄 𝜋 ( 𝐸, 𝑧 ) 𝑄 𝜇 ( 𝐸, 𝑧 ) ª®¬ (17)In order to calculate the time independent neutrino and gamma-ray emissivities, we need, first, to calculate the distributions ofprotons, pions and muons, 𝑁 𝑝 , 𝑁 𝜋 , 𝑁 𝜇 , respectively, from Eq.(17). For these computations we used a code written in the Cprogramming language, mainly following the assumptions of Refs.(Reynoso, Romero & Christiansen 2008; Reynoso & Romero 2009;Vieyro & Romero 2017).In the latter case, the second term of the l.h.s. of the transportequation we replace the escape rate, 𝑡 − 𝑒𝑠𝑐 of Eq. (10) and 𝑏 ( 𝐸 ) = − 𝐸𝑡 − 𝑙𝑜𝑠𝑠 ( 𝐸 ) . The latter quantity the energy loss rate for each particlewhich is given by 𝑏 𝑝 ( 𝐸 ) = − 𝐸𝑡 − 𝑙𝑜𝑠𝑠,𝑝 ( 𝐸 ) = − 𝐸 ( 𝑡 − 𝑠𝑦𝑛 + 𝑡 − 𝑎𝑑 + 𝑡 − 𝑝𝑝 ) , (18) 𝑏 𝜋 ( 𝐸 ) = − 𝐸𝑡 − 𝑙𝑜𝑠𝑠, 𝜋 ( 𝐸 ) = − 𝐸 ( 𝑡 − 𝑠𝑦𝑛 + 𝑡 − 𝑎𝑑 + 𝑡 − 𝜋 𝑝 ) , (19) 𝑏 𝜇 ( 𝐸 ) = − 𝐸𝑡 − 𝑙𝑜𝑠𝑠,𝜇 ( 𝐸 ) = − 𝐸 ( 𝑡 − 𝑠𝑦𝑛 + 𝑡 − 𝑎𝑑 ) (20)It is worth noting that, in the latter equations some additional termsmay appear in cases when the rates of some other processes, ignored MNRAS , 1–11 (2020) in our present work (like e.g. the Inverse Compton Scattering, thepion-muon scattering, etc.), may be considered important for thedescription of other cosmic structures.
In this section, the source functions 𝑄 𝑗 ( 𝐸, r , 𝑡 ) entering the r.h.sof the system of coupled integro-differential equations of transporttype, are discussed. Towards this aim, initially in the system of Eqs.(15), we insert phenomenological expressions obtained from Ref.(Lipari, Lusignoli & Meloni 2009). A mathematical semi-analyticway of solving this system is proposed in Ref. (Kosmas 2020). Also,following Ref. (Tsoulos, Kosmas & Stavrou 2017), a method is de-rived for the numerical solution of the differential system of transportequations. In the assumed p-p mechanism, the density of fast (relativistic) pro-tons injected is more important, while the corresponding density ofslow (thermal) protons is dynamically significant.A usual injection function for the relativistic protons, coming outof the acceleration mechanism, has been taken to be a power-lawwith exponent equal to two, i.e. a function of the energy of the form(Ghisellini, Maraschi & Treves 1985) 𝑄 𝑝 ( 𝐸, 𝑧 ) = 𝑄 (cid:18) 𝑧 𝑧 (cid:19) 𝐸 , (21)where 𝑄 is a normalization constant obtained by specifying thepower in the relativistic protons (Reynoso, Romero & Christiansen2008), see Appendix. The above form is valid for the jet’s frameof reference (the corresponding expression in observer’s system isshown in Appendix).It should be noted that, Eq. (21) enters (through integration) thedetermination of the injection functions of the secondary particles( 𝜋 ± , 𝜇 ± , see below) and this justifies the term systems of coupledintegro-differential equations used in Sect. 3. The latter injectionfunctions are also dependent on the rates and cross sections of thereactions that preceded in the chain processes following the inelasticp-p collision. Under the assumptions discussed before, the system of transportequations (17) can be easily solved and the obtained proton distribu-tion 𝑁 𝑝 ( 𝐸, 𝑧 ) is written as (Reynoso, Romero & Christiansen 2008;Romero & Vila 2008) 𝑁 𝑝 ( 𝐸, 𝑧 ) = ∫ 𝐸 𝑚𝑎𝑥𝑝 𝐸 𝑝 | 𝑏 𝑝 ( 𝐸 ) | − 𝑄 𝑝 ( 𝐸 ′ , 𝑧 ) 𝑒 − 𝜏 𝑝 ( 𝐸,𝐸 ′ ) 𝑑𝐸 ′ (22)where 𝜏 𝑝 ( 𝐸, 𝐸 ′ ) = ∫ 𝐸 ′ 𝐸 𝑡 − 𝑒𝑠𝑐,𝑝 𝑑𝐸 ′′ | 𝑏 𝑝 ( 𝐸 ′′ ) | . (23)The quantity 𝑄 𝑝 ( 𝐸, 𝑧 ) corresponds to the relativistic injection func-tion of protons at the observer’s frame (see Appendix). The minimumenergy of protons is 𝐸 𝑚𝑖𝑛𝑝 = 𝐸 𝑡ℎ = .
22 GeV, while the maximumenergy is assumed to be 𝐸 𝑚𝑎𝑥𝑝 = − GeV. N [ G e V - c m - ] E [GeV]
Proton distribution
M33 X-7SS 433
Figure 4.
Proton energy distribution for M33 X-7 and SS 433 microquasarsystems.
Figure 4 shows the proton energydistribution 𝑁 𝑝 ( 𝐸 ) at the baseof the jet for M33 X-7 (solid line) and SS 433 (dashed line). As canbe seen, the number of protons appears reduced (even by two ordersof magnitude) for the wide half-opening angle system of M33 X-7( 𝜉 = ◦ ) as compared to the narrow of SS 433 ( 𝜉 = . ◦ ). This isbecause, the jet radius 𝑟 ( 𝑧 ) is bigger and hence the magnetic field 𝐵 is smaller in the case of M33 X-7 than that of SS433. For pions produced through the inelastic p-p scattering, as injectionfunction we adopt the one given by (Kelner, Aharonian & Bugayov2006) as 𝑄 𝜋 ( 𝐸, 𝑧 ) = 𝑛 ( 𝑧 ) 𝑐 ∫ 𝜀 𝜎 𝑖𝑛𝑒𝑙𝑝𝑝 (cid:18) 𝐸𝑥 (cid:19) 𝑁 𝑝 (cid:18) 𝐸𝑥 , 𝑧 (cid:19) 𝐹 𝜋 (cid:18) 𝑥, 𝐸𝑥 (cid:19) 𝑑𝑥𝑥 (24)where 𝜀 = 𝐸 / 𝐸 𝑚𝑎𝑥𝑝 , 𝑥 = 𝐸 / 𝐸 𝑝 and 𝐹 𝜋 denotes the distribution ofpions produced per 𝑝 − 𝑝 collision (see Appendix). The steady state energy distribution 𝑁 𝜋 ( 𝐸, 𝑧 ) for pions (createdthrough the scattering of hot protons off the cold ones), resultsthrough a similar way to that of Eqs. (22) and (23). With the re-placement 𝑝 → 𝜋 and 𝑡 − 𝑒𝑠𝑐 → 𝑡 − 𝜋 ( 𝐸, 𝑧 ) , respectively, on the latterequations the corresponding solution for pion energy distribution is,then, written as 𝑁 𝜋 ( 𝐸, 𝑧 ) = ∫ 𝐸 𝑚𝑎𝑥 𝐸 | 𝑏 𝜋 ( 𝐸 ) | − 𝑄 𝜋 ( 𝐸 ′ , 𝑧 ) 𝑒 − 𝜏 𝜋 𝑑𝐸 ′ (25)where the rate 𝑡 − 𝜋 includes contributions from the decay and escaperates [see Eq. (11)]. As mentioned before, in this work in addition to neutrinos comingfrom the decay of charged pions 𝜋 ± , we evaluate neutrino emissivitiesoriginating from the decay of the secondary charged muons 𝜇 ± .For the latter emissivity, we follow Ref. (Lipari, Lusignoli & Meloni2009), in order to take into account properly the muon energy loss.This means that, it is necessary to consider the production of bothleft handed and right handed muons, 𝜇 − 𝐿 and 𝜇 + 𝑅 , separately, because 𝜇 − 𝐿 and 𝜇 + 𝑅 have different decay spectra. MNRAS000
Figure 4 shows the proton energydistribution 𝑁 𝑝 ( 𝐸 ) at the baseof the jet for M33 X-7 (solid line) and SS 433 (dashed line). As canbe seen, the number of protons appears reduced (even by two ordersof magnitude) for the wide half-opening angle system of M33 X-7( 𝜉 = ◦ ) as compared to the narrow of SS 433 ( 𝜉 = . ◦ ). This isbecause, the jet radius 𝑟 ( 𝑧 ) is bigger and hence the magnetic field 𝐵 is smaller in the case of M33 X-7 than that of SS433. For pions produced through the inelastic p-p scattering, as injectionfunction we adopt the one given by (Kelner, Aharonian & Bugayov2006) as 𝑄 𝜋 ( 𝐸, 𝑧 ) = 𝑛 ( 𝑧 ) 𝑐 ∫ 𝜀 𝜎 𝑖𝑛𝑒𝑙𝑝𝑝 (cid:18) 𝐸𝑥 (cid:19) 𝑁 𝑝 (cid:18) 𝐸𝑥 , 𝑧 (cid:19) 𝐹 𝜋 (cid:18) 𝑥, 𝐸𝑥 (cid:19) 𝑑𝑥𝑥 (24)where 𝜀 = 𝐸 / 𝐸 𝑚𝑎𝑥𝑝 , 𝑥 = 𝐸 / 𝐸 𝑝 and 𝐹 𝜋 denotes the distribution ofpions produced per 𝑝 − 𝑝 collision (see Appendix). The steady state energy distribution 𝑁 𝜋 ( 𝐸, 𝑧 ) for pions (createdthrough the scattering of hot protons off the cold ones), resultsthrough a similar way to that of Eqs. (22) and (23). With the re-placement 𝑝 → 𝜋 and 𝑡 − 𝑒𝑠𝑐 → 𝑡 − 𝜋 ( 𝐸, 𝑧 ) , respectively, on the latterequations the corresponding solution for pion energy distribution is,then, written as 𝑁 𝜋 ( 𝐸, 𝑧 ) = ∫ 𝐸 𝑚𝑎𝑥 𝐸 | 𝑏 𝜋 ( 𝐸 ) | − 𝑄 𝜋 ( 𝐸 ′ , 𝑧 ) 𝑒 − 𝜏 𝜋 𝑑𝐸 ′ (25)where the rate 𝑡 − 𝜋 includes contributions from the decay and escaperates [see Eq. (11)]. As mentioned before, in this work in addition to neutrinos comingfrom the decay of charged pions 𝜋 ± , we evaluate neutrino emissivitiesoriginating from the decay of the secondary charged muons 𝜇 ± .For the latter emissivity, we follow Ref. (Lipari, Lusignoli & Meloni2009), in order to take into account properly the muon energy loss.This means that, it is necessary to consider the production of bothleft handed and right handed muons, 𝜇 − 𝐿 and 𝜇 + 𝑅 , separately, because 𝜇 − 𝐿 and 𝜇 + 𝑅 have different decay spectra. MNRAS000 , 1–11 (2020) -6 -4 -2 N [ G e V - c m - ] E [GeV]
Pion distribution (M33 X-7) N π N π ,0 -8 -6 -4 -2 N [ G e V - c m - ] E [GeV]
Muon distribution (M33 X-7) N µ N µ ,0 Figure 5.
Pion (Top) and muon (bottom) distributions for M33 X-7. In bothcase of particles-j, we compare the distribution which considers the importantenergy losses 𝑁 𝑗 ( 𝐸, 𝑧 = 𝑧 ) and with that which neglects energy losses 𝑁 𝑗, ( 𝐸, 𝑧 = 𝑧 ) for 𝑗 = 𝜋, 𝜇 . Both types of particle distributions refer tothe base of the jet, i.e. at 𝑧 = 𝑧 Thus, the injection functions of the left handed and right handedmuons are (Lipari, Lusignoli & Meloni 2009): 𝑄 𝜇 − 𝐿 ,𝜇 + 𝑅 ( 𝐸 𝜇 , 𝑧 ) = ∫ 𝐸 𝑚𝑎𝑥 𝐸 𝜇 𝑡 − 𝑑𝑒𝑐, 𝜋 ( 𝐸 𝜋 ) 𝑁 𝜋 ( 𝐸 𝜋 , 𝑧 ) Θ ( 𝑥 − 𝑟 𝜋 )× 𝑟 𝜋 ( − 𝑥 ) 𝐸 𝜋 𝑥 ( − 𝑟 𝜋 ) 𝑑𝐸 𝜋 (26)and 𝑄 𝜇 − 𝑅 ,𝜇 + 𝐿 ( 𝐸 𝜇 , 𝑧 ) = ∫ 𝐸 𝑚𝑎𝑥 𝐸 𝜇 𝑡 − 𝑑𝑒𝑐, 𝜋 ( 𝐸 𝜋 ) 𝑁 𝜋 ( 𝐸 𝜋 , 𝑧 ) Θ ( 𝑥 − 𝑟 𝜋 )× ( 𝑥 − 𝑟 𝜋 ) 𝐸 𝜋 𝑥 ( − 𝑟 𝜋 ) 𝑑𝐸 𝜋 (27)with 𝑥 = 𝐸 𝜇 / 𝐸 𝜋 and 𝑟 𝜋 = ( 𝑚 𝜇 / 𝑚 𝜋 ) . The corresponding muon energy distributions, coming out of thetransport equation, result from Eq. (25) via the replacement 𝑡 − 𝑒𝑠𝑐 → 𝑡 − 𝜇 ( 𝐸, 𝑧 ) , for muons and are then written as 𝑁 𝜇 𝑖 ( 𝐸, 𝑧 ) = ∫ 𝐸 𝑚𝑎𝑥 𝐸 | 𝑏 𝜇 ( 𝐸 ) | − 𝑄 𝜇 𝑖 ( 𝐸 ′ , 𝑧 ) 𝑒 − 𝜏 𝜇 𝑑𝐸 ′ (28)where the rate 𝑡 − 𝜇 for muons is given in Eq. (11).In Fig. 5, the pion (top) and muon (bottom) energy distributionsfor M33 X-7 are illustrated. In both cases of particles, the distribution 𝑁 𝑗 ( 𝐸, 𝑧 = 𝑧 ) , with 𝑗 = 𝜋 or 𝜇 , which considers the most importantenergy losses, is compared with that which neglects energy losses, 𝑁 𝑗, ( 𝐸, 𝑧 = 𝑧 ) . The solid lines correspond to distributions obtainedby considering energy losses and the dashed lines correspond to thosewhere the energy losses obtained have been neglected. 𝛾 -RAY EMISSIVITIES In this Section, we present and discuss simulated emissivities of highenergy neutrinos and 𝛾 -rays produced in the extragalactic micro- quasar M33 X-7 system which is located in the Messier 33 galaxy,at a distance ∼ −
960 kpc from the Earth (Pietsch et al. 2006).The derived algorithms for this purpose have been tested on thewell-studied Galactic micro-quasars SS 433 and Cyg X-1 system aswell.Because in our previous calculations on neutrino production fromMQs (Smponias & Kosmas 2015, 2017; Kosmas & Smponias 2018),we neglected (due to the complexity and long time consuming) emis-sivity (intensity) of neutrinos originated from the secondary muons( 𝜇 ± ) produced from the charged pion ( 𝜋 ± ) decays, in this section, wepresent contributions also originating from the 𝜇 ± channel [see Eq.(4)].By using the concentrations for protons 𝑁 𝑝 , pions 𝑁 𝜋 and muons 𝑁 𝜇 , the neutrino and 𝛾 -ray intensities are subsequently calculated asdescribed below. After the above discussion, the total emissivity 𝑄 𝜈 ( 𝐸, 𝑧 ) producedfrom a MQ is the sum of contributions from the two sources: (i)the first comes from the direct 𝜋 ± decay (prompt neutrino produc-tion), and (ii) the second comes from the 𝜇 ± decay (delayed neutrinoproduction). Thus, 𝑄 𝜈 ( 𝐸, 𝑧 ) = 𝑄 𝜋 → 𝜈 ( 𝐸, 𝑧 ) + 𝑄 𝜇 → 𝜈 ( 𝐸, 𝑧 ) (29)For pion decays the injection function is given by 𝑄 𝜋 → 𝜈 ( 𝐸, 𝑧 ) = ∫ 𝐸 𝑚𝑎𝑥 𝐸 𝑡 − 𝑑𝑒𝑐, 𝜋 ( 𝐸 𝜋 ) 𝑁 𝜋 ( 𝐸 𝜋 , 𝑧 ) Θ ( − 𝑟 𝜋 − 𝑥 ) 𝐸 𝜋 ( − 𝑟 𝜋 ) 𝑑𝐸 𝜋 (30)with 𝑥 = 𝐸 / 𝐸 𝜋 , while for the four types of muon decays the injectionfunction reads 𝑄 𝜇 → 𝜈 ( 𝐸, 𝑧 ) = Í 𝑖 = ∫ 𝐸 𝑚𝑎𝑥 𝐸 𝑡 − 𝑑𝑒𝑐,𝜇 ( 𝐸 𝜇 ) 𝑁 𝜇 𝑖 ( 𝐸 𝜇 , 𝑧 )× h − 𝑥 + 𝑥 + (cid:16) 𝑥 − − 𝑥 (cid:17) ℎ 𝑖 i 𝑑𝐸 𝜇 𝐸 𝜇 . (31)(with 𝑥 = 𝐸 / 𝐸 𝜇 ).In the last expression, the symbols 𝜇 𝑖 , 𝑖 = , , , 𝜇 { , } = 𝜇 {− , +} 𝐿 , 𝜇 { , } = 𝜇 {− , +} 𝑅 (Lipari, Lusignoli & Meloni2009), while ℎ { , } = − ℎ { , } = −
1. Then, the calculation of eachof the latter two integrals of Eqs. (30) and (31) provides separatelythe partial emissivity of neutrinos for the prompted and the delayedneutrino source, respectively. Needless to note that Earth and spacetelescope are not able to discriminate the two source as it happenswith laboratory neutrino sources.Subsequently, we easily obtain the neutrino intensity (in units
𝐺𝑒𝑉 − 𝑠 − ) by the spatial integration 𝐼 𝜈 ( 𝐸 ) = ∫ 𝑉 𝑄 𝜈 ( 𝐸, 𝑧 ) 𝑑 𝑟 = 𝜋 ( tan 𝜉 ) ∫ 𝑧 𝑚𝑎𝑥 𝑧 𝑄 𝜈 ( 𝐸, 𝑧 ) 𝑧 𝑑𝑧 . (32)Figure 6 (left) shows the neutrino intensity produced at the base ofthe jet from direct decays of secondary pions and muons coming fromp-p collisions in the jets of M33 X-7 and SS 433, respectively. Ascan be seen, the number of produced neutrinos reduces significantlyfor energies 𝐸 > GeV, following the behavior of pion and muondistribution. We can also see that the neutrinos produced in M33 X-7are less than those produced in SS 433, due to the wider-half openingangle 𝜉 as we explained before.It should be noted that, to perform the above calculations we haveupdated and improved our codes to reduce the time consuming to areasonable level, so as the integral involved in Eq. (31) to be easilyobtained, since in going the step from Eq. (30) to Eq. (31) the timeconsuming increases rapidly. MNRAS , 1–11 (2020) I [ G e V - s - ] E [GeV]
Neutrino emission
M33 X-7SS 433 E [GeV] γ - ray emission M33 X-7SS 433
Figure 6.
Neutrino (left) and gamma-ray (right) intensity for M33 X-7 and SS433
As we have discussed in Sect. 2, the p-p collisions inside the jets,produce secondary neutral particles (pions, eta particle, etc.) thatdecay to give 𝛾 -rays. The dominant of these channels goes throughthe reaction 𝜋 → 𝛾 + 𝛾 . (33)For 𝐸 𝛾 > 𝐺𝑒𝑉 , we consider the 𝛾 -ray emissivity at adistance 𝑧 along the jets (in units 𝐺𝑒𝑉 − 𝑠 − ) to be given by(Reynoso, Romero & Christiansen 2008) as 𝑄 𝛾 = 𝑐 ∫ 𝐸 𝛾 / 𝐸 𝑚𝑎𝑥𝑝 𝜎 𝑖𝑛𝑒𝑙𝑝𝑝 (cid:18) 𝐸 𝛾 𝑥 (cid:19) 𝑁 𝑝 (cid:18) 𝐸 𝛾 𝑥 , 𝑧 (cid:19) 𝐹 𝛾 (cid:18) 𝑥, 𝐸 𝛾 𝑥 (cid:19) 𝑑𝑥𝑥 . (34)where 𝐹 𝛾 is the spectrum of the produced 𝛾 -rays(Kelner, Aharonian & Bugayov 2006, 2009) with energy 𝑥 = 𝐸 𝛾 / 𝐸 𝑝 for a primary proton energy 𝐸 𝑝 . The function 𝐹 𝛾 is given in Ap-pendix.Subsequently, the corresponding spectral intensity of 𝛾 -rays, canbe obtained from the following spatial integration over the jet’s vol-ume 𝑉 as 𝐼 𝛾 ( 𝐸 𝛾 ) = ∫ 𝑉 𝑄 𝛾 ( 𝐸 𝛾 , 𝑧 ) 𝑑 𝑟 = 𝜋 ( tan 𝜉 ) ∫ 𝑧 𝑚𝑎𝑥 𝑧 𝑄 𝛾 ( 𝐸 𝛾 , 𝑧 ) 𝑧 𝑑𝑧, (35)Figure 6 (right) shows the 𝛾 -ray intensity for energies 𝐸 𝛾 > 𝛾 − 𝑟𝑎𝑦 intensity reduces rathersteadily following the behavior of the proton distribution in bothcases of MQ systems. Again for the M33 X-7 it is lower mainly dueto the wider half-opening angle (i.e. weaker magnetic field) but alsodue to other parameters.In our ongoing calculations (Kosmas 2020), among others we in-clude results obtained by using the 3-D PLUTO astrophysical hy-drocode where the jets of the studied systems are approximatedas purely relativistic magnetohydrodynamical (RMHD) flows withthe magnetic field (assumed tangled) being dragged with the flow(Smponias & Kosmas 2014). In this work we address neutrino and 𝛾 -ray emissions from mi-croquasars and X-ray binary stars (XRBs) that consist of a stellarmass 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 𝑝 − 𝑝 scat-tering mechanism, i.e. the inelastic collision of relativistic protons ofthe jet with the thermal ones.Such high energy neutrinos and 𝛾 -rays are detectable by operat-ing terrestrial and space telescopes. Among the operating detectorsare the under ice IceCube (at the South Pole), the ANTARES andKM3Net (under the Mediterranean sea), and many others.In the near future, more accurate measurements of gamma-rayand neutrino fluxes by the coming observatories, such as the CTA(Acharya et al.2018), which is expected to shed more light on thenature and the emission sources. From the perspective of neutrinodetection, the addition of more years of data with continuous op-eration of IceCube will improve the sensitivity of the search forGalactic sources of cosmic neutrinos. Furthermore, the next gener-ation detection instrument, IceCube-Gen2, a substantial expansionof IceCube, will be 10 times larger. This next-generation neutrinoobservatory with five times the effective area of IceCube is expectedto improve the neutrino source search sensitivity by the same or-der (Aartsen et al. 2019). With higher neutrino statistics, identifyingGalactic sources will become more promising.On the other hand, theoretically, by modelling the solution ofthe system of coupled transport equations, we were able to performdetailed calculations for various processes taking place inside the jetsof Galactic and extra galactic system M33 X-7 assuming hadroniccontent in their jets. These cooling rates enter the proton, pion andmuon energy distributions through which one obtains neutrino and 𝛾 -ray intensities.Our near future calculations include, among others, multidimen-tional simulations using the PLUTO astrophysical hydrocode whichmay provide both radiative and dynamical description of the as-trophysical outflows within multi-scale and multi-messesger astro-physics. ACKNOWLEDGEMENTS
This research is co-financed (O.T.K.) by Greece and the EuropeanUnion (European Social Fund-ESF) through the Operational Pro-gramme “Human Resources Development, Education and LifelongLearning 2014- 2020" in the context of the project (MIS5047635).D.A.P. wishes to thank Prof. T.S. Kosmas for fruitful discussionsduring his stay in the Dept. of Physics, University of Ioannina.
MNRAS , 1–11 (2020) DATA AVAILABILITY
There is no data used in this paper.
REFERENCES
Aartsen M. G. et al., 2015 ApJ, 805, L5Aartsen M. G. et al., 2016, (IceCube) Phys.Rev.Lett., 117 no.7, 071801.Aartsen M. G. et al., 2018, (IceCube) Science, 361, 147Aartsen M. G. et al., 2019, (IceCube) arXiv:1911.02561v1 [astro-ph.HE]Aartsen M. G. et al., 2020, (IceCube Collaborations), Phys. Rev. Lett., 124,051103, arXiv:1910.08488 [astro-ph.HE].Actis M., Agnetta G., Aharonian F. et al., 2011, Exp. Astron., 32, 193Adrian-Martinez S. et al., 2016, (KM3Net Collaboration) J. Phys. G: Nucl.Part. Phys., 43, 084001, e-Print: arXiv: [astro-ph.IM] 1601.07459Albert A. et al., 2020, (ANTARES and IceCube Collaborations), ApJ, 892,9, arXiv:2001.04412 [astro-ph.HE]Begelman M. C., Hatchett S. P.,McKee C. F., Sarazin C. L., Arons J., 1980,ApJ, 238, 722Bonanos A. Z., Stanek K. Z. et al., 2006, ApJ, 652, 313Bosch-Ramon V., Romero G. E., Paredes J. M., 2006, A&A, 447, 263Charles P. A., Coe M. J., 2006, in Compact Stellar X-ray Sources (eds LewinW. H. G. & van der Klis, M.), Cambridge Univ. Press, Cambridge, UK,215–265Cherepashchuk A.M. et al., 2005, A&A, 437, 561Fabrika S., 2004, Astrophys. Space Phys. Rev., 12, 1Falcke H., Biermann P., 1995, A&A, 293, 665Gaisser T. K., 1990, Cosmic Rays and Particle Physics, Cambridge UniversityPress, CambridgeGhisellini G., Maraschi L., Treves A., 1985, A&A, 146, 204Gies D. R. et al., 2008, ApJ, 678, 1237Ginzburg V. L., Syrovatskii S. I., 1964, The origin of cosmic rays, PergamonPress Ltd.Kelner S. R., Aharonian F. A., Bugayov V. V., 2006 Phys. Rev. D, 74, 034018Kelner S.R., Aharonian F.A., Bugayov V. V., 2009, Phys. Rev. D, 79, 039901Khangulyan D. et al., 2007, MNRAS, 380, 320Körding E. G., Fender R. P., Migliari S., 2006, MNRAS, 369, 1451Kosmas O. T., Smponias T., 2018, Adv. High Energy Phys. 960296,arXiv:1808.00303 [astro-ph.HE].Kosmas O. T., 2020, in preparationLipari P., Lusignoli M., Meloni D., 2009, Phys. Rev. D, 75, 123005Long K. S., Dodorico S., Charles P. A., Dopita M. A., 1981, ApJ, 246, L61Mirabel I.F., Rodríguez L.F., 1999, ARA&A,37, 409Orosz J. A., 2003, in A Massive Star Odyssey: From Main Sequence toSupernova (eds van der Hucht K. A., Herrero A., Esteban C.), Proc. IAUSymp. 212, ASP, San Francisco, 365Orosz J.A., McClintock J.E., Narayan R. et al., 2007, Nature, 449, 872Orosz J.A. et al., 2011, ApJ, 742, 84Papadopoulos D. A., 2020, Study of gamma-ray and neutrino emission frommicroquasar jets, MSc Thesis, University of Ioannina (unpublished)Papadopoulos D.A., Papavasileiou Th.V., Kosmas T.S., 2020,J.Phys.Conf.Ser., to appear, [astro-ph.HE] arXiv:2010.00396.Pietsch W. et al., 2006, ApJ, 646, 420Reid M.J. et al., 2011, ApJ, 742, 83Remillard R. A., McClintock, 2006, ARA&A, 44, 49Reynoso M.M., Romero G.E., Christiansen H.R., 2008, MNRAS, 387, 1745Reynoso M.M., Romero G.E., 2009, A&A, 493, 1Romero G.E., Torres D.F., Kaufman Bernadó M.M., Mirabel I.F., 2003, A&A,410, L1Romero G. E., Vila G. S., 2008, A&A, 485, 623Romero G. E., Boettcher M., Markoff S., Tavecchio F., 2016, Space Sci. Rev.,207, 5Romney J.D. et al., 1987, ApJ, 321, 822Saito T. Y. et al. (for the MAGIC Collaboration), 2009, preprint(arXiv:0907.1017v2)Smponias T., Kosmas T. S., 2011, MNRAS, 412, 1320Smponias T., Kosmas T. S., 2014, MNRAS, 438, 1014 Smponias T., Kosmas O. T., 2015, Adv. High Energy Phys., 921757Smponias T., Kosmas O. T., 2017, Adv. High Energy Phys., 4962741,arXiv:1706.03087 [astro-ph.HE]Spencer R. E., 1984, MNRAS, 209, 869Tsoulos I. G., Kosmas O. T., Stavrou V. N., 2019, Comp. Phys. Commun.236, 237Vieyro F. L., Romero G.E., 2012, A&A, 542, A7, 1
APPENDIX A:A1 The hot proton’s injection function
At the observer’s frame of reference, the injection function of protons, 𝑄 𝑝 ( 𝐸, 𝑧 ) , and is given by (Reynoso, Romero & Christiansen 2008): 𝑄 𝑝 ( 𝐸, 𝑧 ) = (cid:18) 𝑧 𝑧 (cid:19) 𝑄 Γ 𝑏 (cid:16) 𝐸 − 𝛽 𝑏 cos 𝜃 √ 𝐸 − 𝑚 𝑐 (cid:17) × (cid:20) − 𝛽 𝑏 𝐸 𝑐𝑜𝑠𝜃 √ 𝐸 − 𝑚 𝑐 (cid:21) (A1)where Γ 𝑏 is the bulk Lorentz factor of the jet, and 𝜃 denotes the anglebetween the jet’s injection axis and the direction of the line of sight(LOS).Then, the normalization constant 𝑄 is obtained bydetermining the power 𝐿 𝑝 of the relativistic protons(Reynoso, Romero & Christiansen 2008), given by 𝐿 𝑝 = ∫ 𝑉 𝑑 𝑟 ∫ 𝐸 𝑚𝑎𝑥𝑝 𝐸 𝑚𝑖𝑛𝑝 𝐸 𝑝 𝑄 𝑝 ( 𝐸 𝑝 , 𝑧 ) 𝑑𝐸 𝑝 (A2)After performing the latter integration, the normalization constant 𝑄 (see text) reads (Reynoso, Romero & Christiansen 2008) 𝑄 = 𝑐𝑧 𝐾 , 𝐾 = 𝑞 𝑟𝑒𝑙 𝐿 𝑘 𝑐𝑟 𝑙𝑛 (cid:16) 𝐸 ′ 𝑚𝑎𝑥𝑝 / 𝐸 ′ 𝑚𝑖𝑛𝑝 (cid:17) . (A3) A2 The distribution of pions produced per p-p collision
The distribution of pions produced per p-p collision is 𝐹 ( 𝑝𝑝 ) 𝜋 (cid:16) 𝑥. 𝐸𝑥 (cid:17) = 𝛼𝐵 𝜋 𝑥 𝛼 − (cid:16) − 𝑥 𝛼 + 𝑟 𝑥 𝛼 ( − 𝑥 𝛼 ) (cid:17) × (cid:16) − 𝑥 𝛼 + 𝑟 ( − 𝑥 𝛼 ) + 𝑟 𝑥 𝛼 ( − 𝑥 𝛼 ) (cid:17) (cid:16) − 𝑚 𝜋 𝑐 𝑥𝐸 𝑝 (cid:17) (A4)with 𝑥 = 𝐸 / 𝐸 𝑝 , 𝐵 𝜋 = 𝛼 ′ + . 𝛼 ′ = . + . 𝐿 + . 𝐿 , 𝑟 = . /√ 𝛼 ′ , and 𝛼 = . /√ 𝛼 ′ (Kelner, Aharonian & Bugayov2006, 2009). A3 Produced 𝛾 -ray spectrum Following the treatment of (Kelner, Aharonian & Bugayov 2006,2009), the spectrum of produced gamma-rays with energy 𝑥 = 𝐸 𝛾 / 𝐸 𝑝 for a primary proton with energy 𝐸 𝑝 is written as 𝐹 𝛾 ( 𝑥, 𝐸 𝑝 ) = 𝐵 𝛾 𝑙𝑛𝑥𝑥 " − 𝑥 𝛽 𝛾 + 𝑘 𝛾 𝑥 𝛽 𝛾 ( − 𝑥 𝛽 𝛾 ) (A5) × " 𝑙𝑛𝑥 − 𝛽 𝛾 𝑥 𝛽 𝛾 − 𝑥 𝛽 𝛾 − 𝑘 𝛾 𝛽 𝛾 𝑥 𝛽 𝛾 ( − 𝑥 𝛽 𝛾 ) + 𝑘 𝛾 𝑥 𝛽 𝛾 ( − 𝑥 𝛽 𝛾 ) MNRAS , 1–11 (2020) where 𝐵 𝛾 = . + . 𝐿 + . 𝐿 ,𝛽 𝛾 = . 𝐿 + . 𝐿 + . , (A6) 𝑘 𝛾 = . 𝐿 + . 𝐿 + . 𝐿 = 𝑙𝑛 (cid:0) 𝐸 𝑝 /( 𝑇 𝑒𝑉 ) (cid:1) (A7) This paper has been typeset from a TEX/L A TEX file prepared by the author. MNRAS000