Improved photomeson model for interactions of cosmic ray nuclei

TL;DR

This paper introduces an improved, data-driven model for photonuclear interactions of cosmic ray nuclei above the pion production threshold, enhancing the accuracy of neutrino flux predictions in astrophysical scenarios.

Contribution

The authors develop a more accurate, inclusive model for nuclei interactions in the photomeson regime, improving upon the superposition approximation used previously.

Findings

Moderate changes in UHECR mass composition.

Significant reduction in astrophysical neutrino flux predictions.

Model implementation is publicly available.

Abstract

Photon-hadronic interactions are important for the sources and the transport of Ultra-High Energy Cosmic Rays (UHECRs). Current state-of-the-art cosmic ray transport simulations handle nuclear disintegration at energies of the Giant Dipole Resonance at a more sophisticated level, as well as the photohadronic interactions of nucleons in the high-energy regime above the pion production threshold. However, the interactions of nuclei above the pion production threshold are commonly modeled by treating the nucleus as a superposition of free nucleons -- ignoring the effect of the nuclear medium. We construct an improved, inclusive model for the photomeson regime for nuclei with $A \leq 56$ by employing more accurate, data-driven parametrizations of the interaction cross section, the fragmentation of the primary nucleus and the inclusive pion production cross section that directly affects the production of astrophysical neutrinos. We apply our results to two multi-messenger scenarios (Tidal Disruption Events and Gamma-Ray Bursts) in which photonuclear interactions in the photomeson regime are the dominant cooling process for the highest energy cosmic rays. While we find moderate changes to the mass composition of UHECRs, the astrophysical neutrino fluxes exhibit a significant (factor of a few) reduction compared to the na\"ive superposition of free nucleons for sources of UHECR nuclei with a populated cascade. The numerical code implementing the model has been made publicly available, which facilitates the integration of our results in similar frameworks.

Tables

Table 1. Table 1 : Comparison of photomeson models and schematics of the related physical picture. The nuclei and nuclear fragments are represented as collections of circles (nucleons) with color related their role in the interaction (and fraction of the photon energy they receive): blue, spectator nucleons (no extra energy received, boost conservation); purple, non-active participants (some extra energy received); and orange, active participants (receive most of the energy, direct interaction). Green smaller circles represent pions.

		Single Particle Model (SPM)	Residual Decay Model (RDM)	Empirical Model (EM)
Schematics of physical process
${\sigma }_{A\gamma }$		$A{\sigma }_{p\gamma }^{\text{Sophia}}$	$A^{\alpha }{\sigma }_{\mathrm{univ}}$, $\alpha =1$
	${ϵ}_r>1\mathrm{GeV}$		$A^{\alpha }{\sigma }_{\mathrm{univ}}$, $\alpha =\alpha \left({ϵ}_r\right)$
${\sigma }_{A\gamma \to X}^{\mathrm{incl}}$	Fragments produced	n/p and $A-1$	n/p and disintegration products at $E^{\ast }$	$A_{\mathrm{fr}}\in [1-4,$ $A/2\dots A-1]$ (see Appendix 3.3)
	Fragment‘s multiplicity	Always 1 (for any fragment produced)	Combination of disintegration model and Sophia	Empirical formulas and nuclear thermostatistics (see Appendix 3.3)
	${\sigma }_{A\gamma \to \pi }^{\mathrm{incl}}$	$A{\sigma }_{p\gamma \to \pi }^{\mathrm{incl},\text{Sophia}}$		$A^{{\alpha }_{\pi }}{\sigma }_{\pi }^{\mathrm{incl}}$ (see Appendix 3.2)

Table 2. Table 2 : The source parameters shown were taken from the references [ 22 , 29 ] . In the case of the GRB, the parameters for the Optically Thick scenario are shown; the ones for the other scenarios (Empty Cascade and Populated Cascade) are identical except for their luminosities L γ = 10 49 ​ erg / s subscript 𝐿 𝛾 superscript 10 49 erg s L_{\gamma}=10^{49}\rm erg/s and L γ = 10 51 ​ erg / s subscript 𝐿 𝛾 superscript 10 51 erg s L_{\gamma}=10^{51}\rm erg/s , respectively. Parameters critical for the photomeson processes dominating at the highest energies are highlighted in boldface.

Source	TDE [29]	GRB [22]
Gamma factor	$\Gamma =10$	$\Gamma =300$
Redshift	$z=0.001$	$z=2$
Duration	$T={10}^{6}s$	$T=10s$
Luminosity	$L_X={10}^{47}\mathrm{erg}/s$	$L_{\gamma }={10}^{53}\mathrm{erg}/s$
Collision Radius	$R={10}^{9.6}\mathrm{km}$	$R={10}^{8.3}\mathrm{km}$
Injected isotope	${}_{}{}^{14}N$	${}_{}{}^{56}\mathrm{Fe}$
Acceleration efficiency	$\eta =1$	$\eta =1$
Baryonic loading	$\xi =10$	$\xi =10$
Target photon spectrum parameters in the shock rest frame	${\epsilon }_{X,\mathrm{br}}^{\prime }=1\mathrm{keV}$	${\epsilon }_{\gamma ,\mathrm{br}}^{\prime }=1\mathrm{keV}$
	${\epsilon }_{X,\min }^{\prime }={10}^{-6}\mathrm{eV}$	${𝜺}_{𝜸,\mathrm{𝐦𝐢𝐧}}^{\prime }=\mathrm{𝟏𝟎𝟎}\mathrm{𝐞𝐕}$
	${\epsilon }_{X,\max }^{\prime }=300\mathrm{keV}$	${\epsilon }_{\gamma ,\max }^{\prime }=300\mathrm{keV}$
	$𝜶=-\mathrm{𝟐}/\mathrm{𝟑}$, $\beta =-2$	$\alpha =-1$, $\beta =-2$

Equations

$$\Gamma (E)=\int \text{d}\epsilon {\int }_{-1}^{+1}$$

$$Q_{ji}(E_i)=\int d{E}_j{N}_j(E_j){\Gamma }_j(E_j)\frac{\mathrm{d}{n}_{j\to i}}{\mathrm{d}{E}_i}(E_j,E_i),$$

$$\frac{\mathrm{d}{n}_{j\to i}}{\mathrm{d}{E}_i}(E_j,E_i)=\frac{1}{{\sigma }_j(E_j)}\frac{d{\sigma }_{j\to i}^{\mathrm{incl}}}{\mathrm{d}{E}_i}(E_j,E_i)$$

$$M_{j\to i}(E_j)=\frac{1}{{\sigma }_j(E_j)}\int \mathrm{d}{E}_i\frac{d{\sigma }_{j\to i}^{\mathrm{incl}}}{\mathrm{d}{E}_i}(E_j,E_i).$$

$$\frac{\mathrm{d}{\sigma }_{j\to i}^{\mathrm{incl}}}{\mathrm{d}{E}_i}(E_j,E_i)\approx {\sigma }_j(E_j)M_{j\to i}(E_j)\delta \left(E_i-\frac{A_i}{A_j}{E}_j\right).$$

$${\sigma }_{A\gamma }({ϵ}_r)=A_{\mathrm{eff}}({ϵ}_r){\sigma }_{p\gamma }({ϵ}_r)=A^{\alpha ({ϵ}_r)}{\sigma }_{p\gamma }({ϵ}_r).$$

$${\sigma }_{{\pi }^{\pm }}^{\mathrm{incl}}({ϵ}_r)=\frac{M_{N\gamma \to {\pi }^{\pm }}^{\text{textsc}\mathrm{Sophia}}({ϵ}_r)}{M_{N\gamma \to {\pi }^0}^{\text{textsc}\mathrm{Sophia}}({ϵ}_r)}{\sigma }_{{\pi }^0}^{\mathrm{incl}}({ϵ}_r).$$

$${\sigma }_{A\gamma \to \pi }^{\mathrm{incl}}=A^{{\alpha }_{\pi }({ϵ}_r)}{\sigma }_{\pi }^{\mathrm{incl}}=\{\begin{array}{ll}{A}^{{\alpha }_{\pi }({ϵ}_r)}{\sigma }_{N\gamma \to \pi }^{\text{textsc}\mathrm{Sophia}}& \mathrm{for}1\mathrm{GeV}<{ϵ}_r<50\mathrm{GeV}\\ A{\sigma }_{N\gamma \to \pi }^{\text{textsc}\mathrm{Sophia}}& \mathrm{for}{ϵ}_r>50\mathrm{GeV}\end{array},$$

$$R_M^{\pi }({ϵ}_r)=\left(\frac{{\sigma }_{A_X\to \pi }^{\mathrm{incl}}}{A_X}\right)/\left(\frac{{\sigma }_{A_D\to \pi }^{\mathrm{incl}}}{A_D}\right).$$

$$R_M^{\pi }({ϵ}_r)={\left(\frac{A_X}{A_D}\right)}^{{\alpha }_{\pi }({ϵ}_r)-1}.$$

$${\sigma }_{j\to k}^{\mathrm{incl}}({ϵ}_r)={\sigma }_j({ϵ}_r)\left({\delta }_{\mathrm{kN}}+\frac{{\sigma }_{(j-1)\to k}^{\text{textsc}\mathrm{talys}}({ϵ}^{\ast })}{{\sigma }_{(j-1)}^{\text{textsc}\mathrm{talys}}({ϵ}^{\ast })}\right),$$

$$M_{j\to k}=\frac{{\sigma }_{j\to k}^{\mathrm{incl}}}{{\sigma }_j}$$

$${\sigma }_{j\to \mathrm{p}}^{\mathrm{incl}}$$

$${\sigma }_{j\to \mathrm{n}}^{\mathrm{incl}}$$

$${\sigma }_{j\to k}^{\mathrm{incl}}$$

$${\sigma }_{\mathrm{p}}^{\mathrm{dir}}=0.078A^{0.50}\mathrm{mb}.$$

$${\sigma }_{\mathrm{n}}^{\mathrm{dir}}=0.104A^{0.81}\mathrm{mb}.$$

$${\sigma }_{\mathrm{ntot}}^{\mathrm{mul}}$$

$${\sigma }_{x\mathrm{n}}^{\mathrm{mul}}$$

$${\sigma }_{x\mathrm{p}y\mathrm{n}}^{\mathrm{sp}}$$

$$\mathsf{E}$$

$${\sigma }_{\pi }^{\mathrm{prod}}=0.027A^{0.847}\mathrm{mb}.$$

$${\sigma }^{\mathrm{sp}}={\sigma }^{\mathrm{tot}}-\left({\sigma }_{\pi }^{\mathrm{prod}}+{\sigma }_{\mathrm{p}}^{\mathrm{dir}}+{\sigma }_{\mathrm{n}}^{\mathrm{dir}}+{\sigma }_{\mathrm{ntot}}^{\mathrm{mul}}\right).$$

$${\sigma }_{\mathrm{n}}^{\mathrm{mul}}=\sum _{x=2}^{x_{\max }}x{\sigma }_{x\mathrm{n}}^{\mathrm{mul}}.$$

$${\sigma }_k^{\mathrm{sp}}=\{\begin{array}{ll}{\sigma }_l^{\mathrm{sp}}& A_k\in [\mathrm{1...4}]\\ {\sigma }_{x\mathrm{p}y\mathrm{n}}^{\mathrm{sp}}& A_k\in [A/2...A-2]\end{array}.$$

$${\mathrm{C}}_i^r=\left\{{\mathrm{N}}_{i,\mathrm{n}}^r,{\mathrm{N}}_{i,\mathrm{p}}^r,...,{\mathrm{N}}_{i,l}^r,...,{\mathrm{N}}_{i,\alpha }^r\right\}l\in (Z_l,A_l),1\le A_l\le 4,$$

$$x=Z_i=\sum _{{\mathrm{C}}_i^r}{Z}_l{\mathrm{N}}_{i,l}^r,$$

$$x+y=A_i=\sum _{{\mathrm{C}}_i^r}{A}_l{\mathrm{N}}_{i,l}^r.$$

$${\mathrm{N}}_{i,l}=\sum _r{\mathrm{N}}_{i,l}^r{\mathrm{P}}_{i,r}.$$

$${\sigma }_l^{\mathrm{sp}}=\sum _i{\mathrm{N}}_{i,l}{\sigma }_i^{\mathrm{sp}},l\in \mathbb{S}.$$

$${\mathrm{N}}_{i,l}={\mathrm{P}}_i\sum _r{\mathrm{N}}_{i,l}^r,$$

$$A\sum _i{\sigma }_i^{\mathrm{sp}}=\sum _k{A}_k{\sigma }_k^{\mathrm{sp}},$$

$$\mathsf{Z}=\sum _{i,r}{e}^{-\beta ({ϵ}_i-\mu n_r)},$$

Full text

1 Introduction

The mass composition measurements of Ultra-High Energy Cosmic Rays ($E>10^{18}$ eV) from the Pierre Auger Observatory and the Telescope Array are clearly compatible to a mixed composition [1], supporting the presence of nuclei heavier than protons in astrophysical accelerators, such as Active Galactic Nuclei [2], Gamma Ray Bursts [3] or Tidal Disruption Events [4, 5]. Interactions of UHECR occur near the acceleration site for high in-source photon densities and during the transport between source and Earth with the cosmic microwave background (CMB), the infrared and the optical photons. We refer to these interactions as photohadronic interactions of nuclei.

Photohadronic interactions occur at different energy scales [6], where the energy scale typically refers to the photon energy in the nucleus’ rest frame $\epsilon_{r}$. At low energies $\epsilon_{r}\gtrsim 1$ MeV hadronic processes can not occur, hence particles are produced via electromagnetic (Bethe-Heitler) electron-positron pair production. At energies comparable to the nuclear binding energy $\epsilon_{r}\gtrsim 8$ MeV, a resonant motion of the nucleons leads to the disintegration of the nucleus through the excitation of the Giant Dipole Resonance (GDR). Models that parametrize this process for astrophysical applications are based on tables by (for example) Puget-Stecker-Bredekamp (PSB) [7, 8], or are computed from running nuclear event generators, such as Talys [9], GEANT 4 [10] or Fluka [11]. Beyond $\epsilon_{r}\gtrsim 140$ MeV (hadronic scale) hadronic processes start, leading to the production of baryonic resonances that emit pions and photons when they decay. The lowest mass resonances, such as the $\Delta$-resonance, are produced in the $s$-channel for photon projectiles and are thus more efficient in dissipating the available energy into secondary particles compared to proton-proton collisions where they are produced from t-channel processes at threshold. At higher energies $\epsilon_{r}\gtrsim 1$ GeV, the hadronic mass scale, the photon interacts hadron-like as a virtual vector meson. We refer to the latter two energy scales ($\epsilon_{r}\gtrsim 140$ MeV) as “photomeson” regime.

The physics of photomeson production in astrophysical environments has been studied extensively for nucleons (protons and neutrons) as early as [12, 13] and later in Refs. [14, 15, 16, 17, 6, 18]. The Sophia [19] Monte Carlo code is tailored to simulate photon-nucleon interactions around the threshold, including the most important baryonic resonances and direct production ($t$-channel) processes. The model at higher energies is somewhat simpler, but sufficiently precise for modeling multi-pion production in the scope of astrophysical applications. Different parametrizations using Sophia for semi-analytical and numerical transport equation solvers have been developed in Refs. [20, 21, 22]. Nuclei, however, are treated in the superposition approach, in which photons interact with one nucleon and the remainder is left intact. For this single photon-nucleon interaction, the differential cross sections for secondary particles are sampled from Sophia, whereas the remnant nucleus simply has the mass $A-1$ neglecting momentum transfer or recoil, see e.g. Refs. [23, 24, 25, 26, 27]. The effective number of target nucleons seen by the photon projectile is reflected in the “mass scaling” of the cross section $A_{\rm eff}=\sigma_{A\gamma}/\sigma_{p\gamma}\sim A^{\alpha}$. The power alpha lies between two extreme values: $1$ when all nucleons are seen by the photon as targets, and $2/3$ when only nucleons at the surface are active targets and those inside are screened. Refs. [28, 25] make a global scaling assumption for the entire photomeson energy range. In reality, different nuclear scaling powers $\alpha$ are appropriate for the total interaction and the pion production cross sections, with different values for energies close to the threshold and for very high energies.

If photomeson production is the dominant cooling process in a UHECR source, the disintegration/ fragmentation of the UHECR will be described by the interactions above the $\Delta$-resonance threshold. For intergalactic UHECR propagation, the giant dipole resonance (GDR) dominates the photohadronic processes because the energy and high density of the CMB photons constitutes the dominant cooling process. This happens at UHECR energies below the photomeson threshold, hence lowering the importance of the photomeson model for propagation calculations. Under certain circumstances the photomeson regime can dominate the photohadronic interactions in astrophysical sources, if the UHECR does not find interaction partners (photons) to match the GDR energy. Known examples of this are Gamma-Ray Bursts, with a minimal photon energy cutoff [22] which may come from synchrotron self-absorption [24], and jetted Tidal Disruption Events involving massive stars [29] with a spectral photon index such that photon densities of lower energy photons are smaller than usual, making photomeson production the dominant interaction for UHECRs of at the highest energies. We will use these as astrophysical examples in our discussion of the impact of the photomeson model.

Apart from the primary UHECRs, the neutrino production is always related to the photomeson regime: if the neutrinos are mostly produced off nuclei, corrections to the photomeson model are expected to have an impact, whereas if the neutrinos are produced off nucleons, the physics is well described by, for instance, Sophia. The former situation mostly occurs in sources with low radiation fields, which means that the primary nucleus barely disintegrates, whereas the latter situation is typical for sources with strong radiation fields in which nucleons emerging from a nuclear cascade dominate the neutrino production [22].

In this paper we will introduce an improved model for the photomeson production and compare it to the widespread “Single Particle Model” (SPM). The improvements involve modifications relevant for the astrophysical processes: improved cross section descriptions driven by data, and empirical approaches to describe the fragmentation of the nucleus. The impact of these modifications is illustrated at astrophysical simulations of Gamma-Ray Bursts (GRBs) and Tidal Disruption Events (TDEs). The paper is organized as follows: Section 2 discusses the physical differences between photon interactions with nucleons and with nuclei, and introduces the quantities and nomenclature related to multi-messenger astroparticle physics simulations. Section 3 presents the ingredients of the new photomeson model, and Section 4 illustrates the impact of the model in astrophysical source simulations in which the photomeson model is relevant.

2 Astrophysical photohadronic interactions: $\boldsymbol{p\gamma}$

vs $\boldsymbol{A\gamma}$

A convenient reference frame to discuss photonuclear interactions is the nuclear rest frame, in which the photon energy ${\epsilon_{r}=E\,\varepsilon\,(1-\cos\theta)/m_{j}}$ depends on the energy of the relativistic nucleus $E$ and the photon energy $\varepsilon$ in the observer’s or (cosmological) comoving frame. The pitch angle $\theta$ is the angle between incident photon and nucleus such that $\cos\theta=-1$ represents head-on collisions, and $\epsilon_{r}$ is related to the center-of-mass energy by $s=m_{j}^{2}+2m_{j}\epsilon_{r}$ where $m_{j}$ is the mass of the nucleus. The photonuclear interaction rate and the interaction cross section $\sigma$ are related as

[TABLE]

where $n_{\gamma}$ is the photon number density and the rate is expressed in units of inverse length. Depending on the type of source or environment, the photon spectrum can extend from sub-eV up to TeV energies, and its shape can contain peaked (thermal) or power-law (non-thermal) components.

The shape of a typical photonuclear cross section is illustrated in Fig. 1 for $\rm{}^{56}Fe$. The energy range is split in a photodisintegration regime, characterized by the absence of hadron production and negligible momentum transfers or recoils compared to typical cosmic ray energies (boost conservation); and a photomeson regime, with the possibility for production of all sorts of particles (predominantly pions) but with significant momentum transfer. A more detailed discussion about the role of photodisintegration for cosmic ray astrophysics is contained in e.g. [26, 30, 31].

It is instructive to first discuss the photomeson production in the absence of collective effects, i.e. in interactions of photons with free nucleons, which has been extensively studied experimentally and theoretically [35, 36]. The proton and neutron cross sections are different below $\approx 140\,\rm MeV$, since for protons Thomson scattering and pair-production are possible, whereas only the much weaker magnetic moment scattering occurs for the neutron [37]. The production of pions at threshold occurs through the excitation of the lightest baryonic resonance ($\Delta$) in resonant ($s$-channel) scattering. This process has no counterpart in $pp$ scattering since there are no known di-baryon resonant states, and hence no $s$-channel equivalent. Instead, mesons are produced through $t$-channel processes at sufficiently large momentum transfer. Above the pion threshold more channels are available for the production of higher mass resonances and there are small differences between proton and neutron cross sections (see Fig. 2). At high energies above a few GeV, the photon interacts mostly as a virtual vector meson (see for instance LABEL:~[38]) and all phenomena of hadronic interactions can occur.

The transport equations contain re-injection terms that represent the production rate of particles of species $i$ at energy $E_{i}$ converting from the particle species $j$ at energy $E_{j}>E_{i}$ through interactions or decays. The general form for the re-injection rate is

[TABLE]

where $N_{j}$ is the density of $j$ particles, $\Gamma_{j}$ the interaction rate from Eq. (2.1) or a decay rate and

[TABLE]

an energy redistribution function. By integrating the inclusive differential cross section $d\sigma^{\rm incl}_{j\to i}/{{\rm d}E_{i}}\,(E_{j},E_{i})$ we define the multiplicity

[TABLE]

that has the meaning of the average number of particles of species $i$ produced per interaction.

If $i$ and $j$ are both nuclei, the energy redistribution takes a simple form when assuming the conservation of boost $E_{i}/m_{i}=E_{j}/m_{j}$ (or the energy per nucleon)

[TABLE]

For species $i$ other than the remnant nucleus (such as $\pi^{\pm}$, $\pi^{0}$, secondary $\rm p$, $\rm n$ and higher mass hadrons), the redistribution function does not have a simple parametrization and has to be obtained from Monte Carlo simulations (e.g. Sophia).

The simplest extension of the free nucleon interactions to nuclear interactions is referred to as Single Particle Model (SPM) in the cosmic ray astrophysics literature [28, 21]. It assumes that in the photomeson regime the photon interacts always with one nucleon in the nucleus without affecting the rest of the nucleus (quasi-free interaction). The final state particles are the products of the $\gamma N$ interaction and one remnant nucleus with $A-1$ nucleons. The inelastic cross section is frequently assumed to scale with $A$, i.e. $\sigma_{A\gamma}=A\sigma_{p\gamma}$, implying that ${\rm d}\sigma^{\rm incl}_{j\to i}/{\rm d}E_{i}=A_{j}\,{\rm d}\sigma^{\rm incl}_{N\to i}/{\rm d}E_{i}$.

The simplicity of the SPM allows estimating analytically the relative importance of photo-nucleon to photonuclear interactions in astrophysical scenarios. However the mass scaling of the total cross section is reduced with the increase of photon energy $\epsilon_{r}\gtrsim 1$ GeV. Additionally, nuclear medium effects cause differences in the shape of the cross sections. The $\Delta$-resonance in nuclei is broadened and the higher resonances are smeared out as a result of the Fermi motion of nucleons and the in-medium properties of the nucleon resonances and mesons. The production of pions is also strongly affected because of final state interactions (FSI) of the photo-nucleon products. Finally, the disintegration of the nucleus is more important than that considered in the SPM, and multiple fragments can be produced with masses different from $A-1$. The photonuclear interaction can lead to photospallation and fission processes which disintegrate the nucleus into nucleons, deuterons, $\alpha$-particles and larger fragments. The fission processes are more important for higher masses than the ones relevant in this paper, but the model includes it to allow studies where higher masses might be needed. In the next section we implement these mechanisms in a photomeson model with the simplicity characteristic of the SPM.

3 A new model: Empirical photomeson Model (EM)

We improve the SPM in three main aspects:

The absorption cross section: A “universal function” better describes the shape of the cross section per nucleon the near the $\Delta$ resonance energy. At higher energies an energy dependent mass scaling exponent is introduced to account for nuclear shadowing effects. Both modifications are motivated and derived from data. 2. 2.

The pion production cross section: Pion production is known to be strongly influenced by nuclear medium effects [39, 40]. The inclusive pion production cross section is derived from data [41] for different nuclei and parameterized with an additional curve and a mass scaling exponent that differs the absorption cross section. 3. 3.

Nuclear fragmentation: The photonuclear interactions can result in nuclear breakup through different mechanisms not included in the SPM. We implement two alternatives: an Ablation-Abrasion inspired model and an empirical data-driven parametrization. The latter is found in better agreement with detailed simulations and is chosen as our baseline model.

In the following sub-sections, we assess each of the aspects individually. The model described here has been made available in a code [42] and it can be used for reproducing this work, as well as for implementing in other frameworks studying UHECR photomeson interactions.

3.1 Total photonuclear cross section

As discussed in Section 2, the nuclear environment (also referred to as “medium effects”) changes the physics above the pion threshold energies compared to the free nucleon case. Fig. 3 shows the inelastic cross section divided by $A$ as a function of energy. The curves corresponding to protons (SPM curve) and to $\rm{}^{56}Fe$ (EM curve) are compared to photonuclear data compiled from experiments with various target nuclei [43, 44, 33]. At $\epsilon_{r}\lesssim\rm 1\,GeV$, the green circles correspond to light nuclei ($A=2-4$) and the blue circles to $A=7-208$. Within the errors of the data, a scaling with $A$ appears justified. However, the shape of the curve is different compared to $\sigma_{p\gamma}$, with only one pronounced resonance peak in place of the $P_{33}(1232)$ resonance ($\Delta$-resonance), being 20% wider at half height and 30% lower at the peak. The widening has been explained with the $\rm Delta$-hole model [45], where medium effects are taken into account by including the Fermi motion of nucleons, the Pauli blocking restricting decay channels and the $\Delta$-N interactions. The nucleon resonances at energies beyond $500\,{\rm MeV}$ are not visible even for small masses such as Be and C [46, 47].

For $\epsilon_{r}\gtrsim\rm 1\,GeV$ the photonuclear cross section per nucleon is reduced due to the shadowing effect. This reduction is dependent on the nuclear mass, and has been successfully understood within the Vector Meson Dominance model (VMD), which describes the photon’s wave function as a superposition of mesonic states ($\rho$, $\omega$, $\varphi$) [48, 59, 49] that interact hadronically with the nucleus. At higher energies, where the photon can resolve the partons in the nucleus, the parton distribution function is predicted to be high enough that the nucleus becomes opaque to photons, leading to a theoretical limit similar to that of hadron-nucleus interactions $\alpha_{\rm lim}=2/3$ [50]. The mass scaling is typically parametrized in the literature as a power of the nuclear mass:

[TABLE]

In the EM the cross section is parametrized parametrization with the following additional elements:

• •

For $\epsilon_{r}<1$ GeV a “universal curve” (spline fit of data points in Fig. 1) scaled by A ($\alpha=1$) is used. This reflects the universal shape of the cross section per nucleon exhibited by nuclei of a wide range of masses.

• •

For $\epsilon_{r}\geq 1$ GeV the photo-nucleon cross section is calculated with Sophia, scaled by an energy dependent exponent $A_{\rm eff}=A^{\alpha(\epsilon_{r})}$. The energy dependence of the exponent is shown in the insert of Fig. 1.

The data provided in LABEL:~[33] are provided as $f_{A}(\epsilon_{r})=A_{\rm eff}(\epsilon_{r})/A$ for different nuclei (C, Cu, and Pb). Considering the parametrization with mass in Eq. (3.1), the exponent values $\alpha$ are calculated from the data using $\alpha(\epsilon_{r})=1+\ln{f_{A}(\epsilon_{r})}/\ln{A}$ (purple points in the insert of Fig. 3). A linear fit was performed to find the energy dependence, and extrapolated towards the shadowing limit $\alpha_{\rm lim}=2/3$ for high energies (yellow curve in insert). The shaded area contains the 95% confidence interval of the fit in the region where data are available, and has also been extrapolated to higher energies: the upper band with a constant value, the lower band with a transition towards $\alpha_{\rm lim}=2/3$ (Fig. 3).

The mass scaling dependence with energy introduced in th EM at higher energies results in an $A$- dependence of the cross section per nucleon $\sigma_{A\gamma}/A$. The hashed region in Fig. 3 represents this mass dependence for the range of masses $A=4-56$, where the upper line corresponds to the smaller mass and the lower to the larger mass.

3.2 Photoproduction of pions off nuclei

The photoproduction of pions off nuclei is very sensitive to nuclear medium effects near threshold energies [39]. The interaction of the photon can occur via a “quasi-free” process, where one nucleon and the pions produced are ejected while the rest of the nucleus is unaffected. However, the pion(s) produced inside the nucleus interact with the surrounding nucleons with a strong dependence on their kinetic energy. Pions interact weakly with nucleons inside the nucleus for kinetic energies below $\sim 40\,{\rm MeV}$, and exhibit a resonance around kinetic energies $50-100\,{\rm MeV}$ which depends on the nuclear mass [51]. This effect is also evidenced in the photoproduction of pion pairs [52], where the mass scaling of the cross section decreases from $A$ to $A^{2/3}$ as the photon energy becomes large enough to produce pions with kinetic energies above $\sim 40\,{\rm MeV}$. An effective description for the mass scaling [53] can be achieved by separating “surface” and “volume” contributions, but such a description is hard to generalize for a broad range of nuclear masses.

Pion photoproduction off nuclei can be calculated with existing codes using Monte Carlo techniques in a cascade scenario (CRISP [54]), or in terms of transport equations (GiBUU [55]). These simulations can be set up to have a good overall agreement with data. However, a general table for a wide range of nuclear species is computationally expensive and at the same time, the details of pion emission from each species are usually less important in astrophysical simulations compared to the uncertainties of nuclear composition. To keep the EM description of the pion photoproduction off nuclei as simple as in the SPM, the same shape for the inclusive cross section is used for all nuclei with a mass dependence and scaling coefficient $\alpha_{\pi}(\epsilon_{r})$.

Fig. 4 (left) shows the inclusive $\pi^{0}$ photoproduction as a function of incident photon energy $\epsilon_{r}$ for energies below $\sim 1\,{\rm GeV}$. The data were measured for different nuclei and normalized to $A^{2/3}$ (values reproduced from Fig. 9 in LABEL:~[41]). The similarities in the curves point to a dominating “surface”-like mass scaling ($A^{2/3}$) rather than a “volume”-like mass scaling ($A$). The systematic differences between the different curves are due to second order processes which scale with volume like multi-pions produced at low kinetic energies [53]. The data have been fit with a spline (black solid curve), which represents the inclusive pion production cross section fit to all nuclear masses assuming a scaling with $A^{2/3}$. Fig. 4 (right) compares the $\rm{}^{40}Ca$ pion production cross section with different parametrizations. The SPM curve with $A$ scaling (dot-dashed red line) predicts on average more than twice the number of pions compared to data (blue squares), whereas when using $A^{2/3}$ scaling (dashed blue) it leads to less than half. The total cross section per nucleon (universal curve $\sigma_{\rm univ}$) scaled by the pion multiplicity per nucleon $M^{\rm\textsc{Sophia}}_{p\gamma\to\pi^{0}}$ and a $A^{2/3}$ mass scaling (dotted orange) underestimates data, as well. The spline fit to all available data in Fig. 4 (left) scaled by $A^{2/3}$ describes the pion production sufficiently well.

At higher energies heavier mesons start playing a role in the final states, affecting the simplified scaling assumption. For example, studies of photoproduction of $\eta$ mesons up to $2\,{\rm GeV}$ [56] find consistency with a $A^{2/3}$ scaling when the $\eta$ is produced with kinetic energies sufficiently above threshold, demonstrating strong absorption through FSI. These type of effects are currently studied but the theoretical description is not yet complete [57, 58]. However, as the photon energy increases such effects are less important, and have been shown to become irrelevant for photon energies beyond $50\;{\rm GeV}$ [59] where the production of pions is consistent with $A$ scaling.

Our inclusive description with a mass scaling coefficient does not capture the complexity of the exclusive final states. Hence, our model follows the available data that suggests a $A^{2/3}$ scaling near threshold with a common curve for all species of nuclei, and is extended to higher energies by allowing the exponent to increase towards the limit value 1.

The pion photoproduction off nuclei in the EM is, therefore, parametrized as $A^{\alpha_{\pi(\epsilon_{r})}}\sigma^{\rm incl}_{\pi}$ where $\sigma^{\rm incl}_{\pi}(\epsilon_{r})$ is an inclusive cross section for $\pi^{+}$, $\pi^{-}$ and $\pi^{0}$ production, as described below. Below $1\;{\rm GeV}$ $\sigma^{\rm incl}_{\pi}(\epsilon_{r})$ is derived from data: the form of $\sigma^{\rm incl}_{\pi^{0}}$ is obtained by fitting a spline to the experimental values (black curve in Fig. 4). For charged pions, their multiplicities are adjusted proportional to the production ratio off nucleons as taken from Sophia:

[TABLE]

Above $1\;{\rm GeV}$ the EM curve is extended as in the SPM but with a smooth transition from $A^{2/3}$ scaling to $A$ scaling with photon energy (see Fig. 5) at around $50\;{\rm GeV}$:

[TABLE]

satisfying also the relations in Eq. (3.2).

As demonstrated in Fig. 5 the energy dependence of the pion scaling exponent $\alpha_{\pi}(\epsilon_{r})$ is derived from experimental data [60, 61, 62]. The data has been obtained in deep inelastic scattering experiments measuring the semi-inclusive hadron production from the interaction of virtual photons with nuclei. The reported magnitude is the multiplicity ratio, which describes ratio of pion production per nucleon in a nucleus $X$ to that in Deuterium $D$ (see Eq. (3.4)).

[TABLE]

It follows from Eq. (3.3) that the EM the multiplicity ratio is related to the photon energy nucleus mass and the scaling exponent as:

[TABLE]

Below $1\;{\rm GeV}$ the scaling exponent is fixed to $2/3$ as suggested by data in Fig. 4. The multiplicity ratio data, from different nuclear species, show a sharp increase of the scaling exponent in the interval from $1-3\,{\rm GeV}$. A spline fit to the entire energy range (dash-dotted gray curve in Fig. 5) underestimates data from $4-10\,{\rm GeV}$. For this reason the model joins the the two energy ranges below 3 and $4-22\,{\rm GeV}$ using a sigmoid function (dashed black line). This curve applies to all pion species, since they are produced in equal amounts [61].

The energy redistribution of pions (differential cross sections) (see Eq. (2.3)) used are those sampled from Sophia for free nucleons, i.e. any modification to the angular distributions resulting from the nuclear medium is neglected.

3.3 Nuclear fragmentation

In the SPM a photonuclear interaction of a nucleus with mass $A$ always results in the loss of one nucleon per interaction. Hence, the final state always includes a remnant nucleus of mass $A-1$ that does not further disintegrate. In the photomeson regime nuclear breakup is governed by several mechanisms that lead to a more complex final state, including the emission of multiple nucleons, light fragments and a distribution of remnant nucleus masses. We consider two physical approaches to model this effect:

Residual Decay Model (RDM): The initial photon-nucleon interaction is modeled in the same way as in the SPM, assuming a quasi-free target nucleon. The motion of the nucleon through the nuclear medium results in an average energy loss, proportional to the nuclear radius (Abrasion-Ablation hypothesis). This energy left in the residual nucleus with mass $A-1$ is an excitation energy which can been estimated as $\epsilon^{*}=17\,{\rm MeV}(A-1)^{1/3}$ [6]. The subsequent de-excitation proceeds through a re-scattering, or an equivalent interaction a photon with the energy $\epsilon_{r}=\epsilon^{*}$ with the remnant nucleus ($A-1$). The typical energies for the range of masses are close to the GDR regime, for which the multiplicity distributions for final state particles have been obtained from a statistical, Hauser-Feshbach based model Talys [9] that uses combinations different nuclear models to describe the disintegration. These tabulated cross sections have been previously discussed in LABEL:~[26]. The inclusive cross section in the RDM for producing the fragment $k$ from the interacting species $j$ is calculated by normalizing the multiplicity of the fragments produced in the disintegration model at $\epsilon^{*}$ to the total cross section of the interacting species $\sigma_{j}$

[TABLE]

where $\delta_{\rm kN}$ (1 if k is proton or neutron and zero otherwise) accounts for the additional nucleon emitted in the first photon interaction.

Empirical Model (EM): The photonuclear interaction is more complex, and often involves multiple nucleons and intermediate fragments. The production of nuclear fragments is very dependent on the absorption mechanism and the consequent energy transfer to the rest of the nucleus (intra-nuclear cascade). The average production of masses can be modeled using empirical formulas previously obtained in LABEL:~[63] and reproduced in Appendix A.1. The formulas estimate the average cross sections for producing certain fragments in the energy range $0.2\lesssim\epsilon_{r}\lesssim 1\,{\rm GeV}$. Employing those relations and additional considerations from thermostatistics, the inclusive cross sections for a broad range of nuclear fragments have been produced (see Appendix A.1). Schemes of the physical scenarios and a summary of the components of these models are shown in the Tab. 2 in comparison with the SPM.

In Fig. 6 a comparison between of the mass distribution of fragments is shown for two interacting nuclei: $\rm{}^{14}N$ (left) and $\rm{}^{56}Fe$ (right). The multiplicities are calculated from the total cross section of the interacting species $j$ and the inclusive cross section for production of species $k$:

[TABLE]

which is equivalent to Eq. (2.4) and represents the mean number of particles of a certain species produced per interaction. The species are grouped by mass and the multiplicity of the mass group is the sum of the multiplicities of all species with the same mass. The EM values are shown in solid blue, and the RDM in solid green. Additionally, the values obtained with Fluka [11] are shown as reference. In Fluka, detailed Monte Carlo calculations are combined with state of the art theory and nuclear data to simulate the photonuclear interaction and the subsequent production of secondaries. The shape of the EM is in better overall agreement with Fluka than the RDM. The implicit assumptions in the RDM are reflected in the narrower mass distribution which is a typical shape produced in the disintegration at lower energies around the GDR, where the disintegration model for the residual nucleus is sampled. Experimental results for the disintegration at energies corresponding to the photomeson regime point to processes like spallation and fission [65, 66] (see Fig. 11) which are explicitly simulated in Fluka whereas in EM the mass distribution is fixed for all species of the same mass. Fig. 7 shows a comparison of EM and Fluka to data from [64] where the average production over photon energies $\rm 0.3-1\,GeV$ is measured for multiple nuclei. The corresponding values from the EM and Fluka are within a factor 1-3 of the experimental values. For this case, the EM performs similarly to a more elaborate nuclear code. The SPM can not be even compared in the same way, since it produces only one $A-1$ fragment.

The differential cross sections (energy redistributions) of the fragments from the statistical fragmentation in the RDM and EM follow the relation in Eq. (2.5) (boost conservation). The redistribution for secondaries produced directly in the photon-nucleon interaction are sampled from Sophia (see Eq. (2.3)). A fully detailed description of the model’s implementation is located in Appendix A.

4 Impact in astrophysical scenarios

The characteristics of our new model impacts the production of UHECRs and neutrinos, depending on the parameters of the sources. We show in the following two representative example sources where the importance of the photomeson production have been already highlighted: GRBs [22] and TDEs [29]. The examples have been chosen such that photomeson processes dominate the disintegration at the highest energies.

GRBs are the most energetic electromagnetic outburst class. Shells of plasma emitted by the central engine can create internal shocks that become acceleration sites for charged particles. The interactions of those particles with the target photons in the so-called prompt phase can result in the production of a significant number of neutrinos. As a consequence, GRBs are candidate sources for UHECRs and neutrinos, although existing analyses already constrain some regions of the parameter space for GRBs [67] multi-messenger models.

In this work we use an example from LABEL:[22], where the interactions in the source are simulated with the NeuCosmA code for nuclei in the context of a one-zone model, meaning that the collistions of plasma shells cluster at the same collision radius $R$. Provided that baryons are present in the source, the density of the radiation triggers a nuclear cascade in which nuclei lighter than the primary are produced due to photodisintegration and photomeson production, impacting the mass composition of in-source and ejected cosmic rays. In LABEL:[22], several qualitative cases have been distinguished: for a high radiation density in the source, the nuclear cascade strongly develops (“optically thick case”) and produces a high flux of nucleons, which dominate the neutrino production. If the radiation density is low (“empty cascade”), the nuclear cascade does not develop and the neutrino flux is dominated by photomeson production of the primary nuclei. In the intermediate case (“populated cascade”), the nuclear cascade develops and neutrinos are efficiently produced off primary and secondary nuclei. The different cases can be quantitatively distinguished by the optical thickness to photohadronic interactions at the highest energies.

These different cases relating the degree disintegration to the multi-messenger production as introduced in LABEL:~[22] for GRBs, can be applied to other source classes [29]. One example are jetted Tidal Disruption Events (TDEs), which have been proposed as possible UHECR sources [4]. In this scenario, a star is gravitationally disrupted in the vicinity of a black hole by tidal forces, generating a jet in which a nuclear cascade similar the GRB can develop.

The source parameters used in the considered models are listed in Tab. 2. The luminosity $L$ and collision radius $R$ define the photon energy density in the source, together with the Lorentz factor of the shells $\Gamma$. The duration $T$ can vary a lot among different source classes, with consequences for the detection capability of different experiments. The baryonic loading $\xi$ is defined as the energy injected as baryons over the energy injected in photons, and it is fixed to a reference value here. The efficiency of the acceleration controls the maximum energy of acceleration. The spectrum of the GRB is typically described by a broken power law with a spectral break at 1 keV in the observer’s frame and spectral indices $\alpha=-2/3$ $(-1)$ and $\beta=-2$ below and above the break energy for TDEs (GRBs), respectively. Although the photodisintegration through excitation of the GDR is the dominant process for UHECRs with respect to photomeson processes, the interplay of the spectral index of the photons with the GDR and photomeson regimes may favor the photomeson production at the highest energies, as for example in TDEs. The minimal cutoff of the photon spectrum in the source can also drastically reduce the GDR interaction length, as demonstrated already in [22]. The critical parameters influencing photomeson production at the highest energies are marked boldface in the table.

Let us focus on the TDE scenario first, see Fig. 8. The slope of the photon spectrum at low energies reduces the photodisintegration rate at high energies (high energy protons interact with low energy photons), which means that the photomeson regime dominates the photohadronic interactions at the highest energies – as it is shown in the upper left panel of Fig. 8. Here the impact of the cross section systematics is only significant beyond the maximal energy of the cosmic rays, which can be estimated by balancing the acceleration rate with the sum of the energy losses (acceleration is possible only when losses are subdominant). In the upper right panel, the effect of the lower cross section adopted in the EM with respect to the SPM is visible: the injected primary nuclei (blue) is less depleted in the EM, therefore less energy is injected in secondaries. On the other hand the the production of secondaries is distributed over a larger range of masses (see Fig. 6) and the average mass of the secondaries is smaller than in the SPM, which implies a more efficient disintegration. The increase in the nucleons densities is the result of this more efficient disintegration of all the secondaries produced.

The neutrino production, see lower left panel, is mostly affected by the cross section scaling of the pion production; the heavier the primary, the larger the difference between the models. It is noteworthy that the SPM changes the qualitative observation that the neutrino production would be dominated by interactions of nuclei: In the EM, both contributions are similar in magnitude but nuclei dominate slighlty at the highest energies. Overall the neutrino flux is depleted by a factor of $~{}1.5$ if the EM is used with respect to the SPM. The effect of the disintegration channels is appreciable in the $<\ln A>$ where the EM leads to a smaller cascade mass (lower right panel). This happens mainly at the intermediate energies range because of the conservation of the Lorentz factor of the nuclei in the photodisintegration process, and because the secondaries are produced at lower energies in proportion to their mass ratio to progenitor species (due to boost conservation). At the highest energies, where the primary nuclei density dominates, the composition is slightly heavier than in SPM due to the reduced interaction.

For the GRB example, we follow LABEL:[22]. The luminosity and the collision radius reported in Tab. 2 correspond to a high radiation density, and the neutrinos are mainly given by the secondary nucleons produced in the nuclear cascade. In particular, we use the assumptions of Appendix B in LABEL:[22] for the photon spectrum of the source, in which a higher low-energy cutoff of the photon spectrum is investigated resulting in an suppression of photodisintegration rate with respect to the photomeson rate. The upper left panel of Fig. 9 demonstrates how the changes in the cross section affect the interaction rates at energies above the maximum energy of acceleration, which means that the effects on the energy densities (upper right panel) are dominated by the redistribution channels of the secondaries. The impact on the more efficient nuclear disintegration (through a larger number of open channels) is also visible in the $\langle\ln A\rangle$ behavior (bottom right) which overall exhibits a lighter composition, reaching up to three times lower value compared to the SPM. The neutrino production (bottom left) is again dominated by the pion production scaling used in EM ($\rm A^{2/3}$) in contrast with SPM scaling ($\rm A$), and the effect is large because of the heavier primary mass compared to the TDE example.

Fig. 10 shows the impact of the change of the photomeson model in the neutrino fluxes, corresponding to different GRB luminosities identified with different nuclear cascade regimes. Note that, compared to the main text of LABEL:~[22], here the low energy photon cutoff has been applied to all examples. The largest differences are visible in the low and intermediate luminosities, for which the neutrinos are mainly produced off (primary and secondary) nuclei. In the populated-cascade scenario (center), the smaller neutrino yields from the intermediate nuclei in the EM, change the leading production channel of the neutrinos. When using SPM in this intermediate case, the neutrino emission from secondary nuclei dominates the flux, whereas for the EM the contribution from nucleons is comparable to that from secondary nuclei.

5 Summary and conclusions

We have studied photomeson production off nuclei in astrophysical environments, where photomeson interactions refer to photon-nucleus collisions at photon energies $\gtrsim 140$ MeV in the nucleus’ rest frame. Above this threshold, the production of secondary mesons commences, and hence the production of astrophysical neutrinos in the multi-messenger context. Our main focus has been the treatment of the nuclear interactions in astrophysical environments for which the interacting nuclei have large kinetic energies resulting in forward peaked the secondary spectra. We have scrutinized a commonly adopted approach in the literature known as Single Particle Model (SPM), and introduced an improved description: the Empirical Model (EM). The impact of the EM was studied in astrophysical scenarios in which the photomeson processes are known to dominate at the highest energies (certain classes of GRBs and TDEs).

The SPM treats the nucleus as superposition of nucleons scaling the nucleus the total inelastic cross section with $A$ or $A^{2/3}$. We have discussed the disagreements of this model with available data, and improve those in the EM by using: a data-driven parametrization for the total inelastic cross section at low energies and a mass number scaling consistent with data and theory for high energies; a data-driven parametrization of the pion photoproduction cross section resulting in a $A^{2/3}$ scaling at the lower energies that accounts for nuclear medium effects and final state interactions; a nuclear breakup model for the remnant nucleus based on existing empirical formulae for partial photo-emission mechanisms, resulting in a better agreement with distributions obtained from Monte-Carlo simulations, such as Fluka.

These modifications affect the photomeson interactions off nuclei in astrophysical environments, and we used two examples from the literature in which the photomeson regime is known to dominate the nuclear interactions at the highest energies. Our GRB example resembles an effect from synchrotron self-absorption in which low energy target photons are suppressed that would otherwise disintegrate the UHECRs via Giant Dipole Resonance excitation. Our TDE example has a broken power law target photon spectrum with a hard enough (low) power law index such that the photomeson production dominates at the highest energies.

For these cases, we demonstrated that the improved model affects the nuclear cascade in the source resulting in an ejected UHECR composition that is up to three times lighter, and a reduction of the neutrino flux by up to 50%. The nuclear cascade is affected by the cross section (such as in the TDE case) or the additional channels in fragmentation (such as for our GRB case). The impact on the neutrino flux stronger (as in SPM) depends on the mass of the isotope(s) dominating the neutrino production, implying that neutrino emission becomes more sensitive to the choice of the injection composition into an internal-shock GRB model and the degree of the nuclear cascade. For high radiation densities in the optically thick case, neutrino production is dominated by nucleon interactions, and hence the impact from the new model is low. Consequently, the strongest effects will occur for the populated and empty nuclear cascade cases with heavy injection isotopes, which, however, have a smaller neutrino production efficiency compared to cases in which the cascade fully develops.

The impact on the UHECRs is smaller and only expected in environments where the photomeson regime dominates the nuclear disintegration. A prominent counter-example is, for instance, cosmic ray transport in the extragalactic space. In a rigidity-dependent model for the spectra of different nuclear species emitted from the accelerating source, it is found that the UHECR spectrum and composition are best-fitted by a low rigidity cutoff [68, 69, 70], meaning that the photomeson production is not triggered and the nuclear breakup occurs predominantly via photodisintegration off the extragalactic photon backgrounds. The corresponding cosmogenic neutrino flux at the main peak is dominated by interactions of the UHECRs with the cosmic infrared background (because of the low rigidity), and the contribution of nucleon and nuclei interactions in LABEL:~[70] at the peak is about 50-50. This means that the neutrino flux may be even lower (by up to about 25%) for our photomeson model.

We finally note that our model is based on a very small number of parameters and can be easily implemented for cross comparisons in any of the existing codes and calculation frameworks for radiation modeling of High Energy Particle Astrophysics. The tools for reproducing the model are available [42] and can be used for implementation in any framework simulating the interaction of UHECRs.

Acknowledgments

We would like to thank J. Rachen, J. Heinze and A. van Vliet for useful discussions.

This project has received funding from the European Union’s Horizon 2020 research and innovation programme under grant agreement no. 646623.

Appendix A Photomeson model based on empirical formulas

To construct the photomeson model it is necessary to obtain for each nuclear species $j$:

• •

$\sigma_{j}=\sigma_{j}(\epsilon_{r})$ the absorption photonuclear cross section as function of the photon energy in the nucleus rest frame $\epsilon_{r}$.

• •

$\sigma^{\rm incl}_{j\to k}=\sigma^{\rm incl}_{j\to k}(\epsilon_{r})$ the inclusive photonuclear cross section for producing each of the $\rm k$-th possible product particles (pions, nucleons, nuclear fragments).

The form of $\rm\sigma_{j}$ in the EM model is presented in Section 3.1 and summarized in Tab. 1. The expressions for the inclusive cross sections are presented in Section 3.2 for pions and in Section 3.3 for nucleons and nuclear fragments. The following sections detail the calculation of $\sigma^{\rm incl}_{j\to k}$ for nuclear fragments. The energy redistribution of the secondaries is calculated using Eq. (2.5) for the fragments of the nuclear remnant, while for direct interactions of the photon with a nucleon it is obtained from Sophia using Eq. (2.3) (see Appendix A.2).

A.1 Formulas for the inclusive cross sections

The inclusive cross sections for producing nucleons and nuclear fragments are based on the formulas discussed in sections below. The following relations are used in our model:

[TABLE]

where the expressions for $\sigma^{\rm dir}_{\rm p}$, $\sigma^{\rm dir}_{x{\rm n}}$, $\sigma^{\rm mul}_{\rm n}$, and $\sigma^{\rm sp}_{k}$ are obtained using empirical relations, derived from data in the energy range $\rm 0.2-1\,GeV$ in LABEL:~[63]. Here we adapt the naming convetion and remark that these relations do not correspond to the microscopic nuclear processes, but instead give a reasonable representation of data.

For simplicity and due to the sparsity of data, we assume they hold constant for all photon energies in the photomeson regime (see Fig. 1). Below are the processes and their expressions reproduced from LABEL:~[63]:

• •

Direct proton production ($\sigma^{\rm dir}_{\rm p}$) Reactions producing only one proton $(\gamma,\rm p)$ and a residual nucleus of mass $A-1$:

[TABLE]

• •

Direct neutron production ($\sigma^{\rm dir}_{\rm n}$) Reactions producing only one proton $(\gamma,\rm n)$ and a residual nucleus of mass $A-1$:

[TABLE]

• •

Multi-neutron production ($\sigma^{\rm mul}_{x{\rm n}}$) Reactions producing $x>1$ neutrons $(\gamma,x\rm n)$ and the corresponding residual nucleus with mass ${\rm A}-x\rm n$:

[TABLE]

• •

Spallation ($\sigma^{\rm sp}_{x{\rm p}y{\rm n}}$) Spallation reactions where a nominal loss of $x>1$ protons and $y>1$ neutrons occurs $(\gamma,x{\rm p},y{\rm n})$ and the corresponding residual nucleus with mass ${A_{r}=A}-x{\rm p}-y{\rm n}\geq A/2$ is produced:

[TABLE]

• •

Pion production ($\sigma^{\rm prod}_{\pi}$) Reactions producing pions ($\pi^{0},\pi+,\pi^{-}$) and one nucleon $(\gamma,\pi+\mathcal{N})$. This contribution is only used in the EM for normalizing the spallation in Eq. (A.10). The relations for pion production in the EM are discussed in Section 3.2:

[TABLE]

• •

Fission (not considered here) For the nuclei considered in this work (mass up to $\rm{}^{56}Fe$) this process has no values or experimental data [63] but needs to be included for much higher masses.

These relations are represented graphically in Fig. 11, except for the spallation curve which is strongly dependent on the mass and nucleons emitted, thus $\sigma^{\rm sp}$ is represented as the remainder of the total cross section $\sigma^{\rm tot}=0.28A\;{\rm mb}$ (see LABEL:~[63]) after subtracting all contributions

[TABLE]

To obtain the inclusive cross sections used in EM it is necessary to multiply the formulas above by the number of respective particles produced. For example, in the inclusive contributions from direct proton and direct neutron production ($\sigma^{\rm dir}_{\rm p}$ and $\sigma^{\rm dir}_{\rm n}$) only one nucleon and one large fragment ($Z-1$, $A$) or ($Z$, $A-1$) are produced. In the case of the multi-neutron process the cross section $\sigma^{\rm mul}_{x{\rm n}}$ for the emission of $x$ neutrons contributes to the inclusive cross section for neutron emission as $x\sigma^{\rm mul}_{x{\rm n}}$ (see Eq. (A.2)), i.e.:

[TABLE]

In the case of spallation, only one fragment species $k$ with mass $A/2\leq A_{k}\leq A-2$ is produced and the inclusive cross section is just $\sigma^{\rm sp}_{x{\rm p}y{\rm n}}$. We have made additional considerations to group the produced nucleons ($x$ protons and $y$ neutrons) into fragments with masses $A_{k}\in[1..4]$

[TABLE]

Where $\sigma^{\rm sp}_{l}$ stands for any light fragment produced ($\sigma^{\rm sp}_{\rm n}$, $\sigma^{\rm sp}_{\rm p}$, etc.). The method of estimating $\sigma^{\rm sp}_{l}$ is discussed in the following subsection.

A.2 Inclusive cross section of small fragments

In the processes considered in Section A.1 no fragments with masses smaller than $A/2$ are created. We assume for the EM that nucleons produced in the spallation process can be grouped into fragments of no more than four nucleons and estimated the number of those fragments with thermostatistical formulas [71].

For any spallation event $i$ a number of $x$ protons and $y$ neutrons is lost from the interacting nucleus ($Z$, $A$) with an exclusive cross section $\sigma_{x{\rm p}y{\rm n}}$ (which in the following will be labeled as $\sigma^{\rm sp}_{i}$). We consider that the total energy of left from the interaction is taken by the spalled nucleons as their kinetic energy (no internal excitation of products). The spalled nucleons can be configured in a number of small fragments, or a combination $r$:

[TABLE]

where $l$ refers to any of the nuclear species with no more than four nucleons and with a decay half life longer than the relevant timescale of the astrophysical (only common isotopes remain).

The set of all possible $\mathbb{C}_{i}=\{{\rm C}^{r}_{i}\}$ is determined by finding all mixtures of species $l$ in the quantities ${\rm N}^{r}_{i,l}$ such that the proton and nucleon number of the combination matches the spalled numbers:

[TABLE]

An appropriate weight for each combination ${\rm P}_{i,r}$ allows to find the numbers of particles produced in each spallation event $i$ with

[TABLE]

By summing the numbers of particle $l$ over all spallation configurations, the total inclusive cross section of spallation for particle $l$ is obtained

[TABLE]

The following section details the estimation of the combination weights ${\rm P}_{i,r}$.

A.3 Evaluation of the weights of combinations

The simplest assumption is that all combinations are equiprobable, that is they all have the same weight ${\rm P}_{i,r}={\rm P}_{i}$. Then the number of particles of species $l$ (Eq. (A.16)) results

[TABLE]

and the ${\rm P}_{i}$ can be found from the nucleon conservation. This means that the number of nucleons produced in all spallation events is equal to the nucleons present in the interacting nucleus

[TABLE]

where $\sigma^{\rm sp}_{k}$ are those in Eq. (A.12) and Eq. (A.17) with Eq. (A.18).

The equiprobability assumption fails to account for the fact that combinations with more stable nuclei are more likely to occur. In the EM we apply statistical mechanics, assuming that the combinations are possible microstates corresponding to a certain spallation event (macrostates) within the Grand Canonical distribution. A general form of the partition function in this case is [71]

[TABLE]

where $\epsilon_{i}$ is the energy of the macrostate, $\mu$ is the chemical potential of the elements (here only one species is included) and $n_{r}$ the number of elements of a given microstate. The weight of a microstate can be expressed by the probability of the microstate

[TABLE]

We are interested in the relative weights of different microstates with the same energy $\epsilon_{i}$ since we use the normalization condition Eq. (A.19). Hence, the expression for ${\rm P}_{i,n}$ can be rewritten grouping the factors common to a macrostate $i$ into some normalization constant $\mathsf{B}_{i}$

[TABLE]

which allows estimating the weights of the macrostates without calculating the partition function. Considering each nuclear species $l$ as a different constituent of the system, the chemical potential associated will be its ground state energy $\mu_{l}=m_{l}$ (rest mass in units of energy) which can be found in nuclear mass tables. Then the required expression of ${\rm P}_{i,r}$ is

[TABLE]

where the normalization constants $\mathsf{B}_{i}$ are found using Eq. (A.19). The expression chosen for $\beta_{i}$ is commonly used for nuclei [71]

[TABLE]

where the excitation energy $\bar{\epsilon}^{*}$ was chosen as the average photon energy in the range where the spallation formula was obtained

[TABLE]