Theoretical aspects of heavy-ﬂavour production at ultra-high cosmic ray energies

. The main theoretical aspects of heavy-ﬂavour production at ultra-high cosmic ray energies are reviewed, with particular emphasis in the new dynamical effects which are expected to be present in the kinematical range probed by the IceCube and Pierre Auger Observatories. The gluon saturation effects for heavy quark production and the contribution of double parton scattering processes are analysed. Finally, the intrinsic heavy quark hypothesis is presented and some of its phenomenological implications at high energies are discussed.


Introduction
Ultra high energy cosmic rays (UHECRs) remains a puzzle in physics. Although the existence of UHECRs with energies above 10 11 GeV is now a well-established fact, the theoretical understanding of its origin and propagation is a subject of strong interest and intense discussion [1]. In particular, the recent detection of ultra-high energy (UHE) neutrinos by the IceCube Neutrino Observatory [2] starts a new era in the neutrino physics. These and forthcoming data from IceCube and Pierre Auger Observatories may shed light on the Standard Model as well as reveal aspects of new physics [3]. An important subject in cosmic ray physics is the flux of prompt leptons at the Earth which is related to primary interactions at energies that can by far exceed the highest available in accelerators. This flux is directly associated with charmed particle production and its decays, its estimation being strongly dependent on the model used to calculate the charm production cross section and energy spectra [4]. The quantification of the prompt leptons fluxes is essential, for instance, for neutrino physics, since the flux of prompt neutrinos provides the main background of the muon neutrinos from extra-galactic neutrino sources, which are being studied by the IceCube Observatory.
In this contribution I will review some of the main theoretical aspects of heavy-flavour production at ultrahigh cosmic ray energies, with particular emphasis in the new dynamical effects which are expected to be present in the kinematical range probed by the IceCube and Pierre Auger Observatories. Heavy quark production at high energies has been usually described considering the collinear factorization [5], where all partons involved are assumed to be on mass shell, carrying only longitudinal momenta, and their transverse momenta are neglected in the QCD matrix elements (see Fig. 1). The cross sections of the QCD subprocess are usually calculated in a e-mail: barros@ufpel.edu.br leading order (LO), as well as in next-to-leading order (NLO) (see, e.g. Refs. [6][7][8][9]). In particular, the cross sections involving incoming hadrons are given, at all orders, by the convolution of intrinsically non-perturbative (but universal) quantities -the parton densities -with perturbatively calculable hard matrix elements, which are process dependent. The conventional gluon distribution xg(x, µ 2 ), which drives the behaviour of the observables at high energies, corresponds to the density of gluons in the proton having a longitudinal momentum fraction x at the factorization scale µ. This distribution satisfies the DGLAP evolution [10] in µ 2 and does not contain information about the transverse momentum k T of the gluon. On the other hand, in the large energy (small-x) limit, the characteristic scale µ of the hard subprocess of parton scattering is much less than √ s, but greater than the QC D parameter. In this limit, the effects of the finite transverse momenta of the incoming partons become important, and the factorization must be generalized [11][12][13], implying that the cross sections are now k T -factorized into an off-shell partonic cross section and a k T -unintegrated parton density function F(x, k T ), characterizing the k T -factorization approach (see Fig. 2). The function F is obtained as a solution of the evolution equation associated to the dynamics that governs QCD at high energies. Several authors [14][15][16][17][18] have considered the k T -factorization approach in order to analyse some observables and they have obtained a better description of these quantities than the one obtained with the collinear approach. However, the current situation is still not satisfactory, due to the large uncertainty associated to the lack of a complete knowledge of the unintegrated gluon distribution.
The UHECR interactions probe the theory of the strong interactions in a new kinematical range characterized by a center of mass energy of approximately 500 TeV, which is more than one order of magnitude larger than that of the Large Hadron Collider at CERN. In this energy regime, perturbative Quantum Chromodynamics (pQCD) predicts that the small-x gluons in a hadron wave function should EPJ Web of Conferences Collinear factorization: s ≈ µ 2 Λ 2 QCD Figure 1. Heavy quark production at high energies in the collinear factorization formalism. Figure 2. Heavy quark production at high energies in the k T factorization formalism.
form a Colour Glass Condensate (CGC) [19][20][21], which is characterized by the limitation on the maximum phasespace parton density that can be reached in the hadron wave function (parton saturation), with the transition being specified by a typical scale, which is energy dependent and is called saturation scale Q s [22]. In particular, at high energies the nonlinear QCD effects associated to the gluon saturation are expected to contribute significantly and this leads to the breakdown of the twist expansion and of the factorization schemes present in the collinear and k T -factorization approaches [23]. An alternative is to estimate the heavy quark production using the color dipole formalism [24,25] which allows to naturally to include the gluon saturation effects. At high energies color dipoles with a defined transverse separation are eigenstates of the interaction. The main quantity in this formalism is the dipole-target cross section, which is universal and determined by QCD dynamics at high energies. In particular, it provides a unified description of inclusive and diffractive observables in lepton -hadron processes as well as for in Drell-Yan, prompt photon and heavy quark production in hadron-hadron collisions. Furthermore, an important advantage of this formalism is that it is very simple to include nuclear effects. Another possible new effect which can contribute for ultra high energies hadron -hadron interactions is the heavy quark production in multiple gluon -gluon interactions. The high density of gluons in the initial state of hadron -hadron collisions at high energies implies that the probability of these multiple gluongluon interactions within one collision increases. In particular, the probability of having two or more hard interactions in a collision is not significantly suppressed with respect to the single interaction probability. Recent theoretical [26][27][28] and experimental [29] results demonstrate that the contribution of double parton scattering (DPS) process for the heavy quark production is not negligible already for LHC energies. As demonstrated in Ref. [30] its contribution is very important for ultra high energy cosmic ray interactions.
Finally, the study of heavy quark production at ultrahigh energies and forward rapidities allows to test the hypothesis of intrinsic heavy quarks in the hadrons [31]. In the last years several groups have proposed different schemes to determine these distributions considering that the heavy quark component in parton distribution functions (PDF) can be perturbatively generated by gluon splitting (see e.g. [32]). This component is usually denoted extrinsic heavy quark component. Moreover, the possibility of an intrinsic component has been studied in detail and included in the recent versions of the CTEQ parametrization [33,34]. The hypothesis of intrinsic heavy quarks (IHQ) is a natural consequence of the quantum fluctuations inherent to Quantum Chromodynamics (QCD) and amounts to assuming the existence of a QQ (Q = c, b, t) as a nonperturbative component in the hadron wave function. One of the most striking properties of an IHQ state, such as |uud QQ , is that the heavy constituents tend to carry the largest fraction of the momentum of the hadron. Consequently, in contrast to heavy quarks produced through usual perturbative QCD, which emerge with small longitudinal momentum, the intrinsic charm component gives rise to heavy mesons with large fractional momenta relative to the beam particles. Therefore, the existence of an intrinsic component modifies, for instance, the x F and rapidity distribution of charmed particles (see, e.g. [35,36]).
In what follows we present a more detailed discussion about the formalisms and possible new effects in the heavy quark production at ultra-high energies summarized in this Section. In particular, in the next section we present a brief review of the heavy quark production in the collinear, k T -factorization and dipole approaches. In Sect. 3 we discuss the description of the gluon saturation effects and its main phenomenological implications for the heavy quark production in ultra-high cosmic ray interactions. Moreover, in Sect. 4 we discuss the contribution of the double parton scattering processes and in Sect. 5 possible implications of an intrinsic heavy quark component in the hadron wave function for the heavy production at high energies. Finally, in Sect. 6 we summarise the main conclusions.

Heavy quark production
The understanding of the dynamics of heavy quark production has been improved considerably over the last years, in particular by the large amount of precise data from experiments at HERA, Tevatron and LHC. As discussed in the Introduction, there are currently several approaches to treat the heavy quark production. However, the collinear framework [5] is by far the most developed approach. It is based on the factorization theorem, which states that the cross section can be factorized into three different parts (see Fig. 2): the non-perturbative initial Figure 3. The collinear factorization for heavy quark production in hadronic collisions. condition, described by the parton distribution functions f i (x i , µ 2 F ), the hard process itself, perturbative calculable and described by the partonic cross sectionσ i j for the subprocess i j → QQ X, and the subsequent hadronization (also non-perturbative). Thus the total cross section for the heavy quark production is given in the collinear framework by (for an introductory review see Ref. [37]) where µ R and µ F are the renormalization and factorization scales. Predictions for hadronic cross sections crucially depend on the knowledge of bothσ i j and f i/ h 1 (x 1 , µ 2 F ). The short distance cross sectionsσ can be calculated up to a given perturbative order. At leading order (LO), the only processes which can lead to heavy quark production are quark -antiquark annihilation and gluon fusion, as illustrated in Fig. 4. The resulting partonic cross sections are determined by the mass of the heavy quark and by the strong coupling constant. The energy behaviour of the total production cross section for heavy quarks is determined by the gluon distribution, since the gg-channel is dominant at high energies. Moreover, it is finite at p T → 0, which is directly associated to the fact that the minimum virtuality of the t-channel propagator is m 2 Q . It is important to emphasize that the LO calculation gives a rough estimate of the cross section.
In order to obtain a better estimate of the cross section and a rough estimate of its theoretical uncertainty, we need to include the next-to-leading order (NLO) corrections. In  of the produced heavy hadron (hadron containing a heavy quark), the momentum Q 2 in deep inelastic scattering, the mass of the produced heavy hadron and the center-ofmass energy √ s. Different schemes have been developed to obtain predictions from pQCD depending on the specific kinematical region and the relative importance of the relevant scales. In the case of relatively small transverse momentum ( p T ≤ m Q ), the fixed-flavour number scheme (FFNS) [6][7][8][9] is usually applied. Here one assumes that the light quarks and the gluon are the only active flavours within the colliding hadrons and the heavy quarks appears only in the final state. In the complementary kinematical region ( p T m Q ), calculations are usually based on the zero-mass variable-flavour-number scheme (ZM-VFNS) [38]. In the ZM-VFNS, the zero mass parton approximation is applied for the heavy quarks. The heavy quark acts also as an incoming parton with its own parton distribution function. The predictions obtained in this scheme are expected to be reliable only in the region of large p T since all terms of the order m 2 Q / p 2 T are neglected in the hard cross section. In the last years two different approaches, that combines the FFN and ZM-VFN schemes, has been developed: the so-called general-massvariable-flavour-number scheme (GM-VFNS) [39] and the fixed-order plus next-to-leading-logarithm (FONLL) theory [40]. In particular, in the FONLL approach, the cross section for inclusive heavy meson production is calculated from a superposition of the FFN and the ZM-VFN approaches with a p 2 T dependent weight function. This function is fixed in such a way that for p 2 EPJ Web of Conferences cross section is fixed order and for large p 2 T approaches the ZM-VFNS result. In addition, the heavy quark contains perturbative pieces at the starting scale, in order that the ZM-VNS result matches the FFNS in next-to-leading order. Recent analysis [41] show the these approaches are able to describe the current experimental data if the theoretical uncertainties, associated to the choice of the mass, renormalization and factorization scales, fragmentation functions and parton distributions, are taken into account.
The measurement of heavy quark production at high energies in the low transverse momentum region probes the parton distribution functions at small values of the parton fractional momentum x and squared momentum transfer Q 2 . For illustration, using a simplified 2 → 2 kinematics at leading order, we obtains that at LHC energies c quarks with p T ≈ 2 GeV and rapidity y ≈ 0 probe the parton distribution functions at x ≈ 7 × 10 −4 , where the gluon component is dominant. At forward rapidities, smaller values of x will be reached in one of the projectiles. Moreover, in ultra-high energy cosmic ray interactions we obtain that x ≤ 10 −6 for a primary proton with energy larger than 1 PeV. These values are far smaller than anything accessible in the current accelerators, which challenges QCD in a new dynamical region where potentially large higher-order corrections need to be resummed. In particular, in this kinematical range, the characteristic scale µ of the hard subprocess of parton scattering is much less than √ s, but greater than the QC D parameter. In this limit, the effects of the finite transverse momenta of the incoming partons become important, and the collinear factorization must be generalized. In Refs. [11][12][13] the authors have proposed a generalized factorization, denoted k T -factorization, which systematically resums powers of α s ln(s/µ 2 ) present at high energies. In this formalism the Feynman diagrams are calculated taking account of the virtualities and of all possible polarizations of the incident partons. The cross section for the heavy quark production is expressed as a convolution of off-shell matrix elements and k T -unintegrated gluon distribution F(x i , k 2 i T ), which depend on the momentum fraction x i and the momentum transverse k i T of the gluon taking part in the hard scattering. The cross section is given by where y i are the quark rapidities, k i T are the gluon transverse momenta, |M of f g * g * →QQ | 2 is the off-shell matrix element for the production of a heavy quark pair in a gluon -gluon collision and the unintegrated gluon distribution F(x i , k 2 i T ) is solution of the evolution equation associated to the dynamics that governs QCD at high energies (e.g. BFKL, CCFM, BK, ...). As a significant part of the NLO corrections for the heavy quark cross Figure 6. Heavy quark hadroproduction in the color dipole formalism.
section in the collinear factorization are associated with the contribution of nonzero transverse momenta of the initial partons, such corrections have already been included at leading order in the k T -factorization approach. In Ref. [18] the open charm production at the LHC have been studied in detail and obtained that k T -factorization predictions are consistent with recent LHC measurements as well as with the recent FONLL and GM-VNFS predictions if the Kimber -Martin -Ryskin model for the unintegrated gluon distribution is used as input in the calculations. However, as already pointed out in Ref. [18], its predictions are strongly sensitive to the choice of F, which is directly related to the treatment of the QCD dynamics at high energies. The heavy quark hadroproduction can also be described in terms of the color dipole cross section using the color dipole approach [24,25]. In this approach the total heavy quark production cross section is given by [25] where y is the rapidity of the pair, µ F ∼ m Q is the factorization scale, x 1 G(x 1 , µ 2 F ) is the projectile gluon density at scale µ F and longitudinal momentum fraction x 1 and the partonic cross section σ (G N → QQ X) is given by where ρ is the dipole size, α (ᾱ = 1 − α) is the fractional momentum of quark (antiquark) and G→QQ is the pQCD calculated distribution amplitude, which describes the dependence of the |QQ Fock component on transverse separation and fractional momentum. Moreover, σ QQG is the cross section for scattering a color neutral quarkantiquark-gluon system on the proton and is directly related with the dipole cross section as follows The basic idea on this approach is that at high energies a gluon G from the projectile hadron can develop a fluctuation which contains a heavy quark pair (QQ). Interaction with the color field of the target then may release these heavy quarks. The Eq. (3) can be easily interpreted in the target rest frame, where heavy quark production looks like pair creation in the target color field. For a short time, a gluon G from the projectile hadron can develop a fluctuation which contains a heavy quark pair (QQ). The interaction with the color field of the target may then release these heavy quarks. This mechanism corresponds to the gluon-gluon fusion mechanism of heavy quark production in the leading order (LO) parton model. The dipole formulation is therefore applicable only at low x, where the gluon density of the target is much larger than all quark densities. This condition is fulfilled for charm and bottom production at CERN-LHC and ultrahigh cosmic ray interactions.

Gluon saturation effects
The description of QCD dynamics in the high energy limit is still a subject of intense debate (for a recent review see e.g. Ref. [42]). Theoretically, at high energies (small Bjorken-x) one expects the transition of the regime described by the linear dynamics (denoted Low Density QCD in Fig. 7), where only the parton emissions are considered, to a new regime where a high gluon density is present and the physical process of recombination of gluons becomes important in the parton cascade, with the evolution being given by a nonlinear evolution equation. This high density QCD regime is characterized by the limitation on the maximum phase-space parton density that can be reached in the hadron wave function (gluon saturation), with the transition between the low and high density regimes being specified by the saturation scale Q s [22] (see Fig. 7). Currently, the state of art formalism to describe the gluon saturation effects is the Colour Glass Condensate. In this formalism the nonlinear and quantum effects in the hadron wave function are encoded in the imaginary part of the forward amplitude for the scattering between a small dipole (a colorless quark-antiquark pair) and a dense hadron target, at a given rapidity interval Y = ln(1/x), which is represented by N (b, ρ, x). This quantity is directly related to the dipole -hadron cross section, σ di p , which encodes all the information about the hadronic scattering, as follows N (b, ρ, x).
The dipole has transverse size given by the vector ρ = x − y, where x and y are the transverse vectors for the quark and antiquark, respectively, and impact parameter b = (x + y)/2. At high energies the evolution with the rapidity Y of N (ρ, b, Y ) is given by the infinite hierarchy of equations, the so called Balitsky-JIMWLK equations [19,20], which reduces in the mean field approximation to the Balitsky-Kovchegov (BK) equation [19,21]. It is useful to assume the translational invariance approximation, which regards hadron homogeneity in the transverse plane. It implies that the dipole-proton cross section reads where the constant σ 0 , which results from the b integration, sets the normalization. Moreover, the amplitude becomes independent of the impact parameter b and depends only on the dipole size ρ = |ρ|, i.e. N Y (ρ) = N Y (ρ). Although a complete analytical solution of the BK equation is still lacking, its main properties are known: (a) for the interaction of a small dipole (ρ 1/Q s ), N (ρ) ≈ ρ 2 , implying that this system is weakly interacting; (b) for a large dipole (ρ 1/Q s ), the system is strongly absorbed and therefore , is the so called saturation scale. This property is associated to the large density of saturated gluons in the hadron wave function. A numerical solution which considers running coupling corrections to the kernel of BK equation is available in the literature [43]. The calculations using this numerical solution will be denoted as "rcBK" hereafter. Currently, the rcBK model is the most sophisticated saturation model available in the literature. In what follows we will also present the predictions of the phenomenological saturation model proposed in Ref. [44], denoted GBW hereafter, which encodes the main properties of the saturation approaches, with the scattering amplitude parametrized as follows where the saturation scale is given by Q 2 s = Q 2 0 (x 0 /x) λ , x 0 is the value of the Bjorken x in the beginning of the evolution and λ is the saturation exponent. The parameters x 0 and λ are obtained by fitting the ep HERA data. We assume the values obtained in Ref. [45], where the GBW model was updated in order to describe more recent data. Moreover, in order to estimate the magnitude of the gluon saturation effects we will also present the predictions obtained disregarding the nonlinear effects in the dipoleproton scattering amplitude, which implies N (ρ, Y ) = ρ 2 Q 2 s (Y )/4 (denoted GBW linear hereafter). The heavy quark production considering saturation effects was studied in detail in Refs. [28,46]. In [46] we have predicted the energy dependence of the charm and bottom pair production and compared with data 07001-p.5 EPJ Web of Conferences points from UA2, PHENIX and from Cosmic Rays. All these data can be quite well described using the color dipole formalism and an adequate choice of the heavy quark mass. Considering the charm data released by the ALICE Collaboration [47], it was possible to compare in Ref. [28] the color dipole formalism for heavy quark production in hadronic collisions with experimental data at high energies, which is the kinematical range where it is theoretically justified. In Fig. 8 (upper panel) we compare the rcBK, GBW and GBW Linear predictions (upper curves) with the ALICE [47] and PHENIX [48] data considering m c = 1.5 GeV. We can see that the different models for the dipole-target cross sections are able to describe the experimental data. The rcBK and GBW predictions are almost identical in the whole range of energy. When the GBW Linear model is used as input in the calculations, larger values for the charm cross section are predicted. In contrast, for bottom production [see Fig. 8 (lower panel)], the GBW and GBW Linear predictions are identical in the considered energy range and the rcBK one predicts larger values of the cross section. The distinct behaviour observed for charm and bottom production is directly associated to the fact that in the color dipole formalism the contribution of the nonlinear effects is determined by the integrand of the pair separation (ρ) integral, i.e. by the product of the wave function squared and the dipole-target cross section. This integrand has a peak at ρ ≈ 1/m Q . Consequently, for bottom production, the integral is dominated by very small pair separations, probing the linear regime of the dipoletarget cross section. The rcBK prediction is larger than the GBW and GBW Linear ones because its linear regime is associated to the BFKL dynamics, which implies a steeper energy dependence. In contrast, for charm production, we probe larger values of the pair separation, where saturation effects cannot be disregarded. The difference between the rcBK and the GBW predictions is associated to the delayed saturation predicted by the rcBK equation.
In order to estimate the heavy quark contribution for the total cross sections at ultra high energy cosmic ray interactions, in Ref. [30] we have estimated the energy dependence of the ratio where i characterizes the primary cosmic ray, which can be either a photon, neutrino or a proton. Moreover, σ i p H Q is the total heavy quark cross section for the production of a given flavour. In what follows we review our results for the heavy quark production in proton -proton collisions at ultra high energies (for more details see [30]). In Fig. 9(a) we present the predictions for the ratio R H Q [ pp] as a function of the energy E p of the primary proton, obtained assuming that the total pp cross section is given by the parametrization proposed by the COMPETE Collaboration [49], as described in [50]. The results demonstrate that the charm contribution for the total cross sections is larger than 10% at E p ≥ 10 8 GeV and that at very high energies ≈ 10 11 GeV it becomes of the order of the total pp cross section. In order to estimate the magnitude of the gluon saturation effects in the behaviour of the ratio, in Fig. 9(b) we present a comparison between the GBW predictions and those obtained disregarding the nonlinear effects in the dipole -proton scattering amplitude, i.e. with GBW linear model. For bottom production, these effects can be disregarded. In contrast, for charm production, the gluon saturation effects diminish the ratio by ≈ 35 % for E p = 10 12 GeV. In particular, the ratio becomes of order of one at smaller energies if these effects are disregarded. Clearly, these results indicate that new dynamical effects, beyond those included in the dipole factorization (Eq. (3)) and in the QCD dynamics through the gluon saturation effects, should be considered in the kinematical regime probed in UHECR interactions in order to obtain reliable predictions for the charm production. In particular, these results can be interpreted as a signal of the breakdown of the factorization, which are expected at very high energies. In contrast to the standard calculations of the heavy quark contributions for the total cross sections, which are based on the collinear factorization at leading twist and the solution of the linear DGLAP equations, the cross sections obtained in the color dipole approach resums higher-twist corrections beyond the traditional factorization schemes and allow us to include gluon saturation effects in a straightforward and unified way. However, it is important to emphasize that still is not clear if the dipole factorization should not be generalized at larger than the current collider 07001-p.6 energies. Results obtained in Ref. [23] for the quark pair or gluon production in the interaction of two ultra dense systems indicate the breakdown of the factorization at high energies. As the proton gluon density for ultra high cosmic ray energies is expected to be very large, a similar scenario can be present in pp cosmic ray interactions. However, if and for what energy this breakdown occurs are still open questions which deserve more detailed studies. Finally, it is important to emphasize that we have verified that the main conclusions of Ref. [30] are not modified if another phenomenological saturation model or the solution of the BK equation, as given in Ref. [43], are used as input in the calculations.

Double parton scattering
The calculation of the heavy quark cross section in the standard framework assumes that only one hard interaction occurs per collision. This mechanism is called singleparton scattering (SPS), since the Feynman diagram contains one gluon from the hadron target and one gluon from the hadron projectile [see Fig. 10 (upper panel)]. In general, higher order corrections do not change significantly the observables, since the contributions are suppressed by powers of α s . For example, the Q 1Q1 Q 2Q2 (Q i = c or b) production in SPS processes [see Fig. 10 (lower panel)] is suppressed by a factor proportional to α 4 s (for explicit calculations of two heavy quark pairs production in SPS processes see, e. g., Refs. [51,52]). The basic idea, which justifies the SPS approach, is that the probability of a hard interaction in a collision is very small, which makes the probability of having two or more hard interactions in a collision highly suppressed with respect to the single interaction probability. Such an assumption is reasonable in the kinematical regime in which the flux of incoming partons is not very high. However, as already pointed out in Ref. [53], at LHC energies the hadronic cross section appears to be three orders of magnitude higher than the cross section of the partonic subprocess. In this condition, there is a high probability of scattering of more than one pair of partons in the same hadron -hadron collision. This expectation has been confirmed by the LHCb Collaboration [29], which has observed the production of J/ mesons accompanied by open charm and pairs of open charm hadrons in pp collisions at √ s = 7 TeV and verified that the SPS predictions are significantly smaller than the observed cross sections. In [27] it was shown that the DPS contribution (see Fig. 11) to differential cross sections of open charm and charmed meson production is rather significant, being of the same order of magnitude of the SPS differential cross sections. Remarkably, the sum of the SPS with the DPS contribution almost describes the experimental data from ATLAS, LHCb and ALICE on transverse momentum distribution of charmed meson production [27]. The contribution of cccc production in SPS processes to differential cross sections was shown to be negligible when compared with the production in DPS, reaching at most ≈ 10% of the magnitude of this last one, but in most cases being two orders of magnitude smaller [27].
Following the same factorization approximations assumed for processes with a single hard scattering, it is possible to derive the DPS contribution for the heavy quark cross section considering two independent hard parton 07001-p.7  Figure 11. The Q 1Q1 Q 2Q2 pair production in the double parton scattering (DPS) process. sub-processes. It is given by (see, e.g. Ref. [54])

EPJ Web of Conferences
where we assume that the quark-induced sub-processes can be disregarded at high energies, 1 , µ 2 2 ) are the two gluon parton distribution functions which depend on the longitudinal momentum fractions x 1 and x 2 , and on the transverse position b 1 and b 2 of the two gluons undergoing the hard processes at the scales µ 2 1 and µ 2 2 . The functionsσ are the parton level subprocesses cross sections and b is the impact parameter vector connecting the centers of the colliding protons in the transverse plane. Moreover, m/2 is a combinatorial factor which accounts for indistinguishable and distinguishable final states. For Q 1 = Q 2 one has m = 1, while m = 2 for Q 1 = Q 2 . In a rigorous calculation several kinds of correlations between the two gluons in the double gluon distribution function gg h 1 should be taken into account (see, e.g., Ref. [55]). However, a precise estimation of the magnitude of the correlations is very difficult, and in practical calculations most of the correlations are disregarded. Moreover, some of these correlations, as the color correlation and the parton-exchange interference, are Sudakov suppressed at high energies, and can be neglected in this kinematical regime [55]. It is common in the literature to assume that the longitudinal and transverse components of the double parton distributions can be decomposed and that the longitudinal components can be expressed in terms of the product of two independent single parton distributions. The proof of these assumptions in the general case is still an open question (see, e.g. Refs. [54,56]). In the particular case of heavy quark production, in Ref. [26] the authors compared the results of this simple factorized Ansatz with those obtained using double parton distributions with QCD evolution and verified that the predictions are similar if we take into account all uncertainties present in the calculations as, for instance, those associated to the choice of the factorization and renormalization scales.
In Ref. [28] we have estimated the gluon saturation effects for the double parton scattering (DPS) process of the heavy quark production. There we have assumed that the DPS contribution to the heavy quark cross section can be expressed in a simple generic form given by where Q = c or b, σ S P S is the cross section for the heavy quark production considering the mechanism of single-parton scattering (SPS), calculated using Eq. (3), and σ e f f is a normalization cross section representing the effective transverse overlap of partonic interactions that produce the DPS process, which we assumed as being σ e f f = 15 mb. Moreover, 1/2 is a combinatorial factor which accounts for indistinguishable final states. The Eq. (8) expresses the DPS cross section as the product of two individual SPS cross sections assuming that the two SPS sub-processes are uncorrelated and do not interfere (see [28] for more details). In Fig. 7 we present the predictions obtained in Ref. [28] for the energy dependence of the DPS contribution considering different models for the dipole -proton cross section. For charm production the double parton scattering contribution becomes comparable with the single parton scattering one at LHC energies and this result remains valid when one introduces saturation effects in the calculations. This result motivated the study performed in Ref. [30], where the magnitude of this new mechanism for ultra high energy cosmic ray interactions was estimated. In Fig. 12(a) we present the predictions obtained in [30] for the dependence of the ratio R D P S ≡ σ D P S /σ S P S on the energy of the primary proton. The DPS mechanism for charm production becomes of the order of the SPS one at E p = 4 · 10 8 GeV and dominates at higher energies. For bottom production, the ratio only becomes 1 for E p ≈ 10 12 GeV. The magnitude of the gluon saturation effects is estimated in Fig. 12(b). As expected from previous results for the heavy quark production in the SPS process, these effects are negligible for the bottom production and diminishes the ratio R D P S by a factor ≈ 2 for E p = 10 12 GeV. The results presented in Ref. [30] indicate that the DPS mechanism cannot be disregarded at cosmic ray energies and should be included in the airshower simulators in order to obtain reliable predictions for the production of heavy hadrons and consequent flux of prompt leptons at the Earth (for recent studies see Refs. [57,58]). In particular, we expect that this new mechanism could imply an enhancement of this flux at high energies in comparison with previous studies [4], which could have implications at the energies probed at IceCube. This subject deserves more detailed studies in the future.

Intrinsic heavy quark contribution
As discussed in the Introduction, the description of heavy quarks in the hadron wave functions is still an open question. In general it is assumed that the heavy quark content of the nucleon comes mostly from the DGLAP evolution of the initial gluon distribution, as represented in Fig. 13 (left panel). This process is well understood in perturbative QCD and is denoted extrinsic component. However, as pointed out long ago in Ref. [31] (see also 07001-p.8 Ref. [59]) there may be another, higher twist, component which is typical of the hadron where it is created. This component is called intrinsic heavy quark (IHQ). A comprehensive review of the main characteristics of the IHQ models can be found in [60]. A well known model was proposed in Ref. [61]. In this model, the creation of the Q Q pair was studied in detail. It was assumed that the nucleon light cone wave function has higher Fock states, the first one being |qqq Q Q > [see Fig. 13 (right panel)]. The probability of finding the nucleon in this configuration was given by the inverse of the squared invariant mass of the system. Because of the heavy charm mass, this probability as a function of the quark fractional momentum, P(x), is very hard, as compared to the one obtained through the DGLAP evolution. We denote this model by BHPS hereafter. A more dynamical perspective for the intrinsic charm (IC) component is given by the meson cloud model. In this model, the nucleon fluctuates into an intermediate state composed by a charmed baryon plus a charmed meson [62]. The charm is always confined in one hadron and carries the largest part of its momentum. In the hadronic description we can use effective lagrangians to compute the charm splitting functions, which turn out to favor harder charm quarks than the DGLAP ones. The main  Figure 14. Intrinsic charm momentum distribution in the proton with different IC models taken from [33]. difference between the BHPS and meson cloud models is that the latter predicts that the charm and anticharm distributions are different, which generates a difference in the momentum distributions of D + and D − [63].
In spite of its successes, the idea of the intrinsic charm was left aside for some years. In Ref. [33] it was considered as an ingredient in the global fit of DIS data performed by the CTEQ collaboration (for a more recent analysis see Ref. [34]). In [33] the CTEQ group determined the shape and normalization of the IC distribution in the same way as they do for other parton species. In fact they find several IC distributions which are compatible with the world data. Apart from the already mentioned BHPS and meson cloud models, the CTEQ group has tested another model of intrinsic charm, called sea-like IC. It consists basically in assuming that at a very low resolution (before the DGLAP evolution) there is already some charm in the nucleon, which has a typical sea quark momentum distribution ( 1/ √ x) with normalization to be fixed by fitting DIS data. The different parametrizations are presented in Fig. 14. For comparison we also present the no-IC distribution, where the charm content of the nucleon sea comes from the DGLAP evolution (exclusive charm component). We can see that the intrinsic charm distribution is a factor ten (sea like) to twenty (BHPS) larger than the no-IC distribution (the long dashed line) and the peaks of the IC distributions lie in the large x (> 0.1) region. This behaviour is so peculiar that it gives us hope to observe IC experimentally.
Although the existence of IHQ fluctuations in the proton has substantial and growing experimental and theoretical support, as revised in Ref. [60], more definite conclusions are still not possible. One of the shortcomings of previous studies is the use of the QCD factorization at EPJ Web of Conferences large Feynman x F , where it is expected to breakdown (see, e.g., [64]). Basically, at large x F the kinematics is very asymmetric, with the hadrons in the final state emerging from collisions of projectile partons with large light cone momentum fractions (x p → 1) with target partons carrying a very small momentum fraction (x t 1). In this case, we have a scattering of a dilute projectile on a dense target, where the small-x effects coming from CGC physics are expected to show up and the usual factorization formalism is expected to breakdown [23]. This kinematical region has been explored in dAu collisions at forward rapidities at RHIC and the suppression of high p T hadron yields, anticipated on the basis of CGC ideas, was observed. The satisfactory description of these experimental data [65,66] gives us a strong indication that the appropriate framework to calculate particle production at large x F (large rapidities) is one based on CGC physics. Having this in mind, in Ref. [36] we have adapted the hybrid formalism proposed in [67] (to calculate the production of high-x quarks in d Au collisions) to D meson production in pp ( p(d)A) collisions at forward rapidities, including intrinsic charm quarks in the projectile wave function. In this approach the target is treated as a dense system while the projectile proton is taken to be a dilute system in the spirit of standard QCD. In Fig. 15 we present a representation of the scattering process in the hybrid formalism. We assume that the charm quarks are already in the projectile and their densities are those given by the CTEQ parametrization [33]. Moreover, we describe the interaction of these charm quarks in terms of the dipole scattering amplitude proposed in [65] to describe the RHIC data. Finally, we assume that the charm quarks fragment into D mesons according to the fragmentation functions given in Ref. [68]. Our goal in Ref. [36] was to estimate the D-meson spectrum at forward rapidities at RHIC and LHC energies assuming the presence (or absence) of an intrinsic charm component and evaluate how much this quantity is modified.
In Fig. 16 we present the p T -dependence for the ratio between the spectra predicted by the different IC models and the prediction without the IC component (no IC). We consider the D meson production in proton-proton collisions at RHIC energy ( √ s = 200 GeV) and fixed rapidity (y = 4). One can see that the presence of intrinsic charm implies an enhancement of the D meson production in the whole p T range. The behaviour of the different models can be easily understood if we remember that:  From this expression we see that, in the p T range shown in the figure, at RHIC energies and y = 4, the typical values of x p are in the range 0.4 x p 0.95. In this range, all IC models predict a larger charm component than the no-IC model as shown in Fig. 14. This also implies larger values for the cross section calculated with BHPS and meson cloud IC models. The cross section calculated with the sea like IC model is smaller, since it predicts, already at x ≈ 0.6, a very small charm distribution, quite similar to the no-IC one. On the other hand, the BHPS and meson cloud models predict a very large IC distribution in this kinematical range, which implies larger values for the cross section. As expected from Fig. 14, the BHPS prediction for the spectrum is larger than the meson cloud one.

Summary
Heavy quark production in hard collisions of hadrons, leptons, and photons has been considered as a clean test of perturbative QCD. This process provides not only many tests of perturbative QCD, but also some of the most important backgrounds to new physics processes, which have motivated comprehensive phenomenological studies carried out at DESY-HERA, Tevatron and LHC. With the advent of the Pierre Auger and IceCube observatories, the extension of these studies for ultra high energy cosmic ray interactions becomes fundamental. In this contribution we have reviewed the formalisms to treat heavy quark production and new effects which are expected to be present at these energies. We have discussed gluon saturation effects, associated to the physical process of parton recombination. The results indicate that the contribution of heavy quarks at ultra-high cosmic ray interactions is non-negligible and that new dynamical mechanisms, as for instance, double parton scattering processes, should be included in order to obtain reliable predictions.
I am very grateful to R. Engel for the invitation to attend to the International Symposium on Very High Energy Cosmic Ray Interactions 2014. Moreover, I would like to congratulate the organizers for the stimulating and interesting meeting. I also wish