Empirical determination of charm quark energy loss and its consequences for azimuthal anisotropy

The vast amount of data collected at the Relativistic Heavy Ion Collider at Brookhaven National Laboratory, along with those recently collected from the collision of lead nuclei at the Large Hadron Collider at CERN have led to the momentous discovery of quark gluon plasma (QGP). A large part of the effort in the coming years will be devoted to determining the precise values of the transport properties of the QGP. In this context, the energy loss suffered by quarks of different flavours and gluons as they traverse the QGP, undergoing collisions and radiating gluons, is the subject of considerable topical interest. It is most simply demonstrated by a suppressed production of hadrons having large transverse momenta in nucleus–nucleus collisions when compared to appropriately scaled productions in pp collisions at the same centre of mass energy per nucleon.

It is often suggested that heavy quarks may lose a smaller amount of energy per unit length during their passage through the QGP compared to light quarks, due to the 'dead cone effect' [1, 2]. The experimental results, however, show similar suppression for light and heavy mesons [3]. A recent calculation by Abir et al [4] incorporating the generalized distribution of gluons in qc → qcg and gc → gcg processes with a proper accounting of the mass of the heavy quarks, leads to a dE/dx for charm quarks which is similar to those for light quarks at energies ⩾10–15 GeV.

The study of charm quark energy loss provides several unique advantages over those of light partons. The charm mesons easily stand out in the multitude of light mesons. Most of the charm quarks are produced in the initial fusion of quarks and anti-quarks ($q\overline{q}\rightarrow c\overline{c}$) and gluons ($gg\rightarrow c\overline{c}$), though a small additional production is expected [5–9] from multiple scattering between jets, jets and thermalized partons, and thermalized partons. It is not yet clearly established if the charm quarks thermalize in the QGP [10], though it is expected that due to their small numbers their impact on the bulk properties of the QGP would be negligible. One also expects that due to their large mass, charm quarks will not change their direction as they traverse the plasma, though they will slow down, making them excellent probes of azimuthal dependence of the conditions of the interacting system. As charm quarks have a large mass, pQCD remains reasonably valid down to the lowest p_T. There is one additional trait which should help us in getting flavour dependence of energy loss of heavy quarks. While a u or a d or a s quark or a gluon can fragment into one of many mesons or baryons, a charm or a bottom quark would mostly fragment into only a D or a B meson or a charm or bottom baryon, respectively. Thus, a charm quark after losing energy will appear as a D-meson having lower energy. We shall see that this would lead to a characteristic enhancement of charm mesons or single electrons at low p_T.

The paper is organized as follows. Section 2 contains expressions for calculations of the nuclear modification factor-R_AA, azimuthal anisotropy coefficient-v₂, charm production, and our model prescription for energy loss. Section 3 contains a discussion of our results, followed by a summary in section 4.

2.1. Charm production

The energy loss of quarks and gluons is most easily seen via the suppressed production of hadrons measured using the nuclear modification factor, R_AA:

$$\begin{equation} R_{AA}(p_{T},y)=\frac{{\rm d}N_{AA}/{\rm d}^{2}p_{T}{\rm d}y}{\langle T_{AA}\rangle {\rm d}\sigma _{pp}/{\rm d}^{2}p_{T}{\rm d}y}, \end{equation} \tag{ 1 }$$

where N_AA is the hadron production for the nucleus–nucleus system at a given impact parameter, T_AA is the corresponding nuclear thickness, and σ_pp is the cross-section for the production of hadrons at the corresponding centre of mass energy/nucleon in pp collisions.

One can now calculate

$$\begin{eqnarray} &&\fl \frac{{\rm d}\sigma _{pp}}{{\rm d}y_1 {\rm d}y_2 {\rm d}^{2}p_{T}} = 2 x_{a}x_{b}\sum _{ij} \bigg[f^{(a)}_{i}(x_{a},Q^{2})f_{j}^{(b)}(x_{b},Q^{2}) \frac{{\rm d}\hat{\sigma }_{ij}(\hat{s},\hat{t},\hat{u})}{{\rm d}\hat{t}} \bigg.\nonumber \\ &&\hphantom{\fl \frac{{\rm d}\sigma _{pp}}{{\rm d}y_1 {\rm d}y_2 {\rm d}^{2}p_{T}} }+ \bigg.f_{j}^{(a)}(x_{a},Q^{2})f_{i}^{(b)}(x_{b},Q^{2}) \frac{{\rm d}\hat{\sigma }_{ij}(\hat{s},\hat{u},\hat{t})}{{\rm d}\hat{t}}\bigg] \bigg/(1+\delta _{ij}) , \end{eqnarray} \tag{ 2 }$$

where p_T and y_{1, 2} are the momenta and rapidities of produced charm and anti-charm and x_a and x_b are the fractions of the momenta carried by the partons from their interacting parent hadrons. These are given by

$$\begin{equation} x_{a}=\frac{M_{T}}{\sqrt{s}}(e^{y_1}+e^{y_2}) ;\quad x_{b}=\frac{M_{T}}{\sqrt{s}}(e^{-y_1}+e^{-y_2}) , \end{equation} \tag{ 3 }$$

where M_T is the transverse mass, $\sqrt{m_{Q}^{2}+p_{T}^{2}}$, of the produced heavy quark. The subscripts i and j denote the interacting partons, and f_i/j are the partonic distribution functions for the nucleons. The fundamental processes included for LO calculations are:

$$\begin{eqnarray} &&g+g \rightarrow c+\overline{c},\nonumber \\ &&q+\bar{q} \rightarrow c+\overline{c} . \end{eqnarray} \tag{ 4 }$$

We note that the above LO pQCD expression reproduces the NLO results [11] when supplemented with a K-factor ≈ 2 (see [12]).

We have used T_AA= 225 fm⁻² for 0–10% centrality for Au + Au collisions at RHIC, as calculated from the Glauber formalism. For Pb + Pb collisions at LHC, T_AA = 195 fm⁻² for 0–20% centrality has been used. We have used the CTEQ5M structure function along with the EKS98 [13] shadowing function. The factorization, renormalization, and fragmentation scales are chosen as $\sqrt{m_Q^2+p_T^2}$ and the charm quark mass has been taken as 1.5 GeV.

2.2. Energy loss

We propose an empirical model for the energy loss for charm quarks which is inspired by a multiple scattering model used earlier by [14] supplemented with a consideration of Baier et al [15] for partonic energy loss [16]. A straightforward empirical implementation of a similar energy loss procedure for light partons has been found to provide a satisfactory dependence of energy loss on the centrality of the collision [17]. We note once again that the nuclear modification of heavy meson production is similar to that for light mesons.

We perform a Monte Carlo implementation of our model calculations and estimate the momentum loss of charm quarks and nuclear modification of D-meson and single electron production. We assume that the energy loss of heavy quarks proceeds via multiple collisions and that the momentum loss per collision is given by, (see for example [18])

$$\begin{equation} (\Delta p)_i=\alpha \, (p_i)^{\beta }, \end{equation} \tag{ 5 }$$

so that one can write

$$\begin{equation} \frac{{\rm d}p}{{\rm d}x}=-\frac{\Delta p}{\lambda }, \end{equation} \tag{ 6 }$$

where α and β are parameters to be determined and λ is the mean free path of the charm quark, taken as 1 fm, in these initial studies. In what follows, we shall consider charm quarks at central rapidities and thus p = p_T. Thus, the momentum of the charm quark after n collisions will be given by

$$\begin{equation} p_{n+1}=p_n-(\Delta p)_n. \end{equation} \tag{ 7 }$$

The charm quark can continue to lose energy in collisions as long as the resulting momentum remains positive. We estimate the probability for the charm quark to have n collisions, while covering the path length L from a Poisson distribution

$$\begin{equation} P(n,L)=\frac{(L/\lambda )^{n}}{n!}{\rm e}^{-L/\lambda }. \end{equation} \tag{ 8 }$$

Taking a value for the coefficient α and the exponent β, we estimate the largest number of collisions—N, which the charm quark having momentum p_T can undergo. Next we sample the number of collisions n, which the charm undergoes from the distribution

$$\begin{equation} p(n)=P(n,L)\bigg/\sum _{n=1}^N P(n,L), \end{equation} \tag{ 9 }$$

to get the final momentum of the charm quark.

Finally, we fragment the charm quarks into D-mesons. Thus we have,

$$\begin{equation} E\frac{{\rm d}^3\sigma }{{\rm d}^3p}=E_Q\frac{{\rm d}^3\sigma (Q)}{{\rm d}^3p_Q}\otimes D(Q\rightarrow H_Q)\otimes F(H_Q\rightarrow e) , \end{equation} \tag{ 10 }$$

where the fragmentation of the heavy quark Q into the heavy-meson H_Q is described by the fragmentation function D. We have assumed that the shape of D(z), where z = p_D/p_c, is identical for all the D-mesons, [19] and

$$\begin{equation} D^{(c)}_D(z)=\frac{n_D}{z[1-1/z-\epsilon _p/(1-z)]^2} , \end{equation} \tag{ 11 }$$

where _p is the Peterson parameter and

$$\begin{equation} \int _0^1 \, {\rm d}z \, D(z)=1 . \end{equation} \tag{ 12 }$$

We have kept it fixed at _p = 0.13.

Here F(H_Q → e) denotes semi-leptonic decay of D-mesons and the electron distribution is taken from [20].

2.3. Azimuthal anisotropy

Non-central collisions of identical nuclei will lead to an oval overlap zone, whose length in and out of the reaction plane would be different. Thus, charm quarks traversing the QGP in and out of the plane will cover different path lengths and lose differing amounts of energy. This would lead to an azimuthal dependence in the distribution of the resulting charm mesons, whose azimuthal anisotropy could be measured in terms of the v₂ coefficient defined by

$$\begin{equation} v_{2}(p_T)=\frac{\int {{\rm d}\phi \frac{{\rm d}N}{p_{T}{\rm d}p_{T}{\rm d}\phi }\cos {(2\phi )}}}{\int {{\rm d}\phi \frac{{\rm d}N}{p_{T}{\rm d}p_{T}{\rm d}\phi }}}. \end{equation} \tag{ 13 }$$

We have approximated the colliding nuclei as having a uniform density with radius R in these calculations and obtained average path length for the charm quarks along a given ϕ using equation (9) of [17].

Let us now discuss our results for the nuclear modification factor and azimuthal anisotropy. We consider three values for the exponent β; 0, 0.5, 1.0 appearing in equation (5), inspired by the three energy loss mechanisms, namely those applicable in the so-called Bethe–Heitler regime, LPM regime, and complete coherence regimes considered by Baier et al [15, 17] which lead to energy loss per unit length as proportional to energy, square root of the energy, and independent of the energy for light partons. Kampfer et al [21] had earlier used this approach to study the effect of charm quark energy loss on the correlated charm decay.

Next we vary α to get a description of the R_AA for single electrons at RHIC (figure 1) and for D-mesons at LHC (figure 2).

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

In figure 1 we show our results for nuclear modification factor R_AA and azimuthal anisotropy v₂ for non-photonic electrons at the top RHIC energy $\sqrt{s}$ = 200 A GeV [22]. Comparing the results of figures 1, we see that the model, assuming momentum loss per collision as proportional to the momentum, closely follows the shape of the experimentally determined R_AA for single electrons almost over the entire range of p_T under consideration. We add that our theoretical calculations have not included the b → e contribution which is contained in the experimental results which can modify the R_AA for larger p_T by up to 10% as the produced b quarks are much less in number and also lose much less energy [23]. The scenario, where $\Delta p \propto \sqrt{p}$, is only moderately successful in describing the data over a limited p_T range of 2–4 GeV/c (figure 1). While the assumption of a constant momentum loss per collision may bracket the R_AA over the very limited range of 2–3 GeV/c, it does not follow the shape of the p_T dependence (figure 1).

The best values of the α determined from the results in figure 1 are used to estimate v₂ for the single electrons. It is evident that our calculations provide a reasonable description of v₂(p_T) for p_T ⩾ 2 GeV/c and overestimate the results for lower p_T. A relaxation of our assumption of a static medium at a constant temperature and a uniform density of the nuclei may improve this agreement.

Similar results are obtained when we apply the model to the R_AA measured [24] for the D-mesons at the LHC (figure 2). We again see that the model using Δp∝p provides a good description of the data over the entire p_T range, while that using $\Delta p \propto \sqrt{p}$ seems to describe the data for p_T ⩾ 4 GeV/c. The constant momentum transfer collision misses the shape of the p_T distribution completely though it is able to bracket the numerical values over a very narrow p_T range of 3–5 GeV/c (figure 2). The predictions for v₂ are given for reference.

Several factors could affect the value of the energy loss coefficient 'α' for a given mechanism. We have verified that increasing (decreasing) the mean free path, λ, by 0.5 fm results in a increase (decrease ) of the coefficient, α_B such that α_B/λ remains unaltered.

We have kept α_s fixed at 0.3 while estimating the initial charm distribution. Taking the renormalization scale as $C\sqrt{p_{T}^{2}+M_{Q}^{2}}$, with C = 1 or 2 leads to a decrease in the value of R_AA by 7-10 %, which can then be offset by decreasing α_B by about 12%.

It is interesting to note that a Focker–Planck equation given by

$$\begin{equation} \frac{\partial {f}}{\partial {t}}=\frac{\partial }{\partial {p_i}}\bigg[A_{i}({\bf p})f +\frac{\partial }{\partial {p_j}}(B_{ij}({\bf p})f)\bigg] , \end{equation} \tag{ 14 }$$

describes the evolution of the distribution 'f', of charm quarks propagating in quark gluon plasma [25] and losing energy due to multiple soft scatterings with light quarks and gluons. This leads to a drag, A_i(p), and a diffusion B_ij(p) on the momentum of the charm quark. Assuming that A_i(p) depends on momentum only, we have

$$\begin{equation} A_{i}({\bf p})=A(p^2)\,p_i , \end{equation} \tag{ 15 }$$

and the energy loss dE/dx can be related to drag coefficient A(p²) by

$$\begin{equation} \frac{{\rm d}E}{{\rm d}x}=-A(p^2)\, p , \end{equation} \tag{ 16 }$$

where E is the energy of the charm quark, and p its momentum. Considering the average temperature of the plasma attained at the RHIC as ≈220 MeV [26], we can read the drag coefficient from figure 3(a) of [27] as ≈0.02 fm⁻¹ for p_T up to 5 GeV/c. We can rewrite the above equation as

$$\begin{equation} \frac{{\rm d}p}{{\rm d}x}=-A(p^2)\,E . \end{equation} \tag{ 17 }$$

Comparing this with one of the ansatzes used for energy loss per collision in this work, namely Δp = α_Bp, we can write that

$$\begin{equation} \frac{{\rm d}p}{{\rm d}x}=-\frac{\alpha _B}{\lambda }\,p, \end{equation} \tag{ 18 }$$

Thus the effective drag coefficient 'A_eff(p²)', can be written as

$$\begin{equation} A_{\rm {eff}}(p^2)=\frac{\alpha _B}{\lambda } \frac{p}{E}, \end{equation} \tag{ 19 }$$

which reduces to α_B/λ for large values of p.

Zoom In Zoom Out Reset image size

We thus note that the effective drag at RHIC energies is about 0.04 fm⁻¹ compared to 0.02 fm⁻¹ estimated for soft multiple collisions by the authors of [27]. Thus we conclude that at RHIC energies only half of the energy loss could be due to collisions, while the other half could be attributed to radiations of gluons. Similarly estimating the average temperature at LHC as about 270 MeV and using the results for drag due to collisions as ≈0.04 fm⁻¹ from [27] at high momentum, we note that the collisions account for only one-third of the energy loss at 2.76 TeV/nucleon.

It is of interest to compare our results with other studies on medium modification of charm propagation reported in the literature. Thus Moore and Teaney [28] have calculated the diffusion 'D' and drag coefficient 'η_D', (denoted by 'A' here) using LO pQCD as

$$\begin{equation} D\approx \frac{6}{2\pi T}\bigg(\frac{0.5}{\alpha _s}\bigg)^2 \end{equation} \tag{ 20 }$$

and

$$\begin{equation} \eta _D=\frac{T}{M_{Q}D}. \end{equation} \tag{ 21 }$$

Taking α_s ≈0.3, this provides η_D ≈0.06 fm⁻¹ at 220 MeV and η_D ≈0.09 fm⁻¹ at 270 MeV, which are larger than our values at RHIC and smaller than those at 2.76 TeV/nucleon at LHC.

Note however, the results of Bass et al [29] where a description of medium modification as well as v₂ for single electrons at RHIC brackets the diffusion coefficient 'D' between 1.5/2πT and 6.0/2πT, when the flow and contributions of bottom electrons is included. This provides a large value for the drag coefficient between 0.17 and 0.68 fm⁻¹. At first these large values may look surprising. However, from equation (20) we see that these would correspond to values of α_s between ≈1.0 and ≈0.5, which are rather large and make the use of perturbative QCD questionable.

In this connection the studies of Gossiaux et al [30], are also of considerable interest as they suggest that the collisional energy loss could be substantially larger if the Debye mass is replaced by hard thermal loop calculation and a running coupling constant is used. We further note the work of the authors of [31], which suggests that a considerable drag could be produced by the resonant heavy quark-light quark interaction, beyond that determined by LO pQCD interactions.

We have already mentioned that c quarks will materialize mostly as D-mesons. This can have an interesting consequence which is already apparent in the figures 1 and 2. The c quarks, after losing energy, will pile up at lower energies and this would result in a characteristic increase in the production of D-mesons as well as single electrons having lower transverse momenta (figure 3). We note that while there is suppression by almost a factor of 2–4 depending on the incident energy (for single electrons) at larger p_T, there is essentially no suppression at lower p_T at RHIC and even an increase by a few percent at LHC where the energy loss is higher and the momentum spectra of c quarks have less steeper slopes. The increase in the case of D-mesons at the LHC, is rather spectacular. This is different from the normal enhancement of mesons having low p_T due to the so-called Cronin effect, which is expected to be less important at higher incident energies and heavy mesons (see e.g., [32]). We note that a similar enhancement in D-meson production at low p_T has been predicted by considering the drag suffered by the heavy quarks in the plasma and during their hadronization [33]. In a forthcoming paper, we shall show that this can lead to a slight increase in the production of low mass dileptons due to correlated charm decay [34].

At the very minimum the present work describes a simple procedure to implement energy loss of heavy quarks in the relativistic collision of heavy nuclei. It will be of interest to explore the energy and centrality dependence of the momentum loss coefficient α.

We have used a simple model where charm quarks traversing QGP lose energy in multiple collisions. Exploring different parametric forms of energy loss we find that the form where the momentum loss per collision is proportional to the momentum, gives a good description of the nuclear modification of single electron production at RHIC and D-meson production at LHC. Comparing our results with drag coefficients due to collisional energy loss calculated earlier, we find that only a part of the energy loss could be due to collisions. We note a characteristic increase in production of single electrons and D-mesons having low momenta due to energy loss of charms having larger energy.

One of us (MY) acknowledges financial support of the Department of Atomic Energy, Government of India during the course of these studies.