High Energy Neutrino Emission from Astrophysical Jets in the Galaxy

We address simulated neutrino emission originated from astrophysical jets of compact objects within the Galaxy. These neutrinos are of high energies (E] of the order up to a few TeV) and for their observation specialized instruments are in operation, both on Earth and in orbit. Furthermore, some next generation telescopes and detector facilities are in the process of design and construction.The jet flow simulations are performed using themodern PLUTOhydrocode in its relativistic magnetohydrodynamic version. One of the main ingredients of the present work is the presence of a toroidal magnetic field that confines the jet flow and furthermore greatly affects the distribution of the high energy neutrinos.


Introduction
In recent years, a remarkable development of neutrino and gamma-ray astronomy took place.New telescopes and detector facilities were designed and implemented, such as the CTA Cerenkov gamma-ray array [1], Fermi orbital telescope [2], or the IceCube neutrino detector [3].The IceCube neutrino observatory, located at the South Pole, detects neutrinos in a wide energy range, far exceeding the energies of manmade accelerator beams.IceCube features a large detector volume, increasing the possibility to detect neutrinos from individual astrophysical sources.
Potential sites for gamma-ray and neutrino production from jets include galactic sources of X-ray binaries (XRB) [4][5][6], as well as extragalactic sources [7] such as active galactic nuclei (AGN) ( [8]).The above categories may be extended to include a wide range of different phenomena, such as supernova remnants [9,10] and other possible sources of high energy particles.These systems may produce -rays and high energy neutrinos from interactions (collisions) of high energy protons with thermal ones.In general, a large detector volume is needed, since neutrinos are so weakly interacting with matter.
At the same time, -ray emission can also be observed, from such sources, by ground and satellite based gammaray telescopes (e.g., Fermi, CTA).When studying acceleration processes of those sources, it is often useful to compare the neutrino and -ray fluxes emanating from them.
XRB are binary stellar systems comprising a main sequence star and a compact object, emitting in the X-ray band, often also called microquasars [4,11], in relation to the their cousins of galactic scale, quasars.They have relativistic jets that include acceleration sites allowing particles to reach energies up to TeV [12].There, the cooling of TeV electrons is very strong due to the high density of the radiation, implying a hadronic mechanism for the gamma-ray emission.Simultaneous neutrino emission is then possible as well [5,6].
In this work, we aim to simulate the neutrino flux generated from a model microquasar X-ray binary system, assuming a hadronic jet and hadron-related acceleration processes (nonthermal proton acceleration and interactions of high energy protons).Towards this purpose, a power law with an exponential cutoff spectrum is used [5,13].This is an extension of our previous calculations of -ray emission from a simulated relativistic hydrodynamical jet [14].Neutrinos are mainly produced through pion and muon decay with the pions coming from inelastic -scatterings among nonthermal protons and thermal ones within the jet.Pion decay by-products include a muon and a muon neutrino as  ± →  ± + ]  (]  ).The muons can afterwards decay again into an electron or a positron and the associated neutrino.Another pion decay channel leads to two gamma-ray photons [6,13,15]: (  = 1.67 × 10 −24 g and   = 2.38 × 10 −25 g). comprises  0 and  +  − pairs.Subsequently, pion decay leads to muons, -rays, and neutrino production.
In the present work, we concentrate on simulations of neutrino emission from relativistic galactic astrophysical jets.These neutrinos have very high energy (of the order of TeV) and specialized instruments are in operation in order to detect them, both on Earth and in space.Furthermore, some next generation instruments are in the process of design and construction.The jet flow simulations are performed using the modern PLUTO hydrocode in its relativistic magnetohydrodynamic version [16].We assume the existence of a toroidal magnetic field that confines the jet flow (in the region of the compact object) [17], affecting the production of the high energy neutrinos [6].

Description of the Method
The jet is modelled using the relativistic magnetohydrodynamic (RMHD) version of the PLUTO hydrocode [16].The PLUTO RMHD module is employed in order to simulate the jet flow.The Generalized Lagrange Multiplier (GLM) correction method is used, enforcing magnetic divergence suppression through hyperbolic divergence cleaning, while the MUSCL-Hancock scheme is employed as the integrator.A toroidal magnetic field is used that helps constrain the jet, to an extent depending on the field strength.In this paper, we explore the existence of a strong magnetic field that keeps the emitted neutrino flux highly concentrated [17].
The toroidal magnetic field setup employed leads to a rather pronounced jet confinement, due to Lorentz forces acting on the jet matter towards the jet axis.The stellar wind is set to decrease away from the companion star as 1/ 2 , while a corona of 1/ 2 ,  being the jet axis direction, is setup near the compact object, similar to [14].The most important model parameters are shown in Table 1.The IDL suite and the VisIt visualization suite are then employed in order to present the results of the simulations in a graphical manner.
The boundary conditions are outflow at the top and at the sides of the computational domain ("box") and reflective at the bottom, where the jet base is located.The jet emanates from the middle of the bottom plane (- plane) moving upwards, that is, along the model's -axis.

Radiative Transfer and Imaging
The line-of-sight (LOS) code of [18] is used to produce artificial neutrino "images" of the jet-corona system.The particle "emission" is calculated separately from each computational cell, under the twin assumptions that the particle distribution's dominant cooling time is smaller than the hydrocode's time step and also that the mean free path between collisions is smaller than the hydrocode computational cell dimensions.Then each cell can be treated individually from a radiative point of view and its neutrino emission contribution, at a given energy level, is added to its corresponding LOS.Naturally, no absorption is relevant for neutrinos in such systems.
The emission calculation is performed in Mathematica, mainly following the analysis of [5,6].Furthermore, for converting the neutrino emission from the jet reference frame to our rest frame, the calculational procedure can be, for example, that of [19] or that of [20], and we have, respectively, where Ω is the line-of-sight solid angle, Γ denotes the jet beam Lorentz factor (for thermal, slow protons), and  = /, the known ratio.
The difference between the above expressions lies in the powers of the two factors of the denominators on their last fraction.In the present work, we have chosen to apply the expression of [19].
The nonthermal proton distribution suffers synchrotron and adiabatic losses, affecting the balance in the transport between protons and pions.Following the formalism of [5,13], we may then obtain the neutrino emissivity: where  =   /  ,  = 0.427, and the   (for pion distribution) reads The function  (inj)    is given in [13] and   ,   are given in [5].We note the dependency of   on the magnetic field , among others.
On the other hand, adopting the approach of [6], also [15], we have the following line of calculations for the neutrino emission that was actually used in the current work.For the readers convenience we give here a brief presentation of the formalism of [6] (for more details see [6] and references therein).
The pion energy distribution is provided as the solution of a transport equation where Finally, for the emissivity of neutrinos emanating from direct pion decays (prompt neutrinos), we have [6,15] where  = /  and  ,dec is the pion decay time scale.The neutrino emissivity can then be integrated over 3D cell volume (voxel volume) and divided by the surface of a sphere, whose radius is the distance to Earth.The result is a synthetic "neutrino emission observation" of the binary system.By repeating the process for many energies, we can then obtain a SED plot.In order to simplify the computations, the hydrodynamic quantities will first be obtained as an average over the whole of the jet and then used to calculate the point emission from the jet system, at a given energy (see below).component (Figure 1), as opposed to the null RHD case in [14] (Figures 2 and 3).The various emission sites are then concentrated along the jet axis.This effect is even more pronounced because of local turbulence occurring naturally in the jet flow.The jet head advances through the surrounding stellar and corona winds, but its sideways expansion is limited.The magnetic field contributes to the emission mechanism so increased emission is expected from the inner flow funnel of the jet structure.The strong confinement is a result of adoption of a relatively intense field configuration, allowing the beam to remain focused late into the simulation; therefore emission is expected to be present at late stages of the numerical experiment.

Results and Discussion
The LOS code can create a synthetic image and we then add up all of its pixel's intensities for total neutrino intensity from the system, always at a given energy.By repeating the process for many different energies, we can obtain the spectral emission distribution of the model system.Using some kind of normalization in relation to an external factor, such as an energy estimate for the total energy emitted from the system in neutrinos, we could then adopt actual units and directly compare to past and planned observations.Nevertheless, we omit the LOS step and simply add up the intensities of neutrino emission from all volume cells of the system.The reasons are twofold.On the one hand, no absorption occurs for neutrinos; therefore no need exists for a full solution of the equation of radiative transfer.On the other hand, current and near-future observations detect stellar neutrino sources as strictly point ones so no spatial resolution whatsoever exists for them.Combining the above two arguments we arrive at the conclusion that merely adding up emissions from all cells of the system provides an accurate model neutrino emission from a point source system, all the while the system is being modelled internally as a fully dynamical relativistic magnetohydrodynamical jet (RMHD jet).In Figure 1 a snapshot of the 3-dimensional jet magnetic field configuration can be seen, as the jet head advances through the middle of the computational domain.The two main components of the field can be determined and the presence of the toroidal component contributes to the jet confinement through the Lorentz force towards the jet axis.There is a strong dependence of the neutrino emission on the magnetic field, since the synchrotron energy loss mechanism's time scale strongly depends on the value of the field .
In Figure 2 we can see the magnetic field magnitude, roughly halfway into the simulation run.The jet remains well confined, because of the magnetic force towards the jet axis acting due to the toroidal field component.This mechanism intensifies the possible neutrino emission, as it allows the magnetized matter to stay dense and also the field contributes to raise the emission levels (equation ( 4)).
Figure 3 shows the jet density at the middle of the simulation, where the jet is well confined and relatively little mixing has occurred between the jet matter and the surrounding winds.The sideways flow is still of secondary importance and therefore less emission is expected from the sides of the jet.This is opposed to the purely hydrodynamical jet case (e.g., [14]) where the jet expands much more and mixes with the ambient medium, allowing for more dynamical effects to occur over a larger volume at the jet sides.On the other hand, the magnetized jet here demonstrates a denser inner flow that stays focused and does not dissipate into the winds, maintaining stronger emission throughout the simulation.We can see the jet beam advancing along a narrow path through the ambient medium, kept together by the magnetic toroidal field component.This confined jet is then a suitable site for increased high energy -ray and neutrino emission/production, due to maintaining extreme conditions along the jet length.

Neutrino Emission Simulations.
In this work we simulate the neutrino emission from the SS433 microquasar system by extending the models of [5,14] which refer to the emission of -rays from these jets.The differences between the two types of emission necessitate the employment of additional techniques.Towards this aim we have used the following: (i) the formalism of [5,6] which consider the cascade of particles in microquasar jets that lead to the production of neutrinos; (ii) the computationally demanding numerical integration techniques used to reproduce the cascade of particles, at individual hydrocode grid cells, based on the works of [6,13].
Due to the relatively limited availability of computing time the grid resolution of the hydrocode data was reduced to a more manageable size.Nevertheless, we still obtain a better resolution than that currently available from neutrino "observations" of microquasar systems, where only point sources are detected.
We should mention that, in studying the dependence of neutrino intensity originating from our model system, the PLUTO code produces jets whose dynamics directly affect the shape of the curve in the high energy tail (well above 100 GeV).More specifically, the fact that the inner jet has a higher proportion of faster protons and pions, relative to the outer jet, causes the neutrino emissivity there to begin dropping abruptly at higher energies compared to the behavior shown by the results of [6].A more detailed analysis including muonic neutrino emissions and higher hydrodynamical grid resolution is currently performed and will appear in future work.
Figure 4 shows the unnormalized neutrino intensity on Earth,  ] ( ] ) per energy interval, coming out of the jet simulations.For simplicity, only promptly generated neutrinos from pion decay have been considered in this figure (delayed neutrinos [15] are to be included in future work).We concentrate on emissions produced not long after the beginning of the jet ejection event (consequently, our model emission here mainly comes from the inner part of the jet).
As a first approximation, the original hydrocode results of grid size 120 × 200 × 120 have been regridded into a much smaller resolution of 3 × 5 × 3.Then, following [6] the neutrino emission from  interactions was obtained.We can see a general agreement with [6].It should be stressed that the difference between the behaviour of   (  ) in the high energy tail illustrated in Figure 4 and that of [6] is mainly attributed to the strong magnetic confinement assumed in our present work, as compared to that used in [6].
In order to facilitate comparison with observations, the neutrino spectral intensity resulting from the jet is quantified using an energetic argument.At first, we assume that the energy fraction carried by the nonthermal (fast) proton distribution is of the order of 10 −4 [5], while the neutrino energy fraction, of the fast protons, carried away from the system is roughly of the order of 10% [21].In total, neutrinos are considered to carry away from the system around 10 −5 of the total kinetic luminosity of the bulk flow proton stream.The latter's kinetic luminosity is calculated, from the simulation initial conditions to be around 2 × 10 36 erg/sec.
In order to constrain the luminosity  ] ( ] ), we then calculate the area under the curve in Figure 4 equal to 10 −4 (arbitrary units) ×10 5 GeV = 10 (arbitrary units) ×GeV.This may represent the power emitted; therefore we equate this area to 2×10 31 erg/sec or 10 34 GeV/sec, which is the neutrino fraction of the jet kinetic luminosity.This fixes the arbitrary unit to be equal to 10 33 (GeV s) −1 .The intensity curve of Figure 4 flattens at around 10 −4 , which is therefore found to be 10 29 (GeV s) −1 .This value is lower than that of [5] by 1-2 orders of magnitude at lower energies but is higher by a similar amount at higher energies.The reason is that our jet is strongly confined by the toroidal magnetic field component which means that a higher energetic proton population is maintained later on.This is shown before the emission calculation binning and averaging process.Furthermore, this spectral emission distribution appears flatter at first and then drops rapidly.This favors emission at higher neutrino energies.
Detectability by current and upcoming arrays is consequently relatively worse than that in [5] (for the SS433 system's distance), which is taken to be marginal already.This implies that a more energetic jet is needed to provide higher neutrino flux, perhaps a jet of 100 times higher bulk mass flow rate.On the other hand, a more turbulent jet with more shock acceleration sites for nonthermal protons may lead to a  rel of 10 times higher, that is, 0.001.This raises the expected flux by a corresponding factor of 10 as well.We conclude that the model jet may, under favorable conditions, fall within the detection limits of modern detectors discussed in Section 1.

Summary and Conclusions
The neutrino production from a relativistic magnetohydrodynamic model jet was studied, using a hadronic model for proton-proton interactions, leading to pion decay.High energy neutrinos are assumed to be produced in galactic XRB that include a stellar mass compact object.The companion (donor) star is a main sequence one, at an earlier stage of its evolution.The binary system in general emits in many different wavelengths, from radio and IR to high energy gamma-rays and neutrinos.The primary engine driving the jets and their emissions is the gravitational attraction of matter into the compact object.Adiabatic and synchrotron losses were assumed and a steady-state radiative condition was employed at a given hydrocode time step (the radiative process is presumed to occur faster than the cell contents change dynamically).
The toroidal magnetic field configuration of the system confines the jet and creates an environment that favors high energy emission of -rays and production of TeV neutrinos in the jet.A first attempt to model the energetic neutrino production in the system is expected to yield a direct connection between the jet dynamics and the energy spectra of the particles.Further work may accommodate additional effects, such as the full use of hydrocode data into the radiative transfer calculations and the consideration of the relativistic nature of the imaged jet system.Finally, for the numerical solution of the partial differential equations (PDEs) arising from the problem, further, more accurate numerical schemes can be also considered, like the ones of [22].

4. 1 .Figure 1 :
Figure1: A 2-dimensional plot (jet cross-section, depicting a slice cut parallel to the -plane) of the strongest components of the jet magnetic field, in PLUTO simulation units.We can see that the toroidal component (ringwise) constitutes a significant part of the field, while other magnetic lines form part of the poloidal component, especially near the jet axis.

Figure 2 :
Figure 2: A plot of the magnetic field magnitude roughly half way into the simulation.We can see the jet self-confinement due to magnetic forces resulting in a narrow beam.

Figure 3 :
Figure3: A plot of the jet density roughly halfway into the simulation.We can see the jet beam advancing along a narrow path through the ambient medium, kept together by the magnetic toroidal field component.This confined jet is then a suitable site for increased high energy -ray and neutrino emission/production, due to maintaining extreme conditions along the jet length.

Figure 4 :
Figure 4: Unnormalized neutrino intensity,  ] ( ] ), on Earth as a function of the neutrino energy,  ] , per energy interval.This neutrino intensity is coming out of the jet simulations in the presence of the strong toroidal magnetic field that appreciably affects the emission.

Table 1 :
Values of various physical and model parameters pertaining to the -ray emission from the jet, for the simulation run, obtained with the linear method and the Hancock integrator.