Mass composition of cosmic rays above 0.1 EeV by the Yakutsk array data

TL;DR

This study analyzes the mass composition of ultra-high-energy cosmic rays above 0.1 EeV using Cherenkov radiation data from the Yakutsk array over four decades, revealing a change in composition with energy.

Contribution

It provides a detailed analysis of cosmic ray mass composition at ultra-high energies based on extensive Cherenkov radiation observations and models, extending previous measurements.

Findings

Displacement rate dXmax/dE is nonlinear and energy-dependent.

Evidence of changing cosmic ray composition with energy.

Fluctuations of Xmax support composition variation.

Abstract

The paper presents the results of the longitudinal development of extensive air showers (X$_{max}$) of ultra-high energies and mass composition of cosmic rays. The measurements of X$_{max}$ are based on data from observations of the Cherenkov radiation at the Yakutsk array for the period 1974-2014. The cascade curves of individual showers and the depth of maximum X$_{max}$ were reconstructed over the energy range 10$^{16}$-5.7$\cdot$10$^{19}$ eV. It is shown that the displacement rate of the parameter dX$_{max}$ / dE in the atmosphere is nonlinear and depends on the energy. Such a feature indicates a change in mass composition, which is confirmed by fluctuations of X$_{max}$ in this energy region. The composition of cosmic rays was determined by interpolation using the QGSJetII-04 model.

Tables

Table 1. Table 1: Different values for P λ and δ λ subscript 𝛿 𝜆 \delta_{\lambda} in different scales of the atmosphere transparency.

Atmospheric transparency	Pλ=530	${\delta }_{\lambda =530}$	Pλ=369	${\delta }_{\lambda =369}$	Pλ=430
low (2 points)	0.60	0.180	0.35	0.260	0.48
decreased (3 points)	0.70	0.080	0.46	0.120	0.56
normal (3 - 4 points)	0.75	0.060	0.50	0.080	0.60
increased (4 points)	0.80	0.040	0.52	0.060	0.68
high (5 points)	0.85	0.020	0.56	0.030	0.72

Table 2. Table 2: ⟨ X m ​ a ​ x ⟩ delimited-⟨⟩ subscript 𝑋 𝑚 𝑎 𝑥 \left\langle X_{max}\right\rangle , σ 𝜎 \sigma ( X max ) and ⟨ l ​ n ​ A ⟩ delimited-⟨⟩ 𝑙 𝑛 𝐴 \left\langle lnA\right\rangle observed in 1974-2014 at the Yakutsk array.10 16 - 10 18 eV

E0, eV	N, events	Xmax, g/cm2	$\sigma$(Xmax), g/cm2	$⟨lnA⟩$	$⟨lnA⟩$ stat. error
1.74E+16	-	589$\pm$17	-	1.97	0.68
2.62E+16	-	591$\pm$14	-	2.25	0.55
3.50E+16	24	598$\pm$11	53.6$\pm$7.7	2.23	0.43
5.31E+16	50	606$\pm$8	56.7$\pm$5.7	2.35	0.32
7.01E+16	88	608$\pm$6	56.3$\pm$4.2	2.50	0.24
8.85E+16	139	617$\pm$5	58.9$\pm$3.5	2.34	0.20
1.07E+17	158	628$\pm$5	62.8$\pm$3.5	2.17	0.20
1.24E+17	146	631$\pm$5	60.5$\pm$3.5	2.17	0.20
1.42E+17	155	636$\pm$5	62.2$\pm$3.5	2.09	0.20
1.60E+17	158	639$\pm$5	62.9$\pm$3.5	2.07	0.20
1.78E+17	153	646$\pm$5	61.9$\pm$3.5	1.97	0.20
2.68E+17	155	655$\pm$5	62.3$\pm$3.5	1.95	0.25
3.57E+17	156	672$\pm$5	62.4$\pm$3.5	1.48	0.21
5.43E+17	161	688$\pm$5	63.5$\pm$3.5	1.28	0.21
7.18E+17	150	695$\pm$5	61.2$\pm$3.5	1.30	0.22
9.06E+17	155	703$\pm$5	62.4$\pm$3.5	1.18	0.22

Table 3. Table 3: ⟨ X m ​ a ​ x ⟩ delimited-⟨⟩ subscript 𝑋 𝑚 𝑎 𝑥 \left\langle X_{max}\right\rangle , σ 𝜎 \sigma ( X max ) and ⟨ l ​ n ​ A ⟩ delimited-⟨⟩ 𝑙 𝑛 𝐴 \left\langle lnA\right\rangle observed in 1974-2014 at the Yakutsk array. 10 18 ⋅ ⋅ \cdot 5.7 ⋅ ⋅ \cdot 10 19 eV

E0, eV	N, events	Xmax, g/cm2	$\sigma$(Xmax), g/cm2	$⟨lnA⟩$	$⟨lnA⟩$ stat. error
1.09E+18	148	715$\pm$5	60.9$\pm$3.5	0.87	0.22
1.27E+18	143	721$\pm$5	59.8$\pm$3.5	0.73	0.22
1.45E+18	142	725$\pm$5	59.6$\pm$3.5	0.66	0.22
1.64E+18	93	734$\pm$6	58.0$\pm$4.3	0.47	0.25
1.82E+18	92	735$\pm$6	57.4$\pm$4.2	0.47	0.26
2.74E+18	89	749$\pm$6	56.5$\pm$4.2	0.30	0.26
3.66E+18	86	762$\pm$6	55.7$\pm$4.2	0	0.26
5.56E+18	78	769$\pm$6	52.9$\pm$4.2	0.09	0.26
7.35E+18	53	774$\pm$7	50.8$\pm$4.9	0.17	0.30
9.27E+18	47	777$\pm$7	48.1$\pm$5.0	0.34	0.30
1.12E+19	48	784$\pm$7	48.7$\pm$5.0	0.25	0.30
1.28E+19	32	782$\pm$8	45.5$\pm$5.7	0.39	0.35
1.48E+19	35	784$\pm$8	47.2$\pm$5.6	0.48	0.35
1.66E+19	26	782$\pm$9	45.5$\pm$6.3	0.69	0.39
1.86E+19	23	785$\pm$9	42.8$\pm$6.3	0.73	0.39
2.80E+19	11	787$\pm$13	42.9$\pm$9.1	0.99	0.56
3.74E+19	7	783$\pm$15	39.7$\pm$10.6	1.57	0,64
5.69E+19	3	790$\pm$23	41.1$\pm$16.8	1.69	0.97

Equations

$$Q_{exp}={\delta }_Q+{\int }_{X_1}G(R,X/X_2)\cdot N(E_0,X)\cdot K(\lambda ,X)dX$$

$$E_0=E_{ei}+E_{el}+E_{\mu }+E_{hi}+E_{\mu i}+E_{\nu }$$

$$E_{ei}=k(x,P_{\lambda })\cdot \Phi$$

$${\sigma }_{Q_i}^2=\left(0.04+n^2\cdot {\left(\frac{\Delta R}{R}\right)}^2\cdot Q_i^2\right)$$

$$E_0=(1.78\pm 0.44)\cdot 10^{17}\cdot {\left(\frac{Q(200)}{10^{17}}\right)}^{1.01\pm 0.04}$$

$$⟨\sigma (X_{max}{)}_{app}⟩=(38.5\pm 5)-(10\pm 3)\lg \frac{E_0}{10^{18}}$$

$${\sigma }^2(X_{max})={\sigma }^2(X_{max}{)}_{meas}-{\sigma }^2(X_{max}{)}_{app}$$

$$⟨lnA⟩=\frac{X_{max}^{exp}-X_{max}^p}{X_{max}^{Fe}-X_{max}^p}\cdot ln{A}_{Fe}$$

Full text

mass composition \sepcosmic rays \sepCherenkov light \sepultra-high energies \sepYakutsk array

I Introduction

Knowledge of the mass composition (MC) of cosmic rays (CR) in the energy range above 0.1 EeV is very important for astrophysics. It is assumed that a transition from galactic to extra-galactic CR (Berezhko et al., 2012) is within the energy range of 0.1-1 EeV, also irregularities like ”dip-bump” at energies above 5 EeV (Berezinsky et al., 2006), and cutoff of the spectrum of CR associated with the GZK effect at energies $\sim$50 EeV (Greisen, 1966; Zatsepin and Kuzmin, 1966).

One of the crucial techniques in the study of the characteristics of CR by Extensive Air Shower (EAS) method is Vavilov-Cherenkov radiation (Zrelov, 1968; Jelley, 1958) detection, which is produced by shower particles in the visible wavelength. The total flux of Cherenkov light allows determining the energy of the primary particle and the lateral distribution at sea level provides information on the nature of the development of a shower in the atmosphere (Knurenko et al., 2001).

Longitudinal development of the shower in the atmosphere was reconstructed by the lateral distribution function (LDF) of Cherenkov light (Knurenko et al., 2001; Dyakonov et al., 1986) using the inverse problem method (Tikhonov and Arsenin, 1977). To estimate the mass composition, the depth of maximum of air shower development Xmax was used, which is sensitive to the atomic mass of the primary particle.

At the Yakutsk array, in order to estimate air shower energy we used the energy balance method (Knurenko et al., 2006a). The method is based on the energy dissipated by EAS particles in the atmosphere, i.e. full flux of Cherenkov light. In individual showers, the energy was determined by the flux density of Cherenkov light at a distance R = 200 m from the shower axis - Q (200). According to calculations, the classification parameter Q (200) is proportional to the shower energy and have little shower-to-shower fluctuations.

The average depth of the air shower maximum Xmax was used to derive CR MC within the hadronic interaction model QGSJetII-04 (Ostapchenko, 2011) according to (Hörandel, 2006).

The Yakutsk array consists of stations with scintillation and Cherenkov detectors to register elementary particles of different nature: hadrons, electrons and positrons, muons, and Cherenkov photons (Artamonov et al., 1994). Additionally, there are surface electric field and air shower radio emission observations (Knurenko et al., 2017). The array consists of 120 scintillation detectors located on an area of $\sim$13 km2. Three muon stations and 72 channels for integrated and differential Cherenkov detectors are located at the center (Knurenko et al., 1998; Garipov et al., 2001; Ivanov et al., 2015; Egorov et al., 2018) (Fig. 1). The observation energy range of the array is from 1015 eV to 1020 eV, which allows it to study both galactic and extra-galactic CR.

Cherenkov measurements are conducted on clear and moon-less nights. The observation season starts in September and ends in April. The total observation time per season is 450-600 hours, which takes 6-10$\%$ of the annual array operation time. Since optical observations are very sensitive to such atmospheric characteristics as transparency and the aerosol component of the atmosphere (Dyakonov et al., 1999; Zuev, 1970), the atmosphere is monitored daily at the array (Dyakonov et al., 1991; Knurenko et al., 2006b; Knurenko and Petrov, 2014).

II Experimental data and analysis

II.1 Reconstruction of air shower parameters

The depth of maximum Xmax is determined by the cascade curve of Cherenkov light, which is reconstructed from the experimental Cherenkov LDF (Fig. 2) (Knurenko et al., 2001). This algorithm is described in detail in (Knurenko et al., 2001; Dyakonov et al., 1986). The algorithm is based on the Fredholm equation of the first kind (1), which is solved by the adaptive method (Kochnev, 1985):

[TABLE]

In equation (1), $G(R,X/X_{2})$ is the function that represents the lateral-angular distribution of electrons in the electron-photon cascade. In our case, the function $G(R,X/X_{2})$ takes into account the nonequilibrium spectrum of electrons (Belyaev et al., 1980; Dyakonov, 1982), i.e. depends on the age of the shower and is close to the data from Hillas (Hillas, 1982; Patterson and Hillas, 1983) and Giller (Giller et al., 2004). More recent work (Nerling et al., 2006; Smialkowski and Giller, 2018) confirms the correctness of the use of lateral-angular distribution that depends on the age of the shower, for Xmax reconstruction by Cherenkov light of air showers. The use of different models of lateral-angular distribution from the ”equilibrium” approximation and, including calculations of the electron spectra depending on shower age, yielded 3-5$\%$ in the spread of the Xmax estimation. N (E0, X) is a cascade curve, $\delta$Q is the ”noise” level, depends on measurement uncertainty and data processing errors. The choice of lateral-angular distribution was made in (Knurenko, 2003). K ($\lambda$, X) is transmittance of the atmosphere that takes into account the transmission of Cherenkov light throughout the entire thickness of the atmosphere. The coefficient K ($\lambda$, X) was reconstructed by solving the inverse problem using the lateral distribution of Cherenkov light registered under different weather conditions by a photomultiplier with a maximum sensitivity at a wavelength of $\lambda$ = 430 nm (Dyakonov et al., 1999). Reconstructed dependence of Cherenkov light transmission coefficient for different weather conditions on the height is shown in Fig. 3. Approximation of these data is also given there.

For the prompt estimation of atmospheric conditions during the registration of air showers in the optical wavelength at the Yakutsk array, a point system for estimation of the state of the atmosphere was developed. The estimation was tested by instrumental methods, including the rate of triggers for showers with energies of 1015-1016 eV (Dyakonov et al., 1991). The point system for classification of the atmosphere was proposed in the papers (Knurenko et al., 2006b; Guschin, 1988) based on lidar measurements of the spectral transparency Pλ and spectral atmospheric optical thickness $\delta_{\lambda}$. The correlation results for the three wavelengths are given in Table 1.

X1 and X2 - upper and lower limits of the atmosphere for air shower development. The limits were used for reconstruction of the air shower cascade curve.

As seen from equation (1), the method takes into account the physics of air shower electron-photon component development and characteristics of the atmospheric conditions during registration of the Cherenkov radiation (Knurenko et al., 2006b; Knurenko and Petrov, 2014) at the Yakutsk array. Using this technique, average cascade curves were reconstructed (Fig. 4) for energies of $\sim$1018 and $\sim$1019 eV (Tikhonov and Arsenin, 1977) and compared with QGSJetII-04 simulations. As can be seen from Fig. 4, the experimental cascade curve is measured in a wide range of depths due to inclined showers. This circumstance was used to determine the absorption range of cascade particles beyond the maximum development of air showers (Glushkov et al., 1993). On the other hand, cascade curves can be used to test different models of the development of air shower in the atmosphere.

The energy at the Yakutsk array was determined by energy balance method of all air shower particles (Knurenko et al., 2006a). This suggests that the total energy of the shower takes into account the ionization losses of each of the components of air shower: electrons, muons, hadrons and neutrinos (see formula (2)).

[TABLE]

Their share of energy was determined empirically taking into account the experiment conducted at the Yakutsk array. For example, the ionization losses of the electronic component were determined by the total flux of Cherenkov light from air showers, according to formula (3):

[TABLE]

where $\Phi$ is the total flux of Cherenkov light; k (x, Pλ) is the approximation coefficient (calculated value), takes into account the transparency of the real atmosphere, the nature of the longitudinal development of the shower (energy spectrum of secondary particles and its dependence on the age of the shower) and expressed through the depth of the maximum of air shower Xmax, measured at the array (Knurenko et al., 2006a; Ivanov et al., 2015).

As can be seen from Fig. 5, the fraction of Eem is 85-90$\%$, fraction of the remaining components of the shower does not exceed 10$\%$, and their estimation is described in details in papers (Knurenko et al., 2006a; Ivanov et al., 2015).

The experimental data of the Yakutsk array are in good agreement with the calculations by the QGSjetII-03 model in the case of a mixed composition of primary particles. This estimate may change if other models are used, including the QGSjetII-04 (Ostapchenko, 2011). This is stated in our work (Knurenko et al., 2005). In papers (Abbasi et al., 2018a; Barbosa et al., 2004), the calculations of Eem / E0 were performed to estimate the energy of air showers registered at the HiRes and TA arrays. Based on the calculations and the experiment in Yakutsk, it can be said that with an accuracy of 10$\%$, there is agreement on the energy estimates in both arrays.

A test check of the reconstruction method was carried out using QGSjet-01 calculations of the Cherenkov light of air showers. The result of comparing the cascade curve calculated by the model for the energy of 1018 eV and reconstructed by the method described above showed a good agreement Knurenko et al. (2001).

II.2 Air shower reconstruction uncertainty

Measurement errors in our experiment are systematic and random (Knurenko et al., 2015). The main source of systematic uncertainty is the absolute calibration error due to the uncertainty of the transition from the amplitude of the signal at the output of the PMT to the number of Cherenkov light photons. It is determined by the coupling coefficient between the amplitudes of the signals from the reference light source and the plastic scintillator unit, which serves as a calibration light source. The relative error of this coefficient is $\delta_{Q_{cal}}$ = 21$\%$.

Random measurement errors are due to many factors. One of them is the instrumental error of response measurement itself, equal to $\delta_{Q_{app}}$ = 13$\%$. For the short (R $<$ 100 m) distances from the axis, a similar contribution is made by axis determination method, the absolute error of which for showers selected by the central part of the array is $\Delta$R = 10 m for the approximation of Q (R) by a piecewise-power function with value of n $\sigma_{loc}$ = n$\cdot\Delta$ (R / Rthresh.). Along with this, there are Q (R) distortions due to the shower processing procedure. These include the use of different test functions with different boundary conditions for the density of charged particles, the classification of showers according to different measured parameters, the threshold of the Cherenkov light detector. The study of these distorting factors shows that in this case measurement errors of Q (R) increase at small and large (R$>$ 800 m) distances from the shower axis. The influence of the detector threshold begins to affect from a distance of R$>$ 600 m, depending on the energy of the primary particle. It is expressed in the fact that due to the random nature of the detector triggering at the threshold level, zero readings are underestimated and the value of the Cherenkov light flux begins to be systematically overestimated. For E0 = 2$\cdot$1018 eV, this overestimation is 5$\%$ and 20$\%$ at distances of 600 and 800 m, respectively, and this fact had to be taken into account when processing showers. For this reason, showers for the analysis were selected only by the central part of the array, where the symmetry condition of the location of observation stations relative to the air shower axis was met and all showers had approximately the same uncertainty for the reconstruction of the air shower characteristics.

Above mentioned uncertainties were taken into account for reconstruction of the cascade curve according to formula (4):

[TABLE]

Here $\Delta$R - is error of air shower axis determination error, 0.04 - constant associated with the absolute calibration of the light receivers and n - slope value of the power approximation of Cherenkov light LDF.

Uncertainty of the Xmax reconstruction for individual showers was 15-55 g/cm2 (Knurenko, 2003; Knurenko et al., 2001). Uncertainty of energy estimation at the Yakutsk array, taking into account the measurements uncertainties listed above is 23-26$\%$ (Knurenko et al., 2006a; Ivanov et al., 2009).

II.3 Depth of maximum development Xmax

The analysis presented in this paper is based on the Yakutsk array 1974-2014 dataset. The showers selected for analysis are shown in Fig. 6, as the dependence of Xmax on the classification parameter of air showers Q (200) - Cherenkov light flux density at the distance 200 m from the shower axis. The flux density of the Cherenkov light has little shower-to-shower fluctuations and weakly depends on zenith angle (Kalmykov, 1967; Knurenko, 2003), therefore Q(200) was used as the classification parameter for air shower selection. The energy of the shower is derived from Q (200) parameter by formula (5) (Knurenko et al., 2006a):

[TABLE]

In the distribution of air showers given in Fig.6, there are events with Xmax $\geq$ 800 g / cm2, i.e. with the low depth of maximum development compared to iron or even proton induced air showers. These showers can be used to search for neutral particles such as $\gamma$-rays and neutrinos (Knurenko and Petrov, 2018).

The averaged data of the Yakutsk array together with the data of the Auger (Bellido et al., 2018), TA (Abbasi et al., 2018b), HiRes (Abbasi et al., 2010; Ulrich and Cazon, 2012), Tunka (Prosin et al., 2016), LOFAR (Buitink et al., 2016) arrays are shown in Fig. 7 (a) and comparison with different hadronic interaction models QGSJetII-04 (Ostapchenko, 2011), Sibyll 2.3c (Riehn et al., 2015), EPOS LHC (Pierog et al., 2015). The data within the experimental errors are in good agreement with each other, although they were obtained at different experiments and by different methods. The data of the Yakutsk array covers almost four orders of magnitude in energy. Statistics allows with sufficient accuracy to identify the shift of Xmax with increasing energy and to estimate the elongation rate (ER) at different energy ranges. According to our data, the ER is:

E = 1016-1017 [eV]; E.R. = 48$\pm$6 [ g/cm2 ] 2. 2.

E = 1017-1018 [eV]; E.R. = 78$\pm$5 [ g/cm2 ] 3. 3.

E = 1018-1019 [eV]; E.R. = 63$\pm$6 [ g/cm2 ] 4. 4.

E = 1019-1020 [eV]; E.R. = 50$\pm$7 [ g/cm2 ]

The data $\langle X_{max}\rangle$, $\sigma$(Xmax) and other characteristics are given in Tables 2, 3.

In Fig. 7 (b) Xmax fluctuations are presented. These fluctuations were obtained taking into account detectors effect according to the method given below. The accuracy of the Xmax reconstruction was estimated from the simulation of cascade curve reconstruction process by full Monte Carlo method, given the known values of the instrumental uncertainty Q(E0, R, X0) and parameters of the average cascade curve. It was assumed that the measurement uncertainty Q(E0, R, X0) are distributed according to the normal law (Dyakonov et al., 1981, 1983). As follows from (Dyakonov et al., 1981, 1983) $\sigma(X_{max})_{app}$ depends not only on measurements uncertainty Q(E0, R, X0), but also on coordinates of the shower axis relative to the array center (X0, Y0) and the number of stations n participating in the registration of the Cherenkov EAS light. So the average instrumental error at an energy of 1018 eV is equal to $<\sigma(X_{max})_{app}>$ = 38 g / cm2, and its dependence on energy can be expressed as:

[TABLE]

Then physical fluctuations $\sigma(X_{max})$ were found as the difference of measured fluctuations and instrumental fluctuations:

[TABLE]

The data was divided by energy with a step of 1.5, and the value $\sigma$(Xmax) was found for each interval. The figure contains data obtained at the Yakutsk array and data of Auger, TA and HiRes. The figure also shows calculations using modern models of hadron interactions for the primary proton and the iron nucleus.

From Fig. 7 (b) one can see there is a maximum of Xmax fluctuations dependency on energy in the range 2$\cdot$1017-2$\cdot$1018 eV. In this range, the value of $\sigma$(Xmax) is 57-63 g/cm2, which according to model calculations corresponds to the mixed composition of cosmic rays with a high content of p + He and heavier nuclei. However, the composition consisting of a mix of 50% proton p and 50${\%}$ iron nuclei Fe cannot be excluded. A decrease in fluctuations above the energy of 5$\cdot$1018 eV indicates a heavier composition.

II.4 Mass composition

Knowing the depth of the shower maximum Xmax (Fig. 7, a) and average Xmax for the proton and the iron nucleus from QGSJetII-04 simulations, the value of $\langle lnA\rangle$ can be determined by interpolation using formula (8) (Berezhko et al., 2012; Hörandel, 2006):

[TABLE]

where X${}^{exp}_{max}$ is the depth of maximum determined from the experiment, $lnA_{Fe}$ is the natural logarithm of the iron atomic mass.

The mean logarithmic mass values are presented for different energies in Fig. 8. The results are compared with data of Tunka (Prosin et al., 2016), LOFAR Buitink et al. (2016), TA (Abbasi et al., 2019) and Auger (Bellido et al., 2018), derived using the hadronic interaction model QGSJetII-04. There is also the Yakutsk array data for 1994-2010 (Knurenko and Sabourov, 2011)]. This result was obtained for limited data selection with QGSJETII model, which explains discrepancy with the data presented in current paper As can be seen from the figure, there are irregularities of the value of $\langle lnA\rangle$, due to MC change. For example, in the energy region of 1016-1017 eV, the composition is heavier than in the energy region of 1017-1018 eV, where composition mostly consists of protons p and helium He nuclei. Starting from $\sim$1019 eV, the composition becomes heavier, i.e. the fraction of such nuclei as CNO and Fe iron nuclei increases.

III Conclusion

The Yakutsk array is operating continuously for more than 45 years, detecting electrons, muons, Cherenkov light, and radio emission. At the same time, more than 5$\cdot$106 EAS events were recorded in the energy region above 1015 eV. Using a large dataset of air showers with energies of 1016-5.7$\cdot$1019 eV, registered for the period 1974-2014, the longitudinal development of showers was reconstructed and the $X_{max}$ dependence was found. The results showed that elongation rate per decade of energy has nonlinear advancement and takes the values 48$\pm$6, 78$\pm$5, 63$\pm$6, 50$\pm$7 g / cm2 for energy ranges, 1016-1017, 1017-1018, 1018-1019, 1019-1020 eV. As can be seen (Fig. 7), there are kinks at the energies $\sim$1017 eV and $\sim$5$\cdot$1018 eV, i.e. on the ”second knee” area and the beginning of irregularities in the spectrum like ”dip-bump”. As shown by the data of the Yakutsk array (Fig. 8), this is due to a change in the cosmic ray mass composition. At low energies up to 1017 eV, there are noticeably more nuclei with an atomic mass of 4-56, at energies of 1017-1018 eV the proportion of protons increases and reaches a maximum, accounting for 60-80$\%$. In the region of energy greater than 1019 eV, CR predominantly consist of helium nuclei He, CNO, and heavier elements.

Acknowledgement

The reported study was funded by RFBR according to the research project 16-29-13019.