Ultra-High-Energy Cosmic Rays from a Population of Non-identical Sources

. Astrophysical candidates for the sources of ultra-high-energy cosmic rays (UHECRs) exhibit a large diversity in terms of their properties relevant for the acceleration of charged particles, such as luminosity, Lorentz factor, size and magnetic field. Yet, fits of the observed UHECR spectrum and composition often assume identical sources. Here we investigate a population of sources with a power-law distribution of maximum energies. We show that the allowed source-to-source variance of the maximum energy must be small to describe the UHECR data. Even in the most extreme scenario, with a very sharp cuto ff of individual source spectra and negative redshift evolution of the accelerators, the maximum energies of 90% of sources must be identical within a factor of three – in contrast to the variance expected for astrophysical sources. However, the overall population variance can be large when maximum rigidities are distributed as a broken power law, with a steep decline above the break and with hard source spectra. In this scenario, most of the observed UHECR flux is produced by sources near the break.


Introduction
With energies up to and above 10 20 eV, ultra-high-energy cosmic rays (UHECRs) provide exciting prospects for identifying and understanding the most extreme environments in the Universe. However, although tremendous progress has been made since their discovery more than a century ago -and accelerated in the last decades by the commissioning of the latest generation of largescale UHECR observatories, the Pierre Auger Observatory (Auger) and Telescope Array (TA) -central questions about their origin and acceleration mechanism remain to be answered, and individual sources have so far avoided detection at appreciable significance.
A popular approach for constraining the properties of the UHECR sources is based on the combined fit of the observed spectrum and composition at Earth and crossreference with the predictions made by "effective" models that attempt to describe the source population with a limited number of free parameters (e.g. [1][2][3][4][5][6][7]). While these studies have identified some common features present over most analyses, significant assumptions are made in the process. Arguably the strongest of these is that all sources are identical in terms of the parameters relevant for cosmic rays acceleration, such as maximum rigidity 1 , spectral index and emitted cosmic-ray composition.
It is standard to assume an acceleration mechanism universal in rigidity up to some maximum rigidity R max . In this "Peters Cycle" [8] scenario, a remarkably good fit to the observed spectrum and composition can be achieved even with these simplified models if one accepts injection spectra harder than dN/ dE ∼ E −1 and heavy source compositions dominated by CNO-Fe nuclei at the highest energies [3][4][5][6][7]. However, the assumption of identical sources remains questionable.
Motivated by the observational fact that the most likely astrophysical objects where acceleration is predicted to occur, such as active galactic nuclei (AGN), gamma-ray bursts (GRBs) and tidal disruption events (TDEs), exhibit significant diversity in important parameters such as luminosity, magnetic field strength and Lorentz factor, we extend the simplistic model of identical sources by allowing for source diversity in terms of maximum CR rigidities. In particular, we investigate the case of a source population with a power-law distribution of maximum rigidities.
In these proceedings, we summarise our investigation of such non-identical sources and the constraints on source diversity imposed by current UHECR observational data. The extended version of our study and conclusions is available at [9]. We find that observations of UHECR spectrum and composition require their sources to be nearly identical for realistic source parameters.

Population Spectrum of non-identical Sources
Following the Peters cycle prescription, we assume that the cosmic-ray spectrum of an individual source and injected nuclear species i can be described by a power law with cutoff at high rigidities The flux normalisation factors ϕ 0 (Z i ) are the averaged normalisation of the entire population and can already include a hidden variance of the per-source normalisations. They are the same for all sources. The population spectrum, i.e. the cosmic-ray flux produced by the entire population of, possibly non-identical, sources can be obtained from the individual source spectra ϕ src and the distribution of maximum rigidities p(R max ) as The sum over all injection elements is assumed implicitly. If sources are identical, then p(R max ) = δ(R 0 − R max ) and the population spectrum will exhibit the same shape as the individual source spectra. Several different choices can be made for the highrigidity cutoff. We adopt as fiducial scenario a commonlyused exponential cutoff that is expected in some astrophysical contexts [10,11]. It is obtained from a more general super-exponential solution in the limit of λ cut = 1. We also consider other values that will result in less (λ cut < 1) or more (λ cut > 1) intrinsic variance of the maximum rigidity. Assuming a power-law distribution of maximum rigidities and source spectra with super-exponential cutoff, Eq. (2) can be evaluated analytically and the predicted population spectrum reads where γ(. . . ) is the lower incomplete gamma function. The source spectra and resulting population spectra for a range of different cutoff choices are shown in Fig. 1. As evident from the figure and provable analytically, the asymptotic behaviour is entirely independent of the cutoff choice, and we always retrieve

Astrophysical Connection
The simple parameterisation given in the previous section can be connected to observable distributions of specific properties of astrophysical objects to gain further insights into the nature of UHECR sources. In the following, we propose two such connections, one with the Lorentz factor of jetted acceleration environments and one with the observed source luminosity. Figure 1. Top: Cosmic-ray source spectra for different highrigidity cutoffs. R max denotes the maximum rigidity, and the yaxis is scaled to show the ratio to an unmodified power law with spectral index γ src . Bottom: Population spectra resulting from the convolution of a power-law distribution of maximum rigidities above rigidity R 0 (for β pop = 4) and the source spectra displayed in the previous panel.

Lorentz Factor
In the simplified picture where acceleration is limited only by the size of the source -the Hillas criterion [12]for UHECR production in relativistic jets, the maximum rigidity is connected to the bulk Lorentz factor as R max = R 0 Γ jet , with R 0 ∝ l B related to the size l of the source region and the magnetic field B. More generally, this scenario can be expressed as where the Hillas scenario is captured by α = 1 and Espresso acceleration [13][14][15] imposes α ≤ 2.
Assuming a power-law distribution of Lorentz factors dp dΓ jet in agreement with observations of jetted AGN [16,17], this results in the distribution of maximum rigidities p(R max ) = dp dΓ jet dΓ jet dR max (10) Including the effect of energy-boosting on the total output flux of a particular source, the above expression becomes the same as our simple ansatz Eq. (5) after identifying the slope of the R max -distribution as with a time dilation factor ξ = 1 for UHECR production and acceleration inside the jet and ξ = 0 for espresso-style re-acceleration.

Luminosity
Alternatively, we connect the maximum rigidity to the source luminosity by using the Lovelace-Waxman-Blandford argument for the minimum luminosity required for cosmic ray acceleration in expanding flows [18][19][20][21][22] L 0 ≈ 10 45. 5 with β the bulk velocity of the outflow. The maximum rigidity as a function of luminosity is then The additional source diversity introduced by nonidentical outflow speeds is expected to be small and can be neglected. For a power-law luminosity function with slope y 2 , as observed for a large range of potential source classes (e.g. AGN [23,24], TDEs [25,26]), dp dL = our original power-law description of p(R max ) is reproduced exactly, except for an additional normalisation factor κ, for

Methods
We simulate the injection and propagation of UHECRs numerically with CRPropa3 [27] and compare the predictions with the latest available spectrum [28] and composition data [29,30] from the Pierre Auger Observatory. For our fiducial model the free parameters are the minimum maximum rigidity of sources R 0 , spectral index γ src and slope β pop of the R max -distribution. Other parameters, such as the redshift evolution of the source emissivity n(z) and of R 0 , and the shape of the UHE cutoff λ cut are also investigated and mentioned briefly. For a full discussion of these variations, we refer the reader to the primary paper [9]. Sources are assumed to follow a continuous distribution with number density constant in redshift up to a maximum redshift z max = 4. Except for the case of R 0 evolving strongly with redshift, existing multimessenger limits and Figure 2. Source parameter scan for the fiducial model marginalised along all but two axes respectively. The surface plot shows the agreement between prediction and Auger observations in terms of the χ 2 estimator, and the contour lines indicate the one (green) and three (red) sigma confidence intervals for two degrees of freedom. The best fit is marked with a white cross.
observations do not constrain the source parameters in our model, and we omit their discussion for brevity.
The propagation is treated ballistically, and magnetic fields are not included. For the extragalactic background light, we use the model of Gilmore et al. [31]. To limit the complexity of the simulation and the number of adjustable injection fractions, we choose five representative elements ( 1 H, 4 He, 14 N, 28 Si and 56 Fe), which are commonly used for similar fits and provide adequate coverage over the entire mass range.
Auger observations indicate a small variance of the mean shower depth at UHE, implying a relatively pure composition. A spread in maximum rigidity within the source population will increase the mixing between different elements even if the composition at the sources is identical. To establish an upper limit on the allowed population variance, it is necessary to minimise the intrinsic shower variance, which can be achieved by assuming the heaviest UHECR composition compatible with observations since heavy CRs produce less variance in the air showers. To that end, we select Sibyll2.3c [32] as hadronic interaction model for our fiducial scenario and shift the mean and variance of the shower depth to their heaviest realisation, i.e. −1σ syst & + 1σ syst respectively.

Results and Discussion
In contrast to the sizeable diversity of sources expected from the astrophysical connections, we find that sources with power-law distribution in R max must be effectively identical in terms of maximum rigidity if the un-shifted Auger data is considered. Even after including the fiducial shifts of the mass scales, the best-fit population variance is limited, and the softness of the R max -distribution is bounded as (see Fig. 2) This is a consequence of the high-rigidity limiting behaviour of the population spectrum Eq. (6). If source spectra are soft and the distribution of maximum rigidities is relatively flat then ϕ pop is characterised by a large tail toward the highest CR energies even beyond the GZK limit. These extremely-UHE cosmic rays will experience strong interactions and rapid disintegration, resulting in a very mixed cosmic-ray flux at Earth. This is in disagreement with the observed pureness of the UHE composition, and the corresponding source configurations are consequently disfavoured. From the above condition, significant source diversity is only possible for γ src ≳ 1, which is again disfavoured because of the associated, relatively mixed CR flux. At the best fit (Fig. 3), where β pop ≈ 5.2, sources are functionally identical, and population variance exceeding half a decade for 90% of sources is excluded at more than 6σ. Quasi-identical "standard candle" sources (β pop ≳ 10) cannot be rejected.
Varying the redshift evolution of the source number density, n(z) ∝ (1+z) m , we find that the allowed population variance is increased marginally for highly local distributions while sources are required to be even more similar if they are further away on average. This can be explained intuitively by the increased interaction and fragmentation of cosmic rays from more distant sources, which results in more intrinsic mixing between the different mass groups, reducing the permitted inter-source variance. More variance is also possible if the UHE cutoff of the source spectra is steeper than exponential, which again reduces the intrinsic mixing of different nuclei and leaves more room for source-to-source variance. At λ cut = 5.4 +1.7 −2.3 , the best fit is close to a Heaviside-like spectral cutoff.
Even in the most extreme scenario, with Heaviside source spectrum cutoff and negative density evolution (m = −3), the distribution of maximum rigidities cannot be flatter than β pop ≈ 3 − 4; corresponding to a factor of only two to three difference in the maximum rigidity of 90% of sources.
A larger overall population variance is possible if a broken power-law distribution of maximum rigidities is assumed instead. In this case, good agreement with observations can be achieved if the individual source spectra are hard and the R max distribution decays rapidly above the break. The observed UHECR flux is then dominated by sources around the break while the population variance is driven by sources far below the break.

Conclusions
We have performed an analysis of the dispersion in maximum rigidity of a population of cosmic ray sources allowed by current observations of the UHECR spectrum and composition. We have derived analytic expressions for the expected population spectrum for a range of highenergy cutoff functions, assuming a power-law distribution of maximum rigidities.
With these parametrisations as input to our numerical simulation of UHECR propagation, we determine that the sources responsible for the observed UHECR flux must be effectively identical in their maximum rigidity. If the most optimistic, i.e. 'heaviest', mass scale of the observational data is applied, the allowed population variance is increased slightly. However, even in the most extreme scenario of sudden spectral cutoff and predominantly nearby sources, we find that the majority of sources cannot differ by more than a factor of three in maximum rigidity. This is in stark contrast to the population diversity expected for the most probable source classes if a connection between maximum rigidity and jet Lorentz factor or source luminosity is considered.
Only if a broken power-law distribution in maximum rigidity is considered can the total population variance become large, provided the distribution decays rapidly above the break and source spectra are hard.