Numerical evaluation of the nonlinear Gribov-Levin-Ryskin-Mueller-Qiu evolution equations for nuclear parton distribution functions

We numerically study for the first time the nonlinear GLR-MQ evolution equations for nuclear parton distribution function (nPDFs) to next-to-leading order accuracy and quantify the impact of gluon recombination at small $x$. Using the nCTEQ15 nPDFs as input, we confirm the importance of the nonlinear corrections for small $x \lesssim 10^{-3}$, whose magnitude increases with a decrease of $x$ and an increase of the atomic number $A$. We find that at $x=10^{-5}$ and for heavy nuclei, after the upward evolution from $Q_0=2$ GeV to $Q=10$ GeV, the quark singlet $\Omega(x,Q^2)$ and the gluon $G(x,Q^2)$ distributions become reduced by $9-15$%, respectively. The relative effect is much stronger for the downward evolution from $Q_0=10$ GeV to $Q=2$ GeV, where we find that $\Omega(x,Q^2)$ is suppressed by 40%, while $G(x, Q^2)$ is enhanced by 140%. These trends propagate into the $F_2^A(x,Q^2)$ nuclear structure function and the $F_L^A(x,Q^2)$ longitudinal structure function, which after the downward evolution become reduced by 45% and enhanced by 80%, respectively. Our analysis indicates that the nonlinear effects are most pronounced in $F_L^A(x,Q^2)$ and are already quite sizable at $x \sim 10^{-3}$ for heavy nuclei. We have checked that our conclusions very weakly depend on the choice of input nPDFs. In particular, using the EPPS21 nPDFs as input, we obtain quantitatively similar results.