Predictions for nuclear structure functions from the impact-parameter dependent Balitsky-Kovchegov equation

In this work we present dipole scattering amplitudes, including the dependence on the impact-parameter, for a variety of nuclear targets of interest for the electron-ion colliders (EICs) being currently designed. These amplitudes are obtained by numerically solving the Balitsky-Kovchegov equation with the collinearly improved kernel. Two different cases are studied: initial conditions representing the nucleus under consideration and the solutions based on an initial condition representing a proton complemented by a Glauber-Gribov prescription to obtain dipole-nucleus amplitudes. We find that the energy evolution of these two approaches differ. We use the obtained dipole scattering amplitudes to predict ($i$) nuclear structure functions that can be measured in deep-inelastic scattering at EICs and ($ii$) nuclear suppression factors that reveal the energy evolution of shadowing for the different cases we studied. We compare our predictions with the available data.

tions representing the nucleus under consideration and the solutions based on an initial condition representing a proton complemented by a Glauber-Gribov prescription to obtain dipole-nucleus amplitudes. We find that the energy evolution of these two approaches differ. We use the obtained dipole scattering amplitudes to predict (i) nuclear structure functions that can be measured in deep-inelastic scattering at EICs and (ii) nuclear suppression factors that reveal the energy evolution of shadowing for the different cases we studied. We compare our predictions with the available data.

I. INTRODUCTION
Feasibility studies for electron-ion colliders (EICs), like those proposed in the USA [1] or at CERN [2], are an essential ingredient towards the design of these machines. Inclusive measurements are among the most important observables in this context. In particular, the study of the structure function F A 2 (x, Q 2 ) at small Bjorken-x for photons of virtuality Q 2 at a pertubative scale, and for a variety of nuclei A, is expected to yield a new understanding of the high-energy limit of quantum chromodynamics (QCD). Comparison of these measurements with those reported by H1 and ZEUS [3] for the corresponding structure function of the proton, F p 2 (x, Q 2 ), promise to shed new light on the origin of shadowing, the phenomenon that the parton distributions of nucleons bounded in a nucleus are suppressed with respect to those of free nucleons [4].
At small values of x, the dominant parton distribution is that of gluons; thus the case of gluon shadowing has been the focus of attention for theorists since a long time; e.g., [5,6].
A process expected to occur in this kinematic regime is saturation, namely the fact that the density of gluons is so high that they start to interact with each other, even in the domain of perturbative QCD. (For a recent review see [7].) An early equation to describe saturation was introduced in the seminal work [8], while nowadays it is common to use the Balitsky-Kovchegov (BK) equation for this type of studies. The leading order BK equation, discussed in Sec II, was derived in [9] and [10] using two independent approaches. Later on, corrections to account for the running of the coupling [11,12] as well as the resummation of other logarithmic contributions [13,14] were incorporated into this formalism. In the approximation of considering a large homogeneous target, that is, disregarding the impact parameter dependence, this equation has been successfully used to describe the existing F p 2 (x, Q 2 ) data, e.g. in [14,15].
The first attempt at solving the BK equation including the impact-parameter dependence [16] found that the solutions developed so-called Coulomb tails: an unphysical grow of the amplitude at large impact parameters. Nonetheless, using some extra ad hoc corrections it was possible to describe the structure function data of the proton [17,18]. Recently, our group discovered that using the collinearly-improved kernel introduced in [14] the problem of Coulomb tails is tamed such that a successful phenomenology using the BK equation is possible [19,20].
In this article, we solve the BK equation with the collinearly improved kernel for different nuclei of importance for future EICs and predict their structure functions as well as the corresponding nuclear suppression factors, which are a direct measurement of shadowing. We study two cases: solutions obtained from an initial condition representing the nuclei (denoted as b-BK-A below), and solutions of the proton case coupled to a Glauber-Gribov prescription to obtain the nuclear structure functions (denoted as b-BK-GG). Other approaches to the computation of nuclear structure functions can be found for example in [21][22][23][24][25] The rest of the text is organised as follows: Sec. II contains a brief review of the formalism, including the definition of the initial conditions and the values of the parameters used in the computation. Section III discusses the behaviour of the dipole scattering amplitudes obtained by solving the BK equation for the different nuclei and the b-BK-A and b-BK-GG approaches. Section IV presents our predictions for the nuclear structure functions and nuclear suppression factors for all cases under study as well as a comparison with the available data. Finally, in Sec. V we provide a brief summary of the presented work as well as an outlook of future steps.

II. FORMALISM
A. The Balistsky-Kovchegov equation with the collinearly improved kernel The leading order Balitsky-Kovchegov equation [9,10] is It describes the evolution in rapidity Y of the dipole scattering amplitude N( r, b, Y ). Here, the sizes of the mother and daughter dipoles are r ≡ | r |, r 1 ≡ | r 1 |, and r 2 ≡ | r 2 | ≡ | r − r 1 |, respectively. The magnitudes of the corresponding impact parameters between these dipoles and the hadronic target are b ≡ | b|, b 1 ≡ | b 1 |, b 2 ≡ | b 2 |. All these vectors are two-dimensional and live in the impact-parameter plane.
We solve the equation for the case when the following two conditions are fulfilled, (i) the evolution depends only on the magnitude of both the dipole size and the impact-parameter vectors, and (ii) the angle between r and b is fixed to zero: For the kernel we use the recently proposed collinearly improved version [13] K ci (r, r 1 , r 2 ) = α s 2π where (see also [26]) J 1 is the Bessel function, the anomalous dimension is A 1 = 11/12, and L r i r = ln (r 2 i /r 2 ). The sign is positive when the size of the original dipole is smaller than the size of each of the daughter dipoles and negative otherwise. The smallest dipole prescription is used for the running coupling: α s = α s (r min )N c /π, where r min = min(r 1 , r 2 , r). Note that this prescription has also been put forward as the natural scale for the BK equation at next-toleading order [27]. The variable-number-of-flavours scheme is used with the same parameters as in our previous work [19,20].
B. Glauber-Gribov approach to the nuclear dipole amplitude Following [21], one can use the solution of the BK equation for the case of a proton target to obtain the dipole scattering amplitude for a nuclear target by using a Glabuer-Gribov with This approach has been used in other studies, e.g. those reported in [21,22] (see also [24] for a more general approach that reduces to the Glabuer-Gribov case for large nuclei). The nuclear thickness function T A (b) is obtained from a Woods-Saxon distribution for the nuclear matter density  (7), used in the computations reported in this text, and the value of the Q 2 s0 (A) parameter obtained as explained in the text. The Wood-Saxon parameters are taken from [29].
with the x and y coordinates in the impact-parameter plane. It is normalised according to (See for example [28] for full details on the formalism.) The values of the Woods-Saxon parameters are given in Table I. This approach is denoted as b-BK-GG in what follows.

C. Initial conditions for the nuclear targets
To solve the BK equation an initial condition is needed. In our previous work [19,20] we introduced a new functional form for the initial condition given by where Q 2 s0 is a free parameter representing the saturation scale at a zero impact parameter, and b q i are the impact parameters of the quark and anti-quark forming the dipole.
For the case of the proton, we assumed a Gaussian like distribution which leads to The parameter B G was set to 3.2258 GeV −2 , while Q 2 s0 took the value 0.496 GeV 2 . With these values a satisfactory description of HERA and LHC data on the proton structure functions and exclusive production of vector mesons is achieved [19,20].
Here, we follow a similar approach for the nuclear case, but assuming a Woods-Saxon instead of a Gaussian distribution: with where the factor k ensures that kT A (0) = 1. This approach is denoted as b-BK-A in what follows.
As the nuclear parameters are already fixed, the only free parameter is Q 2 s0 (A). We have fixed these parameters using This dipole scattering amplitude at the initial rapidity is used to compute structure functions and to compare them with the predictions obtained using the EPPS16 nuclear parton distributions [30].
In detail, we have varied the value of the Q 2 s0 (A) parameter in order to get a small relative deviation from the structure function F 2 (x = 0.008, Q 2 ) as predicted by the EPPS16 PDFs.
The comparison is done for the following values of the photon virtuality: Q 2 ∈ [3.5, 4.5, 6.5, 8.5, 10, 12, 15, 18, 22, 27] GeV 2 to avoid the nonperturbative region at very low Q 2 and to stay in the region of virtualities where the BK equation is expected to work the best. We have used LHAPDFs [31] to obtain the PDF sets and the APFEL software [32,33] for the computation of the structure function. The values obtained for Q 2 s0 (A) by this procedure are reported in the last column of Table I. Interestingly, this parameter follows a linear behaviour as a function of A 1/3 as shown in  These differences between the dipole scattering amplitudes in the two approaches reflect themselves in one of the most important parameters that can be obtained from these objects: the saturation scale and its evolution. As it is standard, we define the saturation scale at a given rapidity and a fixed impact parameter as the dipole size that produces a scattering amplitude equal to a constant that commonly is chosen to be one half. Figure 4 shows the behaviour of the saturation scale at an impact parameter of 0.01/GeV for two rapidities as a

IV. PREDICTIONS FOR NUCLEAR STRUCTURE FUNCTIONS AND NU-CLEAR SUPPRESSION FACTORS
A. Relation between the dipole scattering amplitude and the structure function Using as input the dipole scattering amplitudes, the structure function F A 2 (x, Q 2 ) is computed as where, following [34],x = x(1 + (4m 2 q i )/Q 2 ) with m q i the mass of the i-quark. The dipoletarget cross section is related to the dipole scattering amplitude by Finally, the wave function representing the probability of a virtual photon splitting into a quark-antiquark dipole, and following the notation of [34], is and for the transverse and longitudinal polarisation of the incoming photon, respectively. The total wave function is In these equations K 0 and K 1 are the MacDonald functions, z is the fraction of the total longitudinal momentum of the photon carried by the quark, e q i is the fractional charge (in units of elementary charge) of quark i, α em = 1/137 and ǫ 2 = z(1 − z)Q 2 + m 2 q i . As in our previous work [19,20] we set the quark masses to 100 MeV/c 2 for light, 1.3 GeV/c 2 for charm, and 4.5 GeV/c 2 for bottom quark. As reported for example in [13] the numerical results do not depend strongly on these choices.

B. Predictions for the nuclear structure function
The nuclear structure functions F A 2 (x, Q 2 ) for Ca and Pb are shown in the upper panels of Fig. 5 as a function of x for two values of the photon virtuality Q 2 . Results in both approaches, b-BK-A and b-BK-GG, are shown in the figure which also shows the ratio of both predicted structure functions. There is a clear difference between both sets of results. Furthermore, the difference shows a dependence on x, on Q 2 , and a striking dependence on A, where the difference between both approaches grows from small to large nuclei.
A measurement of this structure function is expected to be one of the first results of any future EIC. Given the precision expected from these machines, these measurements will select which of the two approaches describes better the data.
The lower panels of Fig. 5 show that the difference between the b-BK-A and b-BK-GG approaches are larger for the charm structure function, F A 2C (x, Q 2 ), than for the inclusive case, and that there is a very soft dependence on kinematic variables and nucleus species.
A measurement of F A 2C (x, Q 2 ) would offer additional stringent constraints to predictions of the structure function of nuclei.

C. Predictions for the nuclear suppression factor
As a final observable we present the nuclear suppression factor, defined as the ratio , which is expected to be unity if the structure of a free nucleon is equal to that of a bounded one. This ratio is the most direct way to observed nuclear shadowing, which for small x is dominated by gluon shadowing and thus may be an important tool to determine the behaviour of saturation across different nuclei.
This factor is shown in the upper panels of Fig. 6 as a function of x. Exiting data at the same Q 2 from [35] is also shown as a cross check of the procedure. For the x dependence of R pA one sees a linear decrease (in logarithmic scale) towards small x for both nuclei, but the linear behaviour is reached later for the lighter nucleus, specially at higher Q 2 scales.
The b-BK-A computation predicts stronger shadowing than the b-BK-GG case with this behaviour seemingly dependent on Q 2 . The same figure also shows, in the bottom panels, the A-dependence of the nuclear suppression factor for two Q 2 scales and for two values of x. As expected, shadowing becomes stronger as the size of the nucleus grows. The different behaviour of shadowing for different nuclei in the b-BK-A and b-BK-GG is clearly seen in this observable. Figure 7 shows the comparison of our predictions with those obtained using EPPS16 which is considered a standard of our current knowledge of nuclear shadowing. The comparison is done for Ca and Pb as middle and large nuclei. The predictions are compared with data from [35]. Note that the predictions are at a Q 2 scale of 2.42 GeV 2 which we considered the lowest we would like to go to stay in a somehow perturbative scale. But the data is measured at a different Q 2 for each x value (as illustrated by the use of empty markers for data at smaller Q 2 ). The Q 2 values are reported in the figure.
Focusing on the predictions, the behaviour at small x is definitely different for the EPPS16 and BK computations. Note that the difference between EPPS16 and b-BK at the initial scale used for the BK evolution have two origins: one, that the parameter of the initial scale shown in Fig. 1 is chosen by comparing with larger values of Q 2 than those shown in Fig. 7, and two that the prediction for the structure of the proton is substantially different for EPPS16 and for the b-BK approach reported in [19,20]. Given that the difference among the approaches goes beyond a normalisation factor and shows a strong x-dependence, data from future EICs are expected to be precise enough to decide which prediction is closer to reality.
Comparing with the currently available data, and taking into account (i) the different Q 2 in data and predictions, and (ii) that for measurements the values are quite low (even below what one would expect to be valid for an approach based on perturbative QCD), the b-BK-A prediction seems to do a reasonable job of describing data. The EPPS16 prediction also does quite well for Pb, but slightly worse for Ca. The b-BK-GG prediction on the other hand is good when comparing with Ca, but it suffers a bit when compared with Pb.

V. SUMMARY AND OUTLOOK
The dipole scattering amplitudes, including the impact parameter dependence, for different nuclei have been obtained by solving the BK equation with the collinearly improved kernel. These amplitudes have been used to predict structure functions and nuclear saturation factors in kinematic ranges of interest for future EICs as those currently planned in the USA and at CERN. We followed two approaches: modelling the target directly as a nucleus using Wood-Saxon parameterisations (denoted as b-BK-A above), and solving for a proton and using a Glauber-Gribov prescription to go to the nuclear level (denoted as b-BK-GG above).
We find sizable differences between these approaches. These differences show a dependence on x, Q 2 and A such that data from a future EIC will be able to select the most appropriate approach for the description of data. We also compared nuclear suppression factors with those predicted using the EPPS16 formalism which is taken as a standard of our current knowledge of nuclear shadowing. We find that all three approaches yield different predictions and that the b-BK-A computation seems to provide a better description of and Q 2 = 7.91 GeV 2 (right) as a function of x for Ca and Pb (upper panels) and as a function of A at different fixed values of x (lower panels). The predictions are compared with data from [35].
existing data.
These studies show that the data expected from a future EIC have the capability of select the best theoretical approach and thus to advance our understanding of the nuclear structure, of shadowing, and of the high energy limit of QCD.
The dipole scattering amplitudes computed in this work are publicly available in the website https://hep.fjfi.cvut.cz/ along with macros and instructions to facilitate their use for anybody interested. Only the data point represented by a solid marker correspond to the Q 2 used for the predictions.
See text for details.