Probing shadowed nuclear sea with massive gauge bosons in the future heavy-ion collisions

The production of the massive bosons Z0\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$Z^0$$\end{document} and W±\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$W^{\pm }$$\end{document} could provide an excellent tool to study cold nuclear matter effects and the modifications of nuclear parton distribution functions (nPDFs) relative to the parton distribution functions (PDFs) of a free proton in high-energy nuclear reactions at the LHC as well as in heavy-ion collisions (HIC) with much higher center-of-mass energies available in the future colliders. In this paper we calculate the rapidity and transverse momentum distributions of the vector boson and their nuclear modification factors in p + Pb collisions at sNN=63\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\sqrt{s_{NN}}=63$$\end{document} TeV and in Pb + Pb collisions at sNN=39\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\sqrt{s_{NN}}=39$$\end{document} TeV in the framework of perturbative QCD by utilizing three parametrization sets of nPDFs: EPS09, DSSZ, and nCTEQ. It is found that in heavy-ion collisions at such high colliding energies, both the rapidity distribution and the transverse momentum spectrum of vector bosons are considerably suppressed in wide kinematic regions with respect to p + p reactions due to large nuclear shadowing effect. We demonstrate that the massive vector boson production processes with sea quarks in the initial state may give more contributions than those with valence quarks in the initial state; therefore in future heavy-ion collisions the isospin effect is less pronounced and the charge asymmetry of the W boson will be reduced significantly as compared to that at the LHC. A large difference between results with nCTEQ and results with EPS09 and DSSZ is observed in nuclear modifications of both rapidity and pT\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$p_T$$\end{document} distributions of Z0\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$Z^0$$\end{document} and W in the future HIC.


Introduction
The production of massive vector bosons Z 0 /W ± in highenergy nuclear collisions has long been regarded as an excellent tool to probe the initial-state cold nuclear matter (CNM) a e-mail: pengru@mail.dlut.edu.cn b e-mail: bwzhang@mail.ccnu.edu.cn effects in the ultra-relativistic heavy-ion collisions (HIC) at the LHC because Z 0 /W ± do not participate in strong interactions and the mean-free-path of massive vector bosons (or more precisely, the leptons decayed from Z 0 /W ± ) is rather long [1][2][3][4][5][6][7][8][9][10][11]. In proposed heavy-ion programs at the Future Circular Collider (FCC) at CERN [12], and the Circular Electron Positron Collider with a subsequent Super proton-proton Collider (CEPC-SPPC) in China [13], the vector gauge boson productions may still play a very important role in investigating the high-density QCD in the initial state and making precise constraints on CNM effects and nuclear parton distribution functions (nPDFs) at very small x. We note that the massive vector bosons are produced very early at time ∼ 1/m Z ,W ∼ 10 −3 fm/c with a decay time smaller than 0.1 fm/c, and a mean-free-path ∼10 fm in the QGP at temperature 1 GeV [14]. Even though hot/dense QCD matter to be created in nucleus-nucleus collisions at the FCC and the CEPC-SPPC may have a much higher initial temperature and longer lifetimes [15], the vector boson with the final state colorless dilepton would still be nearly blind to the medium evolution [14, [16][17][18].
Recent studies have shown that the heavy gauge boson production in nuclear reactions at the LHC energies could shed light on several CNM effects, especially the shadowing and anti-shadowing effects (for Z /W ) and isospin effect (for W ± ) [11,[19][20][21][22][23]. In the conceptual designs of FCC [15] and CEPC-SPPC [24] the center-of-mass energy for protonproton collisions could reach up to 100 TeV, which may give an energy of about √ s N N = 63 TeV for proton-lead collisions and about √ s N N = 39 TeV for lead-lead collisions. With much higher colliding energies available in future HIC, the vector boson would be produced by initial partons with the much smaller momentum fraction x, in which region the shadowing effect will be more pronounced and sea quarks and gluons may play a more important role. It is of great interest to see how the massive gauge boson productions help us understand the CNM effects and impose stringent constraints on nPDFs in future heavy-ion collisions. In this paper we study the nuclear modifications of the vector boson production in the future Pb + Pb collisions at √ s N N = 39 TeV and p + Pb collisions at √ s N N = 63 TeV. The numerical calculations are performed by using the perturbative quantum chromo-dynamics (pQCD) program DYNNLO [25,26] incorporating the parameterized nuclear parton distribution functions (nPDFs) sets EPS09 [27], DSSZ [28], and nCTEQ [29,30].
This paper is organized as follows. In Sect. 2, we briefly introduce the framework of our calculation and then compute the nuclear modification ratios for Z 0 boson rapidity distributions and the charged lepton pseudorapidity dependence of W ± production in Pb + Pb and p + Pb collisions in the future HIC. In Sect. 3, we discuss the nuclear modification of the transverse momentum distributions of Z 0 and W ± bosons in the future HIC. We present our summary in Sect. 4.

Nuclear modification ratio for vector boson rapidity distributions
Before discussing the nuclear modification of vector boson production, we briefly introduce the theoretic framework of this study, which is the same as in Ref. [11]. In hadronic collisions, the differential cross section of vector boson production A + B → V + X → ll + X through the Drell-Yan mechanism [31,32] could be calculated within the QCDimproved parton model [33,34]. In this approach, the cross section could be factorized into the convolution of the parton distribution functions (PDFs) f (x, μ) and partonic cross sectionσ ab→V +X →ll+X , which could be computed with the perturbative expansion, where x a(b) are the momentum fractions of the parton a(b) from the hadron A(B). As the baseline of the study, the hadron-hadron cross sections are numerically calculated with the program DYNNLO [25,26], which provides the O(α s ) and O(α 2 s ) pQCD corrections to the leadingorder (LO) cross section [σ 0 ] qq→V →ll . The O(α s ) calculation includes real corrections and one-loop virtual corrections, and that at O(α 2 s ) includes double-real corrections, real-virtual corrections, as well as two-loop virtual corrections [25,26]. In addition, the parton distribution functions parametrization sets MSTW2008 [35] are used, and the renormalization and factorization scales are set at the boson mass m V .
60 GeV < m ll < 120 GeV p+p √s =7TeV Fig. 2 The normalized Z 0 boson transverse momentum spectrum in p + p collisions at √ s = 7 TeV. The CMS combined data are taken from Ref. [37] We plot in Fig. 1 the Z 0 rapidity distribution, and in Fig. 2 its transverse momentum spectrum, for proton-proton collisions at √ s = 7 TeV. For rapidity dependence, both O(α s ) and O(α 2 s ) pQCD calculations agree well with the ATLAS data [36], which is measured in the invariant mass interval of the lepton pair as 66 GeV < m ll < 116 GeV and also in the lepton transverse momentum region as p l T > 20 GeV. The calculation at O(α 2 s ) also gives a good description of the CMS transverse momentum distribution [37] in the region p Z T > 10 GeV. The CMS measurement is made in the finalstate phase volume defined by the invariant mass, transverse momentum, and pseudorapidity of the dilepton as: 60 GeV < m ll < 120 GeV, p l T > 20 GeV, and |η l | < 2.1. The DYNNLO program has also been shown to give good agreement with the CMS data on the Z 0 rapidity distribution, and with the transverse momentum spectra of Z 0 and W measured by ATLAS [11].
To include several important CNM effects in Pb + Pb and p + Pb collisions, we phenomenologically utilize the parametrization sets of the nuclear parton distribution func-  [27][28][29][30][38][39][40]. In this work, we chose MSTW2008 proton PDFs and multiply them by the flavor and scale dependent factors R f (x, μ) taken from EPS09, DSSZ, and nCTEQ to obtain the parton distribution f p,A (x, μ) of a bound proton in nucleus. With the isospin symmetry being assumed [27][28][29][30], the nuclear parton distribution of a bound neutron f n,A (x, μ) could be obtained. In this paper only minimum-bias (MB) nuclear collisions are considered. Then we focus on the nuclear modification of vector boson production in the future HIC. First we study the Z 0 rapidity distribution in the Pb + Pb collisions at √ s N N = 39 TeV. The Z 0 signal is defined in the invariant mass window of the final state lepton pair as 60 GeV < m ll < 120 GeV, which is the same as the CMS experiments [2,8]. The Z -γ * interference is included in the calculation. We calculate at the order O(α s ), the nuclear modification ratio R A A (y Z ) defined as with N coll the number of binary nucleon-nucleon collisions [41]. In Fig. 3 both the results at the LHC and the future energies are plotted for comparison. Significant differences between R A A (y Z ) at two energies can be seen. The LHC results (left panel) with three nPDFs sets show weak nuclear modifications ( 10 %) in the studied rapidity region. Enhancements with EPS09 and nCTEQ at mid-rapidity can be observed, where suppression is given by DSSZ. However, at the future HIC (right panel), strong suppressions can be found in the whole studied rapidity regime (|y| < 4.0). Especially, nCTEQ nPDFs may suppress more than 30 % of Z 0 yield at the central rapidity, relative to proton-proton collisions. EPS09 and DSSZ give similar nuclear variations for the future HIC. The LO results with EPS09 are also plotted in Fig. 3. It is observed that the higher-order corrections make a negligible alteration on the R A A at the LHC, but give a larger suppression on R A A at the future HIC.
To give a simple picture of the difference between the two collisions energies, we resort to a lowest-order analysis [11]. In the LO pure electro-weak processes, the rapidity of Z 0 is related to the initial quark and anti-quark momentum fractions with the kinematic equation where the narrow width approximation (m ll ≈ m Z ) is used. Equation (3) shows that the momentum fraction of the initial parton is inversely proportional to the colliding energy √ s N N . For example, one has x ∼ 0.033 for mid-rapidity Z 0 production at 2.76 TeV, which falls into the anti-shadowing region of EPS09, whereas at 39 TeV x ∼ 0.0023 is given by Eq. (3), which is small enough to enter the shadowing region of EPS09. The suppression in the case of 39 TeV Pb + Pb collisions is due to the shadowed smaller-x initial partons. To visually display the difference between the two energies, we plot in Fig. 4 the factor R f v,sfs (y Z , μ) written as where R f (x, μ) is the flavor dependent nuclear variation factor given by EPS09, DSSZ, or nCTEQ, and x 1,2 is related to y Z by Eq. (3). The subscripts v and s stand for the valence and sea quarks (no gluon at leading order), respectively. This R f v,sfs (y Z , μ) factor could reflect the nuclear modification on the LO partonic subprocess initiated with f v,s andf s ; thus it could shed light on the R A A (y Z ) [11]. For the four processes plotted in Fig. 4, similar differences between the two collision energies could be observed. Compared to those at the LHC, the PDFs in future Pb + Pb collisions are more depleted in the kinematic region of the mid-rapidity Z 0 production, which results in a larger suppression of R A A (y Z ) in the future HIC than that at the LHC as observed in Fig. 3. We also notice that, although at the future HIC the u v +ū s channel with nCETQ gives less suppression than those with 0.6 us+us ds+ds s+s uv+us LHC: Future: The subprocess dependent factor R fv,sfs (y Z , μ) for Pb + Pb collisions at both √ s N N = 2.76 TeV and 39 TeV with factorization scale set at m Z EPS09 and DSSZ, their contribution combined with relatively large suppressions by other channels with nCTEQ, may give a more suppressed R A A (y Z ) than those with EPS09 and DSSZ. The underlying reason is that the u v +ū s channel gives relatively small contribution to Z 0 yield in the future HIC, and we will discuss this issue in more detail in the following.
We have mentioned that at LO R A A as a function of the vector boson rapidity mainly depends on the nuclear modification on the quark and anti-quark distributions [11], considering the LO process qq → V . Actually anti-quarks only come from the nucleon sea, while quarks could be valence or from the sea. To better understand the contributions to R A A given by some certain flavors, we calculate at LO the contribution ratio of partonic subprocess to the differential cross section. Three primary kinds of processes for both LHC and future HIC are shown in Fig. 5: (1) u v +ū s and d v +d s ; (2) u s +ū s and d s +d s , and (3) s +s and c +c. Processes (2)-(3) completely depend on the sea quarks, while the process (1) contains valence quarks. An obvious difference can be seen between the two energies. At LHC energy, the processes with valence quarks in the initial state dominate the Z 0 production. However, in the future HIC their contribution will be reduced significantly, while the processes with puresea quarks in the initial state become more important. The results demonstrate that sea quarks play an important role in the Z 0 boson production, especially at the future HIC. A main reason is that the sea quark densities increase very fast with the decreasing momentum fraction x, whereas the valence quark densities decrease sharply (more discussion can be found in Appendix 1). We show at the bottom of   at the future HIC is led by the shadowed nucleus sea quarks. It is noted that the nuclear modification of gluon density will also give contributions through high-order corrections.
We note that at the small-x region the uncertainties of nPDFs are rather large, thus the R A A given by different nPDFs may even merge with each other if error bar is included. That means the nuclear effect is very unclear at small-x, and the future Pb + Pb collisions will provide invaluable information to constrain the nPDFs for sea quarks and gluons in the small-x. Especially, the nuclear variations of s and c quarks densities, which have been constrained loosely by vector boson production at the LHC, should become important at the future HIC.
The yield of Z 0 boson in the future p + Pb collisions at √ s N N = 63 TeV is also studied. Figure 6 shows the nuclear modification factor R p A as a function of the Z 0 rapidity in the center-of-mass frame, predicted by EPS09, DSSZ, and nCTEQ nPDFs. Results for LHC energy ( √ s N N = 5.02 TeV) are also shown for comparison. Since the nuclear effects come from the colliding lead nucleus, the asymmetric rapidity dependence of R p A could be observed at both LHC (left panel) and future HIC (right panel) energies. For the LHC energy, enhancement can be seen in the backward region, while suppression is given in the forward region (related to anti-shadowing and shadowing, respectively [11]). For future energy, except the enhancements in the very backward regime (−4.0 < y Z < −3.0), persistent suppressions are given by EPS09, DSSZ, as well as nCTEQ.
The LO analysis is simple for p + Pb collisions, and one could study the R p A from the nPDFs factors R f (y Z , μ) by replacing the momentum fraction x 2 (related to the nuclear parton) with y Z , according to Eq. (3). In Fig. 7  the R f (y Z , μ) factors of certain flavors for both LHC and future HIC. Compared to the LHC energy, the future p + Pb collisions generate more suppression in the rapidity region y Z −2.0 due to the shadowed initial parton from the lead nucleus with smaller-x. Though for s quark density in the LHC backward region, and u v quark density in forward region of both energies, the nuclear modifications with nCTEQ show distinct behavior, their contributions to Z 0 rapidity distribution will be rather small in the corresponding region. To see this clearly, three kinds of contribution of the partonic subprocess initiated by certain nuclear flavors are illustrated in Fig. 8. It is observed that for both LHC and future HIC, the valence quarks' (u v and d v ) contribution decreases with the rapidity (smaller-x), but the sea quarks' increases. Comparison between the two energies implies that the nucleus valence quarks densities give negligible contributions with the increasing energy (smaller-x). The suppression in the future p + Pb collisions is to a great extent the result of the shadowing effect on the sea quarks. The contribution  Fig. 8 The LO rapidity dependent contribution rates for partonic subprocess initiated with various nuclear flavors in p + Pb collisions at both √ s N N = 5.02 and 63 TeV of the s and c quarks will also become non-trivial (more than 30 % in the mid-rapidity at LO). The R A A and R p A as functions of the Z 0 rapidity would become a good probe of the nuclear modification on small-x sea quarks in future heavy-ion collisions.
For completeness, we also study the nuclear modification on the W boson production. Usually, the charged lepton pseudorapidity, instead of the boson rapidity, is measured for W boson production in experiment with the neutrino in the final state (W → lν). We calculate the R A A and R p A as functions of the charged lepton pseudorapidity for W production in Pb + Pb collisions at √ s N N = 39 TeV and p + Pb collisions at √ s N N = 63 TeV, respectively. The transverse momentum region of the charged lepton is chosen to be p l T > 25 GeV [3,6].
For nucleus-nucleus collisions (Fig. 9), asymmetry between R A A (W + ) and R A A (W − ), due to the nuclear isospin effect, can be observed. To be specific, because of the existence of neutrons in nuclei, the nuclear reaction may enhance  the production of W − and reduce the yield of W + , relative to those in the proton-proton collisions. From the panels (d) and (h), one can see that the separation between W + and W − is large at the LHC, and becomes small at the future HIC, especially in the mid-rapidity region. Actually, the isospin effect is related to the parton distribution asymmetry [u(x) = d(x), andū(x) =d(x)], which originates mainly from the asymmetric u/d valence quark distributions in a nucleon [(uud) for proton, and (udd) for neutron]. As is discussed before, at the future energy the initial parton tends to come from the smaller-x regime, where the sea quarks dominate the valence quarks. Thus, at smaller-x the isospin effect becomes weaker (more detail can be found in 1), and nuclear modifications on W + and W − become more symmetric. It is also noted that the nuclear effect of the total W (= W + + W − ) production is very similar to that of the Z 0 rapidity dependence at the same collision energy. The differences among the R A A with three nPDFs sets are also consistent, for various particles (W , W + , W − , and Z 0 ) at the same collision energy. Three nPDFs sets predict similar nuclear modifications at the LHC [panels (a-c)], whereas nCTEQ gives much stronger suppressions in the future HIC [panels (e-g)]. Since isospin effect is very weak at the central rapidity, the suppression on the W boson production at the future HIC is mainly due to the shadowing effect on the nuclear sea quarks. For proton-nucleus collisions (Fig. 10), the nuclear variation on W ± is again the result of the isospin effect versus (anti-)shadowing effect. At the LHC, one can observe that the isospin effect results in a separation between W + and W − in the backward direction (larger x of the nucleus), and shadowing effect suppresses their production in the forward region (smaller x of the nucleus). However, the isospin effect is rather weak at the future energy, and suppressions on both W + and W − in a wide range can be observed.
The charge asymmetry observable A of W ± production, defined as is also studied for p + p, p + Pb and Pb + Pb collisions at both √ s N N = 5.02 TeV and √ s N N = 39 TeV, and results are plotted in Fig. 11 as a function of the charged lepton pseudorapidity. One can observe the charge asymmetries at the LHC show quite distinct behaviors in three different colliding systems, whereas the three curves at the future energy lie near the horizon and are close to each other. It has been found that the nuclear variation of the charge asymmetry observable A of W ± is mainly due to the isospin effect, and not sensitive to the other nuclear effects [11]. Therefore it is understandable that the isospin effect will become weak in the future HIC, as shown in Fig. 11.

Nuclear modification ratio for vector boson transverse momentum distributions
In this section we investigate the CNM effects on vector boson transverse momentum spectrum, which may be used to get access to quite different kinematic regions from that by the boson rapidity production [11]. In this paper we mainly focus on the massive boson production with p V T > 10 GeV, where DYNNLO provides an excellent baseline description for p V T distribution in elementary hadron-hadron collisions. First, we study the nuclear modification on Z 0 boson transverse momentum distribution in future Pb + Pb collisions at   √ s N N = 39 TeV. Besides the same Z 0 invariant mass interval as used in the last section, the rapidity range |y Z | < 2.0 is also chosen. The factor R A A , as a function of Z 0 boson transverse momentum p Z T , is calculated at the order O(α s ). Although the differential cross section will be changed a little by higher-order corrections, R A A defined as the ratio of the cross section of A+A to that of p + p would not be very sensitive to the higher-order corrections or the variation of hard scales [11,19]. Predictions by EPS09, DSSZ, and nCTEQ, for both LHC and future HIC, are shown in Fig. 12. Obvious difference between the two collision energies can be seen. Appreciable enhancements are given by three nPDFs sets in the region 10 GeV p Z T 250 GeV at the LHC, whereas considerable suppressions are given by them in the region 10 GeV p Z T 150 GeV for the future Pb + Pb collisions. The nCTEQ nPDFs supports stronger nuclear modifications than the other two in the studied p T regime.
To well understand the differences among the nPDFs sets and also between the two collision energies, one can again perform the analysis at LO. At this order, we obtain the kinematic relation at mid-rapidity with the narrow width approximation. Then according to Eq. (6), we can replace the initial parton momentum fraction x with p Z T and plot the flavor dependent factor R f ( p Z T , μ) of nPDFs in Fig. 13. For both LHC and future HIC, those factors of EPS09, DSSZ, and nCTEQ are shown. Compared to the nuclear effect at the LHC, higher suppression on the parton distribution could be seen in the region p Z T 200 GeV for the future nucleus-nucleus collisions, except for the nCTEQ u valence quark distribution. Because the related initial parton momentum fractions are different (smaller at future energy), the shadowing effect dominates the future Pb + Pb collisions, whereas the antishadowing effect shows up at the LHC.
To see clearly how the nuclear modification on each flavor enters into the observable R A A ( p Z T ), we again calculate the contribution rates of the various partonic subprocesses for both the LHC and the future HIC. Three kinds of contribution are plotted in Fig. 14: (1) partonic subprocesses with valence quarks in the initial state, such as u v +ū s or d v + g; (2) pure u, d sea initiated subprocesses as u s +ū s and d s +d s ; and (3) gluon initiated partonic processes with On the other hand, processes of Type (2) and Type (3) and other sea quark initiated processes dominate in the futureenergy collisions. In previous work [11,22] on the nuclear effects on p T spectra of gauge bosons at LHC, it has been shown that the gluon's contribution is predominant. Here we demonstrate that the impact of the nuclear modification of the sea quark distribution on the R A A ( p Z T ) is comparable to that of the gluon, in the future nucleus-nucleus collisions. For example, the large enhancement in the region 100 GeV p Z T 300 GeV at the LHC with nCTEQ is mainly due to the gluon nuclear modification but not the s quark (also enhanced a lot by nCTEQ), while at the future HIC the shadowing of both gluon density and sea quark (u, d, s, c) distributions in nCTEQ gives significant suppressions to the p T distributions of Z 0 .
The nuclear modification factor of the Z 0 boson transverse momentum distribution in the future p + Pb collisions at √ s N N = 63 TeV is also calculated at O(α s ), as shown in Fig. 15. Similar to that in Pb + Pb collisions, the suppression due to the shadowing effect could be observed at the p + Pb in the future HIC. The nuclear modification in p + Pb collisions is weaker than that in Pb + Pb collisions because in p + Pb only one colliding object is a lead nucleus. The partonic subprocess contribution rates are calculated for both LHC and future p + Pb reactions. Figure 16 shows the rates for three kinds of process, each of which is initiated by certain flavors in the lead nucleus. One can see that for both of the two energies, processes initiated with sea quarks or gluons give the main contribution. Especially at future energy the valence quarks' contribution is marginal; the nuclear effects are brought about by the sea quarks and gluon.
Next we turn to the study of the nuclear modification on the W boson transverse momentum spectrum. The W boson production is calculated in the final-state phase space defined  by the center-of-mass frame pseudorapidity and transverse momentum of the charged lepton as |η l | < 2.5 and p l T < 25 GeV [6]. In Fig. 17 x). Although W ± production obviously depends on the isospin effect, that of the total W receives a Z 0 -like nuclear variation mainly because of the (anti-)shadowing effect. One can observe that the R A A of the W − boson are very different for two collision energies. It is because the isospin effect strongly enhances the W − production at the LHC, while the shadowing effect is more pronounced than the isospin effect in the region p W T < 100 GeV at the future HIC. It should also be noted that the suppression on W + at the LHC is due to the isospin effect, but that at the future energy is caused by the shadowing effect.
The nuclear modification ratio R p A ( p W T ) is also calculated for the future p + Pb collisions at √ s N N = 63 TeV, and the results are shown in Fig. 18. Results for the LHC energy are also shown for comparison. Same trends as in Pb + Pb collisions could be seen in p + Pb. It should be noted that although both isospin effect and (anti-)shadowing effect are weaker in p + Pb as compared to Pb + Pb, the shadowing effect still surpasses the isospin effect in the p W T < 100 GeV region and thus suppresses the W − yield in the future HIC.
Next we study the charge asymmetry observable A as a function of the W transverse momentum. The results for p + p, p + Pb and Pb + Pb collisions at both √ s N N = 5.02 TeV and √ s N N = 39 TeV are plotted in Fig. 19. One can observe that the nuclear modification of A due to the isospin effect increases with p W T and decreases with the collision energy. Actually A 's nuclear modification also increases with its value in proton-proton collisions, which is related to the parton distribution asymmetry ratios r ud (x) and rūd (x) (definitions can be found in Appendix A). Figure 19 demonstrates the isospin effect will be reduced in the future HIC.
From the numerical results we note that the differences among various nowadays nPDFs sets are considerable. As an example, we calculate the ratio of the R A A ( p W T ) at the LHC ( √ s N N = 2.76 TeV) to that at future energy ( √ s N N = 39 TeV) for total W production, as The results for p W T 0 = 90 GeV are shown in Fig. 20. One can observe that EPS09 and DSSZ predict similar ratios of R A A at the LHC and future HIC, with a very flat p W T dependence. However, nCTEQ gives a very different result: the ratio is obviously larger than those with the other two nPDFs, and it increases faster with p W T . The measurement on the nuclear modification of the vector boson production in the future heavy-ion collisions will provide very useful constraints on nPDFs.

Summary
We investigate in this paper the initial cold nuclear matter effects on the heavy vector boson production in the future nuclear collisions. Within the framework of perturbative QCD, we numerically calculate the nuclear modification ratios of the Z 0 boson rapidity distribution, the charged lepton pseudorapidity dependence for W production, as well as Z 0 /W ± transverse momentum spectra. Both Pb + Pb collisions at √ s N N = 39 TeV and p + Pb collisions at √ s N N = 63 TeV are studied, where the CNM effects are included by utilizing the parameterized nuclear parton distribution functions sets (EPS09, DSSZ, and nCTEQ).
Different nuclear modifications from those at the LHC energies could be found. First, the nuclear shadowing effect would play a more important role in massive gauge boson production in the future HIC. As consequences, both the (pseudo-) rapidity dependence and the transverse momentum distribution of the boson would receive suppressions in a wide regime. Second, the nuclear isospin effect would be rather weak in the future nuclear reactions, which could be clearly observed in the W ± production. For example, the W − production in the future HIC is even suppressed in the p W T < 100 GeV region, though at the LHC it is enhanced by the isospin effect.
These very distinct nuclear modifications at the future HIC compared with those at the LHC result mainly from the smaller momentum fraction x carried by the initial partons at the fixed rapidity or p T . At the LHC, the momentum fraction of the initial-state partons may fall in the anti-shadowing region (for EPS09), however, in the future HIC it may be seen in the kinematic region of the shadowing effect. Furthermore, the calculations of the contribution rates for various partonic subprocesses show that the valence quarks gradually lose their relevance with increasing colliding energy because of their sharply depleted densities in the smaller-x region. With the highly suppressed u and d valence quark distributions, the valence quark density asymmetry is minor, and the isospin effects would become much weaker in the future HIC.
In addition, large differences among results with several parametrization sets of nPDFs could also be seen at the future HIC. For example, results with nCTEQ predict more sup-pression for the vector boson yields than that with EPS09 and DSSZ in future HIC. The theoretical understanding of nuclear modifications of PDFs in nuclear collisions is far from satisfactory, and the massive gauge boson production in the future HIC could provide an excellent probe of the shadowing effect on the nuclear sea quark distribution as well as gluon density by imposing new and precise constraints on the nPDFs at rather small-x.  Panels (c, d) the parton distribution asymmetry ratio r u,d (x) and rū ,d (x). The factorization scale is set at the Z 0 boson mass m Z . The marks at the x-axis are corresponding to the different colliding energies be emphasized. Actually, the parton distribution asymmetry is to a large extent caused by the valence quark (u v , d v ) density asymmetry, as a result the tenuous valence quark density would give a much weaker r ud (or rūd ) in the small-x regime.