Intrinsic charm in the nucleon and charm production at large rapidities in collinear, hybrid and $k_T$-factorization approaches

We discuss the role of intrinsic charm (IC) in the nucleon for forward production of $c$-quark (or $\bar c$-antiquark) in proton-proton collisions for low and high energies. The calculations are performed in collinear-factorization approach with on-shell partons, $k_T$-factorization approach with off-shell partons as well as in a hybrid approach using collinear charm distributions and unintegrated (transverse momentum dependent) gluon distributions. For the collinear-factorization approach we use matrix elements for both massless and massive charm quarks/antiquarks. The distributions in rapidity and transverse momentum of charm quark/antiquark are shown for a few different models of IC. Forward charm production is dominated by $gc$-fusion processes. The IC contribution dominates over the standard pQCD (extrinsic) $gg$-fusion mechanism of $c\bar c$-pair production at large rapidities or Feynman-$x_F$. We perform similar calculations within leading-order and next-to-leading order $k_T$-factorization approach. The $k_T$-factorization approach leads to much larger cross sections than the LO collinear approach. At high energies and large rapidities of $c$-quark or $\bar c$-antiquark one tests gluon distributions at extremely small $x$. The IC contribution has important consequences for high-energy neutrino production in the Ice-Cube experiment and can be, to some extent, tested at the LHC by the SHIP and FASER experiments by studies of the $\nu_{\tau}$ neutrino production.


I. INTRODUCTION
The text-book proton consists of uud valence quarks. This picture is by far too simplified. In fact there is strong evidence for internal strangeness and somewhat smaller for internal charm content of the nonperturbative proton. Different pictures of nonperturbative cc content were proposed in the past. A first example is relatively old BHPS model [1] which assumes uudcc 5-parton Fock configurations (see also Refs. [2,3]). Another picture proposed in the literature is a meson cloud model (MCM) [4][5][6][7][8], where the p →D 0 Λ c or (1.1) It is rather upper limit but this value depends somewhat on the model of charm content of a proton. In general, for sea-like models the probability can be slightly larger than for the BHPS one. In the sea-like case the charm is concentrated at lower values of x. Very recent lattice study of charm quark electromagnetic form factors suggested asymmetry of c-quark andc-antiquark distributions [9].
Recently there is a renewed interest in the intrinsic charm (IC) which is related to experiments being performed at the LHC [10][11][12][13]. The intrinsic charm is often included in global parton analyses of world experimental data [14][15][16].
The highly energetic neutrino experiments, such as IceCube, could put further constraints on the intrinsic charm [17][18][19]. Here, however, the IC contribution may compete with a concept of the subleading fragmentation [20]. Similarly, future LHC high and low energy forward experiments like FASER and SHIP could also be very helpful in this context (see e.g. Ref. [21] and Ref. [22], respectively). Also the LHCb experiment in its fixed-target mode could be sensitive to the contributions coming from intrinsic charm in a proton, especially in the case of open charm production [23], where some problems with a satisfactory theoretical description of the experimental data were reported (see also discussion in Ref. [46]).
In this paper we concentrate therefore on forward production of charm quarks/antiquarks. There were already some studies performed with color glass con-densate approach and compared to the dipole approach at forward directions [25,26]. In our approach we will use instead collinear, hybrid and k T -factorization approach. The latter two were not studied so far in the context of IC and forward production of charm.

II. MODELS OF INTRINSIC CHARM IN A NUCLEON
In the five-quark Fock component uudcc heavy quark/antiquark carries rather large fraction of the mother proton. In the BHPS model, after some approximations the probability to find c orc (the same for both) can be expressed via a simple formula: The normalization constant A depends on integrated probability for cc component and is 6 for 1 % probability. Please note that the quark mass is not explicit in this simplified formula.
In the meson cloud models c is in the baryon-like object andc in the meson-like object.
Then the probabilistic distribution can be obtained as 3) The f B and f M functions, the probability to find meson or baryon in proton, can be calculated from corresponding Lagrangians supplemented by a somewhat arbitrary and poorly known vertex form factors and can be found e.g. in Ref. [8]. In general, such an approach leads to c(x) =c(x).
In practice both models give rather similar distributions as will be shown in the following, so using one of them as an example is representative and sufficient. These are models of large-x components of IC. In principle, the IC may have also small-x component known under the name of sea-like, however, only simple ad hoc parametrizations were used in the literature.
There is another category of processes leading to sea-like IC (see Fig. 1 where an example of corresponding dynamical processes is shown). Using intrinsic glue in the nucleon (see e.g. Ref. [27]) one can generate intrinsic charm sea. The intrinsic gluon distribution x dy 1 y P qg x y g(y) , (2.5) where g is the intrinsic gluon distribution. With the model from Ref. [27] the c andc distributions are integrable, concentrated at x ∼ 0.1-0.2 and corresponding probability is 2.7 %. It should be less for massive quarks/antiquarks. In Fig. 2 we show x-distribution of the IC for the BHPS and for the sea-like model described above. In this calculation BHPS model with 1% probability was used for illustration.
In the GRV approach [28] the charm contribution is calculated fully radiatively as the convolution of gluon distribution with appropriate mass-dependent splitting function: The explicit formula for C c g,2 and a, including mass of quarks/antiquarks, can be found in Ref. [28]. In the following calculations we will use more modern gluon distributions.
In this paper we concentrate on large-x component and completely ignore the sea-like component(s). Charm can also be generated by evolution equations via g → cc transition (splitting). Often it was included in the evolution as a massless parton with zero as initial condition at the starting scale µ 2 ∼ m 2 c . In a dedicated fits, the intrinsic charm distribution is used as initial condition for DGLAP evolved charm distributions (see e.g. Ref. [29]). In the right panel Fig. 4 we show charm distribution in a proton without (dashed line) and with (solid line) the IC distribution taken as initial condition of the evolution.

A. The collinear approach
In the present study we discuss production of the final states with one charm quark or charm antiquark. In the collinear approach [30] the final state charm must be associated with at least one additional gluon or (light) quark. Typical leading-order mechanisms for charm production initiated by charm quark in a initial state are shown in Fig. 9. The diagrams correspond to the gc → gc (or gc → gc) subprocesses that are expected to be dominant at high energies, however, the qc → qc andqc →qc (or qc → qc and qc →qc) mechanisms with q = u, d, s are also possible and will be taken into account in the following numerical calculations. In the collinear approach the differential cross section for forward charm production within the gc → gc mechanism 1 can be calculated as where M gc→gc is the on-shell matrix element for gc → gc subprocesses and g(x 1 , µ 2 ) and c(x 2 , µ 2 ) are the collinear gluon and charm quark PDFs evaluated at longitudinal momentum fractions x and factorization scale µ 2 .
Including the mass of charm quark the on-shell matrix element takes the following form: where g 2 s = 4πα s (µ). In the massless limit m c → 0 one recovers the known textbook formula: A role of the charm quark mass in the matrix element will be discussed when presenting numerical results. x. Here different models for initial intrinsic charm quark distributions are shown (left panel) and a comparison between charm quark distributions obtained with and without concept of intrinsic charm in the proton. 1 Here and in the following we concentrate only on the forward production mechanisms (with charm quark having positive-rapidity), but, the formalism for symmetric backward configuration is the same.
In the numerical calculations below the intrinsic charm PDFs are taken at the initial scale m c = 1.3 GeV, so the perturbative charm contribution is intentionally not taken into account. We apply four different grids of the intrinsic charm distributions from the CT14nnloIC PDF [16] that correspond to the BHPS 1% and BHPS 3.5% as well as the sealike LS (low-strength) and sea-like HS (high-strength) models for initial intrinsic charm distribution. The distributions are compared with each other in the left panel of Fig. 4.
In the right panel we present in addition the difference between the CT14nnloIC charm PDF obtained with and without intrinsic-charm concept.
On the other hand the collinear gluon PDFs g(x, µ 2 ) are taken at the running factorization scale related to the averaged transverse momentum of the outgoing particles, i.e.
The charm quark mass m c = 1.3 GeV plays here a role of the minimal scale and ensures that we are not going beyond the fitted PDF grids where unconstrained extrapolation procedures are applied. We keep the charm quark mass here even when the massless matrix element and/or kinematics are used.
As will be shown later, the numerical results strongly depend on how the longitudinal momentum fractions x 1 and x 2 (arguments of parton distributions) are calculated. In the massive scheme of the calculations the quantities are defined as follows: In this equations p t1 is transverse momentum of the outgoing gluon (or light quark/antiquark) and the m t2 is c quark (c antiquark) transverse mass defined as m t = p 2 t + m 2 c . As will be discussed further it is crucial to include in kinematics the mass of the final charm, while the initial charm can be considered massless. In the following numerical studies all the calculations in the massless limit with massless matrix elements will be done within the kinematics corrected in the above manner. The effect of the correction will be also explicitly shown.
Considering forward production of charm at the LHC energies one is exploring asymmetric kinematical regions where x 1 is very small (down to 10 −5 ) and x 2 is rather large (about 10 −1 ). Thus in this reaction small-x gluon PDF and intrinsic large-x charm content of the proton are probed simultaneously. As it is shown in Fig. 5 both distributions are not well constrained by the global experimental data. In the left panel we show the central fits for intrinsic charm distribution from the CT14nnloIC and the NNPDF30nloIC PDF sets [31] together with 1σ standard deviation. In the right panel we compare gluon PDF fits from different collaborations, including MMHT2014nlo [32], JR14NLO08FF [33] and CT14lo/nnlo sets. Clearly the current level of knowledge of both distributions is rather limited and the large uncertainties prevent definite conclusions. In principle, a study of far-forward production of charm may improve the situation by exploring unconstrained areas.       In the present study we go beyond the leading-order mechanisms and include also higher-order processes that are expected to play important role. We take into account all 2 → 3 and 2 → 4 processes at tree-level that lead to a production of charm quark or  Having massless partons (minijets) in the final states considered in the present work it is necessary to regularize the cross section that has a singularity in the p t → 0 limit.
We follow here the known prescription adopted in PYTHIA where a special suppression factor is introduced at the cross section level [35]: for each of the outgoing massless partons with transverse momentum p t , where p t0 is a free parameter of the form factor.

B. The hybrid model
Within the asymmetric kinematic situation x 1 ≪ x 2 described above the cross section for the processes under consideration can be also expressed in the so-called hybrid factorization model motivated by the works in Refs. [36,37]. In this framework the small-x gluon is taken to be off mass shell and the differential cross section e.g. for pp → gcX via g * c → gc mechanism reads: where F g * (x 1 , k 2 t , µ 2 ) is the unintegrated gluon distribution in one proton and c(x 2 , µ 2 ) a collinear PDF in the second one. The dσ g * c→gc is the hard partonic cross section obtained from a gauge invariant tree-level off-shell amplitude. In the present paper we shall not discuss the validity of the hybrid model on the theoretical level and concentrate only on its phenomenological application in forward production. A derivation of the hybrid factorization from the dilute limit of the Color Glass Condensate approach can be found in Ref. [38].
The gluon uPDF depends on gluon longitudinal momentum fraction x, transverse momentum squared k 2 t of the gluons entering the hard process, and in general also on a (factorization) scale of the hard process µ 2 . In the numerical calculations we take different prescription [42].
All of the models, except the PB-NLO-set1, are constructed in the way that allows for resummation of extra hard emissions from the uPDFs. It means that in the hybrid model already at leading-order some part of radiative higher-order corrections can be effectively included via uPDFs. However, it is true only for those uPDF models in which extra emissions of soft and even hard partons are encoded, including k 2 t > µ 2 configurations. Then, when calculating the charm production cross section via e.g. the g * c → gc mechanism one could expect to effectively include contributions related to an additional extra partonic emission (i.e. g * c → ggc) which in some sense plays a role of the initial state parton shower. In Fig. 7 we plot the gluon transverse momentum dependence of the different gluon uPDFs from the literature. At the small x-values and low scales the differences between the model are quite significant. There are ongoing intensive works on construction of the full NLO Monte Carlo generator for off-shell initial state partons that are expected to be finished in near future [43].
This framework seems to be necessary in phenomenological studies that are based on the PB uPDFs [44]. The extra hard emissions from the DGLAP-based uPDFs are usually strongly suppressed which leaves a room for higher-order terms. Therefore, in this case one needs to include usual leading order subprocesses properly matched with a number of additional higher-order radiative corrections at the level of hard matrix elements. In the moment, it can be done only at tree-level.
In consequence, the numerical calculations with the PB-NLO-set1 uPDFs, are done including in addition all 2 → 3 and 2 → 4 channels of partonic subprocesses that lead to a production of charm quark or antiquark and are driven by the g * c and q * c (orq * c) initial state interactions, similarly as in the collinear case. Here we follow a dedicated matching procedure to avoid double-counting as introduced in Ref. [44], and further used in Refs. [45,46].

C. The k T -factorization
Another possible theoretical approach to perform the calculations for the processes considered here is the k T -factorization [47]. This framework extends the hybrid model formalism and includes in addition effects related to off-shellness of the initial state charm quark. In principle, it allows to study intrinsic charm contribution to charm production via mechanisms where both incident partons are off mass shell.
A topology of possible diagrams present in the k T -factorization in the case of intrinsic charm studies is not the same as in the collinear case. Here one can follow two different ways of calculation and consider: • g * c * → gc (and/or q * c * → qc,q * c * →qc ) mechanism, The second one is not present in other approaches and can be treated as leading-order.
The first mechanism directly corresponds to the scheme of the calculations applied in the hybrid model and can be classified as higher-order. However, their mutual coincidence is not clear and strongly depends on the model of unintegrated PDFs used in the numerical calculations.

The 2 → 2 partonic mechanism
The k T -factorization cross section for the pp → gcX reaction driven by the typical 2 → 2 mechanisms, e.g. like the g * c * → gc, can be expressed as follows: where F g (x 1 , k 2 1,t , µ 2 ) and F c (x 2 , k 2 2,t , µ 2 ) are the gluon and intrinsic charm quark uPDFs, respectively, for both colliding hadrons and M off−shell g * c * →gc is the off-shell matrix element for the hard subprocess. Here the Feynmann diagrams are the same as shown in Fig. 9.
The extra integration is over transverse momenta of the initial partons. Here, one keeps exact kinematics from the very beginning and additional hard dynamics coming from transverse momenta of incident partons. Explicit treatment of the transverse momenta makes the approach very efficient in studies of correlation observables.
Considering forward production of charm one should not expect that the initial state (intrinsic) charm quark could have large transverse momenta. Rather small deviations from the collinear limit are more physically motivated here. Therefore, for the unintegrated charm distribution F c (x, k 2 t , µ 2 ) we will assume Gaussian distributions with rather small smearing parameter σ 0 . The unintegrated c (c) distributions are constructed as: is a standard two-dimensional Gaussian distribution and σ 0 is in principle a free parameter which governs the nonperturbative effects in the proton wave function. The factor π is because of our normalization of unintegrated parton distributions: The hard off-shell matrix element M off−shell g * c * →gc is known only in the massless limit. Within this limit the relevant calculations can be done in the KATIE Monte Carlo code, where the matrix element is computed numerically. Its analytic form can be obtained according to Parton-Reggeization-Approach (PRA) and was published in Ref. [48]. For the higher-order tree-level diagrams with off-shell initial state partons and with extra partonic legs in the final state one can also use KATIE, which is very efficient in this type of calculations and which was very recently equipped with tools that allow in addition for generation of analytic form of matrix-elements for a given hard multileg processes [34].

The 2 → 1 partonic mechanism
In the k T -factorization framework the charm quark/antiquark can be created at oneorder higher approach. A relevant formalism was used previously for production of forward pions in Ref. [49]. In Fig. 8 we show basic graphs for charm quark production within the 2 → 1 mechanisms. The emitted charm-quark (or antiquark) momentum-space distribution can be written as: In the formula above f 's are unintegrated gluon or charm quark/antiquark distributions.
For unintegrated gluon distributions we will take the ones used recently in the literature in the context of η c or χ c production [50,51] where the kinematics is similar. For Ω 2 we The longitudinal momentum fractions are calculated as

IV. RESULTS
We divide the section with numerical results to four subsections. First three of them are devoted to numerical calculations obtained with the collinear-, hybrid-and the k Tfactorization approach, respectively. The last subsection contains explicit predictions for impact of intrinsic charm mechanism on forward production of charm in different experiments, including low energy LHC experiments like fixed-target LHCb and SHIP, as well as high energy FCC and LHC experiments, like proposed recently LHC-FASER.

A. The collinear approach
We start presentation of numerical predictions with the results for pp → gcX reaction driven by the gc → gc leading-order mechanism calculated in the collinear framework within massive matrix element and kinematics for the energy √ s = 7 TeV. Here we take the gluon and the intrinsic charm distributions as encoded in the CT14nnloIC collinear PDFs. The three different lines in Fig. 9 correspond to a different choice of the p t0 parameter used for the regularization of the cross section. We see that the predictions for charm quark transverse momentum (left panel) and rapidity (right panel) distributions are very sensitive to the choice of this parameter, especially, at small charm quark transverse momenta, which also affects the rapidity spectrum. In the numerical studies below p t0 = 1.0 GeV will be taken as a default choice which leads to a central value of the uncertainty related to the choice of the parameter.
In Fig. 10 we present again collinear results for the leading-order gc → gc mechanism but here we applied four different sets of the intrinsic charm distribution in a proton at initial scale µ = 1.3 GeV as incorporated in the CT14nnloIC PDFs. Again, we show the differential cross sections as a function of the charm quark transverse momentum (left panel) and rapidity (right panel). The solid, long-dashed, dotted and dash-dotted lines correspond to the BHPS 1%, BHPS 3.5%, sea-like LS and sea-like HS models, respectively.
The sea-like models lead to a larger cross section than in the case of the BHPS model in the midrapidty region. On the other hand, a larger cross section in the forward direction is obtained within the BHPS models. Clearly, large uncertainties due to the intrinsic charm input are found.  The intrinsic charm component in the proton is not the only source of uncertainties related to the collinear PDFs. As it is shown in Fig. 11 the gluon PDF also leads to a sig- Having discussed the dominant leading-order mechanism we wish to move beyond and consider importance of higher-order corrections for the charm quark forward production mechanisms with intrinsic charm in the initial state. In Fig. 13 we compare our collinear predictions for the leading-order 2 → 2 mechanisms, both gc → gc (dotted histograms) and qc → qc (short-dashed histograms) shown separately, and for the higherorder 2 → 3 (long-dashed histograms) and 2 → 4 (dash-dotted histograms) mechanisms calculated at tree-level. A sum of the four different components denoted as 2 → 2 + 3 + 4 is also shown but it does not follow any merging procedure here 2 . For the higher-order contributions the partonic subprocesses with gc and qc initial states are added together.
We report a huge contribution to the cross section coming from the higher-order mechanisms (more than order of magnitude). It clearly shows that the leading-order mechanisms are not enough in order to get reasonable predictions for the impact of intrinsic charm concept on forward charm quark production. Full NLO and even NNLO frameworks are required for precise studies of the subject within the collinear parton model.
The situation in the case for other approaches, like the hybrid-and the k T -factorization is quite different than in the collinear case what will be discussed in next two subsections.

B. The hybrid model
Now we wish to start presentation of our numerical results obtained in the hybrid model. Here the incident small-x parton is assumed to be off-mass shell in contrast to the large-x intrinsic charm which is kept on-shell. In Fig. 14 we show theoretical predictions for charm quark transverse momentum (left panel) and rapidity (right panel) distributions for forward charm production within the leading-order g * c → gc and the q * c → qc mechanisms. Here the KMR-CT14lo gluon and light quark/antiquark uPDFs are used.
We observe that here much larger cross sections are obtained than in the analogous calculations done in the collinear framework (see two lowest histograms in Fig. 13). Especially, in the hybrid model the gluonic component is much bigger than its collinear counterpart.
Significant effects related to the off-shellness of the incident gluons are found. Considering far forward rapidities of the produced charm quark the transverse momenta (virtualities) of the incident small-x gluons start to play a very important role and lead to a sizeable enhancement of the predicted cross section with respect to the leading-order collinear calculations.
Since in the hybrid model the leading-order quark component q * c → qc is found to be negligible one can safely concentrate on the gluonic g * c → gc channel only. In Fig. 15 we show the relevant predictions for different unintegrated gluon densities from the litera- Above, we have used those gluon uPDF models that are assumed to allow for an effective resummation of extra real emissions (real higher-order terms). Therefore, their can be successfully used in phenomenological studies based even on leading-order matrix elements (see a discussion in Refs. [44,46]). Here we wish to present results obtained within the scheme of the calculations where the higher-order corrections are not resummed in the uPDF but are taken into account via the hard-matrix elements. This procedure can be tested with the help of the DGLAP-based Parton-Branching uPDFs as was proposed in Ref. [44] and further applied in Refs. [45,46]. In Fig. 16 we show predictions of the hybrid model for the 2 → 2, 2 → 3 and 2 → 4 mechanisms, as well as for their sum 2 → 2 + 3 + 4 obtained using a dedicated merging procedure. The results are calculated with the PB-NLO-set1 quark and gluon uPDFs. For the leading-order 2 → 2 mechanisms we show g * c and q * c channels separately while for the higher-order components we plot sum of all possible gluonic and quark channels. As in the collinear case, the higher-order mechanisms are found to be very important also here.
For a better transparency in Fig. 17 we compare the hybrid model results obtained with the KMR-CT14lo (solid histograms) with the PB-NLO-set1 (dashed histograms) uPDFs, that correspond to the two different hybrid calculation schemes, together with the results obtained in the collinear approach (dotted histograms). Both types of the hybrid model calculations seem to lead to a very similar predictions. It seems to justify the proposed 2 → 2 + 3 + 4 hybrid calculation scheme with the PB uPDFs and with the applied merging in a qualitative way. On the other hand, the collinear 2 → 2 + 3 + 4 results seem to be larger by a factor of 2 than their hybrid model counterpart. However, this might be related to a lack of a relevant merging procedure in the collinear case.

C. The k T -factorization approach
Now we wish to present results obtained within the k T -factorization approach. So here we take into account also effects related to off-shellness of c quark of the intrinsic charm in the proton. The transverse momentum dependent intrinsic charm uPDF is obtained by Gaussian smearing of the collinear PDF. Rather small smearing parameter is used that do not allow for a large transverse momenta of the intrinsic charm. It seems to be appropriate for the case of the forward production of charm. For the unintegrated gluon density the KMR-CT14lo model is used.
In Fig. 18 we show predictions for the g * c * → gc mechanism with both initial state partons being off-shell. Here three different values of the smearing parameter in the calculation of the intrinsic charm uPDF are used: σ 0 = 0.5 GeV (solid histograms), 3.5 GeV (dotted histograms) and 7.0 GeV (dashed histograms). The larger σ 0 is taken the smaller cross section at small outgoing charm quark transverse momenta is obtained (left panel). The same is true for the rapidity spectrum in the forward region (right panel). In Fig. 19 we illustrate mutual relations between the results obtained with the hybrid and the k T -factorization frameworks. When the smearing parameter in the calculation of the intrinsic charm uPDF is small, e.g. σ 0 = 0.5 GeV, the hybrid model g * c → gc results coincide with the g * c * → gc results obtained within the full k T -factorization approach.
Finally, we wish to present results of the k T -factorization approach for the g * c * → c mechanism. In Fig. 20  to quite different results. The discrepancies between the uPDF models obtained here seem to be larger than in the corresponding case of the g * c → gc calculations within the hybrid model.

D. Predictions for future experiments
Before we go to predictions for different future and present experiments we wish to summarize the conclusions drawn in the previous subsection by a direct comparison of the results corresponding to the approaches discussed above. In Fig. 21 we compare predictions of the three different approaches used in the previous subsections: the 2 → 2 + 3 + 4 collinear (dashed histograms), the hybrid g * c → gc (solid histograms) and the k T -factorization g * c * → c (solid lines) calculations. Different models lead to a very different results with more than one order of magnitude difference between the lowest and the highest predicted cross section. Huge cross section for gc → c or cg → c may be partly due to ignoring other emissions than c orc in the evolution of x 1 and x 2 . These large uncertainties of the predictions can be reduced only by a forward experiments at forward directions. Forward charm production data sets that will be dominated by the contribution from intrinsic charm are necessary to draw definite conclusions about the level of applicability of the different theoretical approaches.
Therefore, now we wish to present results of the study of the impact of the intrinsic charm component on the forward charm particle production in already existing or future experiments at different energies. We start with predictions for the high energy experiments at the LHC and the FCC, at √ s = 13 and 50 TeV, respectively (top and bottom panels in Fig. 22). In the LHC case the considered kinematics correspond to the planned FASER experiment. Here we compare predictions of the k T -factorization approach for the g * g * → cc mechanism which is known to give a very good description of the LHC open charm data [44], and predictions of the g * c → gc mechanism (dashed) within the hybrid model. In both cases the charm production cross section starts to be dominated by the intrinsic charm component at very forward rapidities, i.e. y ≥ 7. In this far-forward region, the transverse momentum distribution of charm quark is also dominated by the contributions of the intrinsic charm. The predicted enhancement of the charm cross section could certainly be examined by the FASER experiment dedicated to a measurement of forward neutrinos originating from semileptonic decays of D mesons. The actual predictions for neutrinos will be presented elsewhere. In addition, we also analysed a possibility of experimental study of the intrinsic charm concept at lower energies. In Fig. 23  We observe that also at relatively small energies the intrinsic charm contributions could be identified experimentally. It seems that the already existing data set on open charm meson production in the fixed-target LHC mode [23]

V. CONCLUSIONS
In this paper we have discussed the effect of intrinsic charm in the proton on forward production of c quark orc antiquark at different energies. Three different approaches: collinear, hybrid and k T -factorization have been used with modern collinear and unintegrated parton distribution functions.
The production mechanism of c-quarks andc-antiquarks originating from intrinsic charm in the nucleon is concentrated in forward/backward directions, but details depend on collision energy. The absolute normalization strongly depends on the approach used. The leading-order (LO) collinear framework leads to the smallest cross section. The cross section becomes much bigger in the k T -factorization or in the hybrid model which effectively include higher-order corrections. The next-to-leading (NLO) and even nextto-next-to-leading (NNLO) tree-level corrections are found to be very important here.
Therefore, the k T -factorization or the hybrid model will give stringent limits on the intrinsic charm which cannot be constrained at present from first principles.
We have shown that in the collinear approach the LO calculations of the intrinsic charm component are insufficient. We have included the NLO and NNLO components at tree-level which were found to significantly contribute to the cross section.
Working in the hybrid model or in the k T -factorization approach we have shown that the effects related to the off-shellness of the incident partons (especially gluons) are large.
In both cases higher-order corrections are effectively included already within the basic gc → gc mechanism. We have used different models for gluon unintegrated parton distribution functions (uPDFs) from the literature. We obtained different results for different gluon uPDFs. The forward charm production was recognized as a useful testing ground for the small-x behaviour of the gluon uPDFs. We have shown in addition that the final results are also sensitive to the concept of gluon saturation in a proton. Unintegrated gluon densities derived from linear and non-linear evolution equations lead to a quite different results. We have performed also leading-order calculations within k T -factorization approach where the basic process is either g + c → c or c + g → c as done for forward production of charm quarks.
We have shown that the intrinsic charm component dominates over the standard pQCD (extrinsic) mechanism of cc-pair production at forward (or far-forward) rapidities neutrino production such as SHIP sand FASER are an interesting alternative in next few years.
In the present study we intentionally limited to the production of charm quarks/antiquarks. The production of charmed mesons or baryons is currently uder discussion and a new fragmentation scheme was proposed [58] very recently. We leave the predictions for production of charmed hadrons and their semileptonic decays for a separate study. However, the consequences for high-energy neutrino production have been discussed shortly in the context of the IceCube experiment and experiments proposed at the LHC (SHIP and FASER).