# Heavy Quarks in Turbulent QCD Plasmas

###### Abstract

The quark-gluon plasma, which is produced at an early stage of ultrarelativistic heavy-ion collisions, is expected to be initially strongly populated with chromodynamic fields. We address the question how heavy quarks interact with such a turbulent plasma in comparison with an equilibrated one of the same energy density. For this purpose we derive a Fokker-Planck transport equation of heavy quarks embedded in a plasma of light quarks and gluons. We first discuss the equilibrium plasma and then the turbulent one applying the same approach, where the heavy quarks interact not with the plasma constituents but rather with the long wavelength classical fields. We first consider the three schematic models of isotropic trubulent plasma and then the simplified model of glasma with the chromodynamic fields only along the beam direction. The momentum broadening and collisional energy loss of a test heavy quark are computed and compared to those of equilibrium plasma of the same energy density.

## I Introduction

The early stage of relativistic heavy-ion collisions is least known because there are hardly any experimentally accessible signals of the phase. Nevertheless one expects that the quark-gluon plasma, which is produced in the collisions, is initially strongly populated with chromodynamic fields. Within the framework of the Color Glass Condensate (CGC) approach, see e.g. the review Gelis:2012ri (), color charges of partons confined in the colliding nuclei act as sources of long wavelength chromodynamic fields which can be treated classically because of large occupation numbers of the soft modes. Since the density of color charges per transverse area of heavy nuclei is large, the corresponding momentum scale is expected to be significantly bigger than the QCD scale parameter . Consequently, the coupling constant is presumably sufficiently small and perturbative methods are applicable. The system, however, is rather strongly interacting because of the high-amplitude fields present in the system.

A momentum anisotropy of the early stage quark-gluon plasma makes it unstable with respect to chromomagnetic modes which in turn cause a spontaneous generation of the fields, as explained at length in the review article Mrowczynski:2016etf (). Therefore, the effect of strong fields is further enhanced. Following the terminology of electromagnetic plasma, we call such a nonequilibrium system of fields as the turbulent plasma meaning that numerous modes are excited in the system. In the CGC approach the non-equilibrium system of fields from the early stage of relativistic heavy-ion collisions is called glasma Gelis:2012ri () and it can be treated as a specific realization of the turbulent QCD plasma. Leaving aside the mechanism of field generation and its detailed structure, one asks what are the transport properties of turbulent plasmas. We are specifically interested how heavy quarks - charm or beauty - behave in such a system when compared to the equilibrium plasma of the same energy density.

Heavy quarks are often treated as a probe of strongly interacting matter created in relativistic heavy-ion collisions, see e.g. the review Prino:2016cni (). Thanks to their large masses the quarks are produced only at the earliest stage of the collision due to hard interactions of partons from incoming nuclei. Later on they propagate through a surrounding medium testing the entire history of the system. It has been long believed that the interaction of heavy quarks is significantly weaker than that of light quarks or gluons but experimental data clearly contradict the expectation. As discussed in the review Prino:2016cni (), the behavior of mesons containing a heavy quark is rather similar to that of light mesons at both small and large transverse momenta. The problem is not fully resolved.

The medium created in relativistic heavy-ion collisions evolves fast towards the locally equilibrated quark-gluon plasma which expands hydrodynamically and ultimately is converted into a hadron gas. Final momentum spectra of heavy quarks are mostly shaped in the long-lasting equilibrium phase which is relatively well understood. An effect of a pre-equilibrium phase is often entirely ignored but this transient phase can significantly influence heavy-quark spectra because of its high density. Non-equilibrium calculations recently performed in a framework of kinetic theory Das:2017dsh () confirm the suggestion. However, we are interested even in the earlier phase when the medium is not described in terms of quasi-particles, as in a kinetic theory, but rather as a system dominated by classical fields which is the turbulent plasma.

A simple parametric estimate suggests that the interaction of heavy quarks in a turbulent plasma is much stronger than in the equilibrium one, if the coupling constant is small. The momentum broadening parameter , for example, is of order in equilibrium plasmas. Since the quark of interest actually interacts with soft gluons emitted by plasma constituents, one can think that the factor is composed of two pieces of . The first one is related to the gluon emission and the second one to the gluon absorption. If soft chromodynamic fields are present in the plasma, the interaction of the quark should be rather of order than . In Sec. VI we argue that is indeed not of the order , not even but presumably of the order in a turbulent plasma.

Because of their big masses, relaxation times of heavy quarks, which are produced in relativistic heavy-ion collisions, are expected to be significantly longer than that of light quarks and gluons. When an equilibrium or, more generally, a stationary state is reached by light quarks and gluons, heavy quarks need some extra time to adjust to the state of the plasma. Such a situation is naturally described in terms of the Fokker-Planck transport equation which was indeed repeatedly applied to heavy quarks in Moore:2004tg (); Svetitsky:1987gq (); vanHees:2004gq (); Mustafa:2004dr (). The equation is usually derived from the Boltzmann equation by applying the so-called diffusion approximation to the collision term LP81 (). The approximation assumes that the momentum transfer to the heavy quark in every collision is much smaller than the quark momentum.

The aim of this paper is threefold. In the first part we rederive the Fokker-Planck equation of heavy quarks which do not interact with plasma constituents but rather with soft classical fields present in the plasma. Specifically, we apply the so-called quasi-linear theory known from the electromagnetic plasma Ved61 (); LP81 (). The theory assumes that the distribution function can be decomposed into a large but slowly varying regular part and a small fluctuating or turbulent one which oscillates fast. The average over a statistical ensemble of the turbulent part is assumed to vanish and thus the average of the distribution function equals its regular part. The turbulent contribution to the distribution function obeys the collisionless transport equation while the transport equation of the regular part is determined by the fluctuation spectra which provide the collision term. The derivation presented here closely follows the procedure which was developed for QCD in Mrowczynski:2009gf (), where, however, only the longitudinal chromoelectric field was taken into account and here the complete chromodynamic field is considered. The equilibrium correlation functions of chromodynamic fields, which are needed to obtain the quasi-linear transport equations, were derived in Mrowczynski:2008ae ().

Our second aim is to confront the equilibrium plasma with the turbulent one. Therefore, we consider three models of isotropic turbulent plasma in the second part of the paper. Postulating a form of the correlation functions of chromodynamic fields, we derive the coefficients of the Fokker-Planck equation which can be related to the energy loss, momentum broadening and diffusion coefficient of heavy quarks in the plasma. The transport coefficients of turbulent plasma are compared to those of the equilibrium one at the same energy density.

The third aim is to study an evolution of heavy quarks at the earliest stage of relativistic have-ion collisions. Since the chromoelectric and chromomagnetic fields spanned between the receding nuclei are initially mostly parallel to the beam direction, we model the glasma with the boost invariant correlation functions of longitudinal fields. The energy loss and momentum broadening of heavy quarks are computed, assuming that all energy of the glasma is accumulated in the longitudinal fields.

At the end of the introductory remarks we note that an approach similar to ours, which was also inspired by the electromagnetic plasma studies Dupree:1966 (), was formulated in Asakawa:2006jn (), see Majumder:2007zh (); Chandra:2015gma () as well. We also mention an attempt Song:2015jmn () to study transport of heavy quarks in a plasma populated by strong chromodynamic fields. Unfortunately, the paper is flawed as the framework of an isotropic Langevin approach is applied to anisotropic plasmas.

Throughout the paper we use the natural system of units with ; our choice of the signature of the metric tensor is . Lorentz indices are denoted with and label the Cartesian coordinates . The color indices of the adjoint representation of gauge group are .

## Ii Derivation of Fokker-Planck equation

Our derivation of the Fokker-Planck equation of heavy quarks embedded in quark-gluon plasma starts with the transport equation of the Vlasov form

(1) |

where the distribution function of heavy quarks is the hermitian matrix which belongs to the fundamental representation of the SU() group. The distribution function depends on the time (), position () and momentum () variables. There is no explicit dependence on the time-like component of the four-momentum as the distribution function is assumed to be non-zero only for momenta obeying the mass-shell constraint that is . The quark velocity equals and with being the chromodynamic potential in the fundamental representation. The mean-field term of the transport equation (1) is expressed through the color Lorentz force with the chromoelectric and chromomagnetic fields also belonging to the fundamental representation. The symbol denotes the anticommutator. The derivation of the transport equation (1) is discussed in detail in the review Mrowczynski:2016etf ().

Further on we assume that the chromodynamic fields and the distribution function which enter the transport equation (1) can be decomposed into a regular and fluctuating or turbulent component. The distribution function is thus written down as

(2) |

where denotes ensemble average; is called the regular part while is called the fluctuating or turbulent one. It directly follows from Eq. (2) that . The regular contribution is assumed to be color neutral or white, and it is expressed as

(3) |

where is the unit matrix in color space. Since the distribution function transforms under gauge transformations as , where is the transformation matrix, the regular contribution of the form (3) is gauge independent. We also assume that

(4) |

but at the same time

(5) |

What concerns the chromodynamic fields, we assume in accordance with Eq. (3) that their regular parts vanish and thus

(6) |

We substitute the distribution function (2) into the transport equation (1) and linearize the equations in the fluctuating contributions. Thus we get the equation

(7) |

where is the substantial or material derivative.

Now we substitute the distribution functions (2) into the transport equations (1) but instead of linearizing the equation in the fluctuating contributions, we take the ensemble average of the resulting equation and trace over the color indices. Thus we get

(8) |

Since the regular part of distribution function is assumed to be color neutral, see Eq. (3), the term vanishes because the fields are traceless. The trace over color indices also cancels the terms originating from covariant derivatives like . We finally note that the trace is gauge independent as the regular distribution function is.

Now, we are going to write down the transport equation (8) in the Fokker-Planck form. For this purpose we observe that due to the condition (5), the space-time dependence of the regular distribution function can be neglected in the linearized transport equation (7) and then, the equation becomes easily solvable. We solve it with the initial condition

(9) |

using the one-sided Fourier transformation defined as

(10) |

The inverse transformation is

(11) |

where the real parameter is chosen in such a way that the integral over is taken along a straight line in the complex plane, parallel to the real axis, above all singularities of .

The linearized transport equation (7), which is converted into the algebraic equation by means of the one-sided Fourier transformation, is solved as

(12) |

We stress that although we have ignored the (weak) frequency and wave number dependence of the regular distribution , the fields retain their full frequency and wave number dependence in the expression (12). Inverting the one-sided Fourier transformation, one finds the solution of the linearized transport equation as

(13) |

where we assumed that and are analytic functions of .

With the help of the solution (13), the force term in the transport equation (8) becomes

(14) |

The second term in the r.h.s. of Eq. (14) can be manipulated to the form

(15) |

which is effectively the definition of the vector . We also introduce the tensor

(16) |

and we note that, as explained in the subsequent sections, and become time independent for a sufficiently long . Then, the transport equation (8) can be written as the Fokker-Planck equation

(17) |

Since the distribution function carries no information about color degrees of freedom, the function is gauge invariant, and consequently and should be gauge invariant as well. However, one observes that and as defined by Eqs. (15) and (16) are gauge dependent because the traces are of nonlocal quantities in the definitions (15) and (16). The starting transport equation (1) is gauge covariant but the linearization procedure breaks the covariance because the covariant derivative is replaced by the normal one. Consequently, the solution (13) is not gauge covariant – the right-hand side of Eq. (13) transforms differently under local gauge transformations than the left-hand side. To cure the problem, one modifies the solution (13) by means of the link operator which is also called the gauge parallel transporter, see e.g. Sec. IIIE of the review article Mrowczynski:2016etf (). Then, the modified solution obeys Eq. (8) with the covariant derivative instead of the normal one. Let us briefly discuss the procedure in a context of the Fokker-Planck equation (17).

According to Eq. (16), the tensor is determined by the traces of the field correlation functions like where chromodynamic fields are written in the adjoint representation of the group which is used further on. The trace becomes gauge invariant under the replacement

(18) |

where is the link operator defined as

(19) |

Here is the adjoint representation generator of the group and denotes the ordering along the path connecting the points and . Since the fields transform as vectors under the local gauge transformation and the link transforms as

(20) |

one checks that the trace of the correlation function which includes the link is indeed gauge invariant. Consequently, the tensor is gauge invariant. Analogously one achieves the gauge invariance of the vector . Further on, whenever any nonlocal correlation function shows up a presence of the link operator is implicitly assumed even so it is not explicitly written.

## Iii Fokker-Planck equation

Although this is a textbook material we briefly discuss here the Fokker-Planck equation (17). We first note that in the isotropic plasma the tensor and vector both depend on a single vector that is the heavy-quark velocity . Therefore, they can be written as

(21) | |||||

(22) |

where and the coefficients and are equal to

(23) |

The equilibrium distribution function of the form

(24) |

with being the temperature of the plasma of light quarks and gluons, where heavy quarks are embedded, is expected to solve the transport equation (17). This is indeed the case if the coefficients and obey the condition

(25) |

which in the isotropic plasma reads

(26) |

When the plasma is isotropic and the coefficients and are equal to each other and independent of , the Fokker-Planck equation reads

(27) |

where .

The quantities and have a clear physical meaning. As discussed in e.g. the classical monograph Kampen:1987 (), the average momentum change per unit time and the correlation of momentum changes per unit time are given as

(28) | |||||

(29) |

Using the formulas (28) and (29), and can be related to the collisional energy loss and transverse momentum broadening of a heavy-quark in the quark-gluon plasma, which play an important role in a theoretical description of the jet quenching phenomenon. The parameter controls the radiative energy loss in a plasma medium Baier:1996sk (). One easily finds that

(30) |

Since , the energy loss per unit path equals

(31) |

which in isotropic plasmas reads

(32) |

The coefficient , which is the broadening per unit path of the distribution of the test parton’s momentum transverse to the initial parton’s momentum, is immediately found as

(33) |

and in an isotropic plasma it equals

(34) |

When we deal with an equilibrium plasma and the coefficients are equal to each other and independent of , the Fokker-Planck equation can be related to the nonrelativistic Langevin equation Kampen:1987 (). Then, the diffusion constant can be expressed as

(35) |

We are mostly interested in turbulent QCD plasmas populated with strong chromodynamic fields but we start with the equilibrium system where the fields are only at a level of thermal noise. We rederive the known Fokker-Planck equation Svetitsky:1987gq () in a different way to demonstrate the reliability of our approach.

## Iv Equilibrium plasma

We assume that the quark-gluon plasma, in which heavy quarks are embedded, is in thermodynamical equilibrium and we first derive in this section explicit expressions for the coefficients , and which enter the Fokker-Planck equation.

### iv.1 Computation of and

As the formula (16) shows, the quantity is given by the correlations functions , , , and which were studied in detail in Mrowczynski:2008ae (). The explicit expressions are collected in Appendix A. Since the correlation functions are of the structure (109), the tensor is written as

(36) |

where the chromodynamic fields are expressed in the adjoint representation of the group and is the fluctuation spectrum. For a translationally invariant system it is defined as

(37) |

A more general definition is discussed in Appendix B. Combining the equilibrium fluctuation spectra (110), (111) and (112), one finds

(38) | |||

where and are chromodielectric functions which for the equilibrium plasma of massless particles are also given in Appendix A.

After performing the elementary time integration in Eq. (36), one is left with the integral over the four-vector . Taking into account only the terms of the integrand which are even as a function of and give nonzero contributions, one obtains

(39) |

where . In the limit , we have

(40) |

and thus we get

(41) |

One can show that the expression (41) properly approximates the formula (39) if the spectrum weakly changes as a function of in the interval . When the time grows the condition is easier and easier to fulfill.

Since the plasma under consideration is isotropic, the tensor is fully determined by the two functions and introduced in Eq. (21). Substituting the explicit form of the fluctuation spectrum (38) into Eq. (41), one finds the coefficients and as

(42) | |||||

where is the Casimir invariant.

The coefficient , which is determined by the correlations functions , and can be derived directly from the formula (15). Such a derivation for a simplified case of only longitudinal electric field present in the plasma can be found in Mrowczynski:2009gf () where it is shown that the condition (25) or (26) is indeed satisfied. Here instead we refer to the condition (25) to obtain .

Once the coefficients is given by Eq. (21) together with Eqs. (42), (IV.1) and by Eqs. (22) and (26), the Fokker-Planck equation (17) is fully specified. We note that and depend on heavy quark momentum and its mass only through the velocity . Therefore, and become independent of when the heavy quarks of interest are truly relativistic and . Although, the coefficients and are independent of the quark mass, the Fokker-Planck equation (17) does depend on which is evident when the momentum derivatives are replaced by the velocity derivatives.

Since the coefficients and depend on the quark mass only through the velocity, one might think that the corresponding Fokker-Planck equation is applicable to quarks of any mass. However, it is not true. As mentioned in the introduction, the typical momentum transfer in a single collision, which is of order , must be much smaller than the quark momentum that is . And here the quark mass matters.

### iv.2 Limit of small

In the limit of small velocities of heavy quarks, the equilibrium Fokker-Planck equation gets a simpler form and the coefficients can be estimated analytically. Indeed, when heavy quarks are nonrelativistic (), we have and one can use the approximate formulas (116) and (117). Then, Eqs. (42) and (IV.1) read

(44) | |||||

(45) |

where the contributions originating from appear to vanish. Introducing spherical coordinates with the axis along the vector , the integrals (44) and (45) are rewritten as

(46) | |||||

(47) |

where the trivial azimuthal integrals are performed and with being the angle between and .

When , the integrals (46) and (47) can be estimated as follows. One first observes that the dominant contribution comes from . Assuming that , the integrals over are easily computed and one obtains

(48) |

Approximating the integrand as , we finally get

(49) |

As one sees in Eq. (48) or (49), and are independent of and equal to each other. The formula (35) is therefore applicable and the inverse diffusion constant equals

(50) |

which agrees with Eq. (12) of the study Moore:2004tg () when only the logarithmic term is taken into account.

### iv.3 Numerical results

Here we show some numerical results for the equilibrium QGP of and . The Debye mass is computed according to the formula (115). Fig. 2 presents the coefficients and , which are obtained directly from Eqs. (42) and (IV.1), as functions of the velocity . The coupling constant and the temperature are and . In Fig. 2 we show how and depend on the temperature . The coupling constant is again and the velocity equals .

We have checked that the values of and , which we obtained numerically, agree rather well with those computed by Svetitsky Svetitsky:1987gq () except in the domain of small velocities . The agreement is not trivial because the coefficients of the Fokker-Planck equation were derived in Svetitsky:1987gq () from the matrix elements of heavy quark binary interactions with plasma constituents. To remove infrared divergences of the matrix elements, a mass parameter, corresponding to the Debye mass, was included in the gluon propagator. Since the procedure is not very accurate, it presumably explains the difference with our and in the domain of small velocities. On the other hand our approach does not treat properly the interactions with a momentum transfer exceeding the Debye mass. Nevertheless the results of both approaches are numerically rather similar.

## V Turbulent QGP

In this section we consider a Fokker-Planck equation of heavy quarks in a turbulent QGP which is populated with strong chromodynamic fields. The plasma is assumed to be isotropic and translationally invariant both in time and space. The tensor , which enters the Fokker-Planck equation (17), is given by Eq. (41). The method to obtain , which is used in Mrowczynski:2009gf (), works only for equilibrium plasmas. Therefore, we will refer to the relation (25). However, it implicitly assumes that in the long time limit the system of heavy quarks described by the Fokker-Planck equation reaches a state of thermal equilibrium with temperature . First of all, the value of is, in principle, unknown and one needs additional arguments to determine it. There is also a more important problem - it is unclear what are the properties of for which the assumption of equilibrium makes sense. If, for example, there are only magnetic fields in the plasma, is purely transverse and . Consequently, the Fokker-Planck equation reads

(51) |

and any stationary isotropic function solves the equation because . Therefore, pure magnetic fields do not drive the system to the thermal equilibrium, as expected. In spite of these concerns, we will use the relation (25) to get .

### v.1 Gaussian correlation functions of independent and fields

We start with a simple model proposed in Asakawa:2006jn () where the correlation functions of electric and magnetic fields are chosen in the following Gaussian form

(52) | |||||

(53) | |||||

(54) |

The correlation lengths and the parameters of dimension mass to the fourth power will be discussed later on. The fluctuation spectrum obtained from the correlation function (52) is also Gaussian and it equals

(55) |

Substituting the correlation functions (52), (53), and (54) into Eq. (41), the tensor is found as

(56) |

which gives

(57) | |||||

(58) |

When and , the coefficients and are, as in the equilibrium case, equal to each other and independent of , and

(59) |

The Fokker-Planck equation is then of the form (27).

### v.2 Gaussian correlation function of vector potentials

Since the electric and magnetic fields are, in general, coupled to each other, the functions , , , and are not fully independent from each other. The electric and magnetic fields are automatically related to each other if one postulates the correlation function of the four-potential and then computes the correlation functions of the and fields. In this section we follow this path. Specifically, we assume the Gaussian correlation function of the vector potential

(60) |

The parameter of the dimension mass squared will be discussed later on. The fluctuation spectrum of the potential equals

(61) |

We further choose the radiation gauge

(62) |

and the electric and magnetic fields are obtained in the linear regime as

(63) |