LBL-37688 Thermal Equilibration in an Expanding Parton Plasma ∗ 6 9 H. Heiselberg 9 1 NORDITA, Blegdamsvej 17, DK-2100 Copenhagen Ø., Denmark n a Xin-Nian Wang J Nuclear Science Division, Mailstop 70A-3307 1 1 Lawrence Berkeley National Laboratory University of California, Berkeley, CA 94720 USA 1 v 7 4 2 1 0 6 9 / h p Abstract - p e Thermalization in an expanding parton plasma is studied within the frame- h work of Boltzmann equation in the absence of any mean fields. In particular, : v we study the time-dependence of the relaxation time to the lowest order in i X finite temperature QCD and how such time-dependence affects the thermaliza- r tion of an expanding parton plasma. Because of Debye screening and Landau a damping at finite temperature, the relaxation time (or transport rates) is free of infrared divergencies in both longitudinal and transverse interactions. The resultantrelaxation timedecreaseswithtimeinanexpandingplasmalike1/τβ, with β < 1. We prove in this case that thermal equilibrium will eventually be established given a long life-time of the system. However, a fixed momentum cut-off in the calculation of the relaxation time gives rise to a much stronger time dependence which will slow down thermal equilibrium. It is also demon- strated that the “memory effect” of the initial condition affects the approach to thermal equilibrium and the final entropy production. PACS numbers:12.38.Mh,25.75+r,12.38.Bx ∗ ThisworkwassupportedbytheDirector,OfficeofEnergyResearch,DivisionofNuclearPhysics ofthe Office ofHighEnergyandNuclear Physicsofthe U.S. DepartmentofEnergyunder Contract No. DE-AC03-76SF00098. 1 Introduction Perturbative-QCD-based models developed in the last few years predict that nucleus- nucleus collisions at future collider energies are dominated by hard or semihard pro- cesses [1, 2, 3, 4]. These processes happen during the very early stage of the collisions and they produce a rather large number of semihard partons which essentially form a hot and undersaturated parton gas [4, 5]. However, this parton gas is initially far away from thermal and chemical equilibrium [6]. Secondary parton scatterings in the gas may eventually lead to local thermal and chemical equilibrium if the parton interactions are sufficiently strong. Transport calculations based on a semiclassical parton cascade model [4] indicate that thermal equilibrium could be established within a rather short time of about 1 fm/c. However, the complexity of the Monte Carlo simulations makes it difficult to obtain a lucid understanding of the dependence of the thermalization time on the many parameters employed in the model. One such parameter is the cut-off of mo- mentum transfer in binary parton scatterings. The cut-off was first introduced to regularize the infrared divergency of the cross section between two massless partons in high-energy pp and pp¯ collisions [7]. The value of the momentum cut-off is deter- mined phenomenologically to reproduce the measured total cross sections of pp and pp¯ collisions. However, this cut-off is not necessary anymore in a high-temperature quark-gluon plasma, since the Debye screening and Landau damping provide natural regularizations of the infrared divergency. Since transport times depend sensitively on the screening masses which in turn depend on the temperature, the introduc- tion of an artificial cut-off could give rise to a completely different behavior of the thermalization time and consequently the approach to thermal equilibrium. The approach to thermal equilibrium in relativistic heavy ion collisions is dictated by the competition between expansion and parton interactions [8]. If the expansion is much rapid than the typical collision time among partons, e.g., shortly after partons are initially produced, the expansion is closer to free-streaming than hydrodynamic expansion. Only at times in the order of the collision time may the parton gas reach local thermal equilibrium and expand hydrodynamically. Furthermore, the time dependence of the collision time (or the relaxation time) will determine whether the system can eventually reach local thermal equilibrium because of the competition between expansionandpartoninteractions. Ifthecollisiontimeincreases rapidlywith time, the parton system may never thermalize, leading only to a free-streaming limit. The collision time, therefore, is a very important quantity which in turn depends sensitively on the infrared behavior of parton interactions. In QCD, parton scattering cross sections exhibit a quadratic infrared singularity duetotheexchange ofamassless gluon. Theinfraredbehavior canbeimproved byin- cluding corrections from hard thermal loops to the gluon propagators. Resummation of these thermal loops gives rise to an effective gluon propagator which screens long range interactions (Debye screening). Braaten and Pisarski [9] have developed this resummation technique systematically and used it to calculate the damping rate of a soft gluon (p gT) which is gauge invariant and complete to the leading order in the ∼ QCD coupling constant g [10]. For a fast particle (p>T), the exchanged gluons probe ∼ the static limit of the magnetic interactions which by the transversality condition are not screened. One thus has to introduce a nonperturbative magnetic screening mass to regularize the logarithmic infrared singularity in the static limit [11, 12, 13]. As a result, the damping rate for an energetic particle to the leading order in g is Γ T [α ln(1/α )+ (α )] , (1) s s s ∼ O where α = g2/4π. However, as we will argue, damping rates do not determine how s fastasystem approachlocalthermalequilibrium. Whatreallydeterminethethermal- ization processes are the transport rates which are free of the logarithmic divergency after the resummation of thermal loops [14, 15]. This is because thermalization is achieved to the leading order mainly through momentum changes in elastic scatter- ings. Thus, the effective cross section should be weighted by the momentum transfer and the dynamic screening due to the Landau damping of the gluons is sufficient to regularize the logarithmic singularity in the transverse interactions [16, 17, 18]. The resultant transport times for a system near thermal equilibrium behave like 1 Tα2ln(1/α ) (2) τ ∼ s s tr to the leading order in α . s For a system near local thermal equilibrium, the time dependence of the transport times is through the temperature according to Eq. (2). This dependence is in general slower than 1/τ and thus can lead to local thermal equilibrium according to our earlier argument based on the relaxation time approximation [8]. However, if one introducesanartificialcut-offforthemomentumtransfersofelasticpartonscatterings as in the numerical simulation of a classical parton cascade [4], the time dependence will be much stronger. Consequently, as we will demonstrate in this paper, the systemwillapproachlocalthermalequilibriummuchslowly. Wewillalsodemonstrate that inclusion of the screening effects is the key to a slower time dependence of the relaxation time, therefore a faster approach to thermal equilibrium. Thispaperisorganizedasfollows. InSectionII,wefirstre-examinetheBoltzmann equation and the evaluation of the damping rate and the relaxation time to the lowest order, including only 2 2 processes. We will also discuss the time dependence of ↔ the relaxation time in different scenarios. In Section III, we will solve the Boltzmann equation inthe relaxationtime approximation anddemonstrate how timedependence of the relaxation time will affect the approach to thermal equilibrium. We also show how initial conditions of a system affect the thermalization processes and the final total entropy production (or “memory effect”) in Section IV. Finally in Section V we give a summary and an outlook, especially of the numerical simulations of parton thermalization, taking into account of the Debye screening and Landau damping effects without double counting. 2 Time dependence of thermalization time In a system with two-components as, e.g., the quark and gluon plasma, the one that interacts the strongest will thermalize faster than the other. Subsequently there will be momentum and energy transfer between them. Since numerical simulations indicate that the initially produced partons are mostly gluon, we consider here a gluon gas only for simplicity. The Debye screening of color fields in the presence of semihard gluons [19] will also allow us to neglect the effect of mean fields. We furthermore assume that the spatial variation of the system is small on the scale of a collision length so that we can approximate the evolution of the system by the Boltzmann equation [20], 1 v ∂f (p ) = ν dp dp dp F [f] M 2(2π)4δ4(P +P P P ),(3) 1 1 1 2 2 3 4 1234 12→34 1 2 3 4 · − 2| | − − Z F [f] = f f (1 f )(1 f ) f f (1 f )(1 f ), (4) 1234 1 2 3 4 3 4 1 2 ± ± − ± ± where P = ( p ,p ) are the four-momenta of massless partons and dp d3p/(2π)3. i i i i | | ≡ To keep the formula general, are used for bosons (gluons) and fermions (quarks ± and anti-quarks), respectively. The statistical factor ν is 2(N2 1) = 16 for gluons 2 c − and 12N for N flavors of quarks and anti-quarks with N = 3 colors. The squared f f c matrix element, M 2 2/(16E E E E ) is summed over final states 12→34 12→34 1 2 3 4 | | ≡ |M | and averaged over initial states. For gluon-gluon scatterings, su st tu 2 = C 4g4 3 , (5) |M12→34| gg − t2 − u2 − s2 (cid:18) (cid:19) where C = N2/(N2 1) = 9/8 is the color factor of gluon-gluon scatterings, s, gg c c − t and u are the Mandelstam variables. There is clearly a quadratic singularity for small energy ω and momentum q transfers because of the long range interactions mediated by the massless gauge bosons. Because the final state has two identical particles, su/t2 should contribute equally as st/u2 in Eq. (5). We thus can approxi- mate 2 8C g4s2/t2 for small-angle gluon scatterings. Since the collisional 12→34 gg |M | ≈ integral is dominated by contributions from near the singularity, we can assume a small angle scattering approximation, i.e., ω,q E ,E . Then energy-momentum 1 2 ≪ conservation leads to p = p +q, p = p q , 3 1 4 2 − E = E +ω, E = E ω , 3 1 4 2 − ω v1 q v2 q . (6) ≈ · ≈ · The integration over p and p can be rewritten as 3 4 1 q (2π)4Z dp3dp4δ4(P1+P2−P3−P4) = (2π)2 Z d3qZ−qdωδ(ω−v1·q)δ(ω−v2·q) . (7) In a medium, one can use a resummation technique to include an infinite number of loop corrections to the gluon exchange. Using Dyson’s equation, this amounts to an effective gluon propagator. One can use this effective propagator to obtain the effective matrix element squared for forward gluon scatterings (see Appendix A), 1 (1 x2)cosφ 2 M 2 2C g4 − , (8) | gg| ≈ gg (cid:12)q2 +πL(x) − q2(1 x2)+πT(x)(cid:12) (cid:12) − (cid:12) (cid:12) (cid:12) (cid:12) (cid:12) where cosφ = (v qˆ) (v qˆ(cid:12)) and x = ω/q. The scaled self-e(cid:12)nergies in the long 1 1 × · × wavelength limit are given by [21] x 1+x π π (x) = q2 1 ln +i x , (9) L D − 2 1 x 2 (cid:20) (cid:18) − (cid:19) (cid:21) x2 x 1+x π π (x) = q2 + (1 x2)ln i x(1 x2) , (10) T D" 2 4 − (cid:18)1−x(cid:19)− 4 − # where q2 = g2(N + N /2)T2/3 is the Debye screening mass in thermal QCD. The D c f imaginarypartsprovideLandaudampingtopartoninteractionsinathermalmedium. We see that the longitudinal interactions are screened by thermal interactions. How- ever, the transverse interactions still have a logarithmic singularity in the static limit. This singularity can only be regularized by introducing a nonperturbative magnetic screening mass in the calculation of the damping rate of a fast parton. We can use the definition of parton interaction rates, ν 1 Γloss(p1) = (2π2)5 d3p2d3qdωf2(1±f3)(1±f4)2|M12→34|2δ(ω −v1·q)δ(ω−v2·(q1)1,) Z ν 1 Γgain(p1) = (2π2)5 d3p2d3qdωf3f4(1±f2)2|M12→34|2δ(ω −v1·q)δ(ω −v2·q) , (12) Z and rewrite the Boltzmann equation as v ∂f = f Γ (p )+(1 f )Γ (p ) 1 1 1 loss 1 1 gain 1 · − ± ˜ = Γ(p )(f f) , (13) 1 1 − − where Γ(p) = Γ (p) Γ (p) (14) loss gain ∓ ˜ is usually referred to as the damping rate of a particle (or a quasiparticle), and f is defined as Γ (p) ˜ gain f(p) = . (15) Γ (p) Γ (p) loss gain ∓ For a system in local thermal equilibrium, one can relate the damping rate to the imaginary part of the gluon self-energy [22], ImΠµ(p) = 2(p u)Γ(p) , (16) µ − · where the factor 2 comes from our definition of the interaction rates among identical particles. One can also show that, if f(p) takes the local equilibrium form feq(p) = (exp(p u/T) 1)−1, · ∓ Γ (p) loss = ep·u/T, (17) Γ (p) gain usingtheenergyandmomentumconservationandtheidentity1 feq(p) = feq(p)exp(p ± · u/T). Therefore, by definition, f˜(p) becomes feq(p). Thus, the global equilibrium distribution feq(p) is a solution to the Boltzmann equation if the flow velocity u is independent of space and time. We can complete the angular integrations in Eqs. (11) and (12). Making approx- imations f(p ) f(p ) and f(p ) f(p ), we can also complete the integration over 3 1 4 2 ≈ ≈ p by cutting off the integration over q at q 3T [17]. We then obtain the gluon 2 max ≈ damping rate, g4 C ν 1 q2 1 1 (1 x2)2 Γgg = gg gT3 dx maxdq2 + − , (18) 4π 12 Z−1 Z0 (|q2 +πL(x)|2 2|q2(1−x2)+πT(x)|2) including only gluon-gluon scatterings. The contribution from longitudinal interactions is finite and proportional to g2T due to the Debye screening. However, Debye screening is absent in the transverse interactions in the static limit. There is a logarithmic divergency even if Landau damping is taken into account. One solution to this problem is to add a nonpertur- bative magnetic mass m g2T to the transverse self-energy π (x). In this case, mag T ∼ the dominant contribution of the transverse interactions comes from m <q<q mag D ∼ ∼ and thus is proportional to g2T ln(q /m ) which is independent of q q . D mag max D ≫ With π (x) and π (x) given by Eqs. (9) and (10), we can complete the numerical L T integration. A fit to the numerical result gives us g4 C ν T3 q2 m2 q2 q2 Γgg = gg g (ln D 1.0+2.0 mag 0.32 D )+1.1 max , (19) 4π 6 q2 " m2 − q2 − q2 q2 +q2 # D mag D max max D where the first term comes from the transverse interactions while the second from the longitudinal ones. Using the estimate of m 0.255 N /2g2T from Ref. [23] and mag c ≈ neglecting the quark contribution to the Debye screeninqg mass, we have Γgg N α T [ln(1/α ) 0.1+ (α )] . (20) c s s s ≈ − O Note that contributions to the order α2 in Eq. (19) have been neglected, since they s are not complete in our calculation. In order to have a complete calculation of such higher order corrections, one has to include thermal vertex and vacuum corrections which should depend on the renormalization scale. The final result to this order should be invariant under the renormalization group. This result agrees with previous calculations[11,12,13]totheleadingorderofα whichdependsonlyontheimaginary s part of the transverse self-energy, π (x) i(π/4)q2 x at small x. Inclusion of the T ≈ − D full expression of the self-energy only contribute to the next order corrections. The increase in scattering by including quarks is exactly compensated by the increase in Debye screening due to quarks [12]. For a soft gluon (p gT), one can not neglect its thermal mass anymore. The ∼ damping rate for a gluon at rest will not have the logarithmic divergency in the transverse interaction, since the exchanged gluon must carry nonzero momentum and energy at least of order of gT and thus never approach to the static limit. In addition the transverse (magnetic) interactions are reduced by velocity factors which for the massive partons are smaller than the speed of light. The damping rate in this case was found by Braaten and Pisarski [10] to the leading order as Γgg(0) 1.1N α T . (21) c s ≈ Apparently, the damping rate has a nontrivial momentum dependence [24]. As we have mentioned, feq(p u/T) is a solution to the Boltzmann equation as far · as the flow velocity is uniform in space and time. We should emphasize here that the damping rate does not determine how rapidly a system near equilibrium approaches it as one would naively think. The thermalization time is actually related to the transport rates [14, 15, 16, 17, 18]. The easiest way to prove this is to check that the logarithmic divergency that has plagued the calculation of the damping rate of a fast gluon does not appear in the Boltzmann equation. To check this, we make the following expansion: ω2 f(p ) f(p ) ωf′(p )+ f′′(p ) , (22) 3,4 1,2 1,2 1,2 ≈ ± 2 for small angle scatterings. The function F [f] in Boltzmann equation becomes 1234 q2x2 F [f] = [f (1 f )f′′ +f (1 f )f′′ 2(1 f f )f′f′] , (23) 1234 − 2 1 ± 1 2 2 ± 2 1 − ± 1 ± 2 1 2 which is proportional to q2. Here we have dropped terms linear in x since they vanish afterintegrationoverx. Onecanverifythatthisfunctionafterexpansionstillvanishes for the equilibrium distribution, F [feq] = 0. For a system away from equilibrium, 1234 the collisional integral in the Boltzmann equation is nonzero but finite despite the logarithmicsingularityinthetransversepartofthematrixelementsquared, M 2 12→34 | | whenq,ω 0. Becauseofthefactorq2 inF [f],Landaudampingintheself-energy 1234 → of the exchanged gluon is sufficient to give a finite value of the collisional integral. In other words, thermalization not only depends on the parton interaction rates but also on the efficiency of transferring momentum in each interaction. Those interactions with zero energy and momentum transfers do not contribute to the thermalization process, thoughtheir crosssectionsareinfinitely large. This iswhy thethermalization time and other transport coefficients do not suffer from the infrared divergency as pointed out in a number of papers [12, 16, 17, 18]. Forasystem nearequilibrium, onecancharacterize thedeviationfromequilibrium by δf = f feq = θ(p)p∂f/(pu) in a relaxation time approximation. Equivalently, − − · · one has 1 f feq p ∂f = − . (24) p u · − θ(p) · In general the relaxation time θ(p) depends on momentum p and in principle can be obtained by solving the linearized Boltzmann equation. In this paper, we neglect the momentum dependence of the relaxation time. For a pure gluonic gas near local thermal equilibrium where thermalization is achieved through viscous relaxation, the relaxation time is (Appendix B) 1 1.6 0.92N2Tα2ln . (25) θ ≃ c s N α (cid:18) c s(cid:19) Again, because of the extra factor q2 appearing in the transport rate, the Debye screening and Landau damping are sufficient to regularize the effective transport cross section. The dominant contribution comes from interactions with q <q<q , D max ∼ ∼ leading to a logarithmic factor ln(q /q ) as compared to ln(q /m ) in the gluon max D D mag damping rate. Therefore, the dependence of the relaxation time on the (weak) cou- pling constant and the color dimension N is quite different from the gluon damping c rate. If the system is close to thermal equilibrium, the hydrodynamic equations from energy-momentum conservation to the zeroth order of δf can give us the time evolu- tion of the temperature T. For an ideal gluon gas with one-dimensional expansion, T decreases like T/T = (τ /τ)1/3. Therefore, the relaxation time θ increases with time 0 0 with a power of 1/3, (τ/τ )1/3 0 θ = . (26) 0.92N2T α2ln(1.6/N α ) c 0 s c s A more general time dependence of the relaxation time can have a power-law form, θ = θ (τ/τ )β , (27) 0 0 which also covers both the constant (β = 0) and the linear (β = 1) cases as have been studied by several authors [25, 26, 27, 28]. The latter case arises when a constant scattering cross section, σ, is assumed for the relaxation time, θ (σn)−1, and with ∼ a density decreasing as n τ−1 due to one-dimensional expansion. For a system ∼ far away from thermal equilibrium, the time dependence may differ from Eq. (26). Our earlier calculations [8] show that an initially free-streaming system has only a logarithmic time dependence. This is because the phase space for small angle scatterings (q q ) opens up quadratically with time in the free-streaming case and D ∼ it balances the decrease in parton density. We would like to emphasize that the weak time dependence of the relaxation time in Eq. (26) depends very sensitively on the Debye screening of the small-angle parton scatterings which restricts the momentum transfer to q>q = gT. Smaller Debye D ∼ screening mass due to the decrease of the temperature, gives a larger interaction rate which then compensates the decrease of the parton density and thus gives the weak time dependence of the relaxation time. If we use a fixed momentum cut-off q , as cut in most of the numerical simulations of parton production [3] and cascade [4], instead of a time dependent Debye screening mass, the effective transport cross section will remain constant, proportional to N2α2/T2ln(T2/q2 ). The resultant relaxation time c s cut 1 N2α2T ln(T/q ) (28) θ ∝ c s cut will increase more rapidly with time. If we have to include transverse expansion later in the evolution of a system, then the temperature will decrease faster, like τ /τ, than in the one-dimensional case 0 assuming hydrodynamic expansion. The relaxation time even with the inclusion of Debye screening will increase linearly with time. The relaxation time with a momen- tum cut-off applied to parton interactions will increase faster than linear with time, which will only lead the system into free-streaming. At this point we should emphasize that we have only considered the lowest or- der contribution from 2 2 processes in our calculation of the relaxation time. In ↔ principle, higher order processes, like 2 2+n, should also contribute to the ther- ↔ malization. Such processes can be included by considering high order thermal vertex corrections. For a complete calculation, one should also include vacuum corrections and the result should depend on the renormalization scale and obey the renormaliza- tion group equation. In general, contributions from 2 2+n processes should have ↔ a form [29], Γ Γ [α ln(q2/q2)]n, (29) n ∼ 0 s 0 where q2 is the momentum scale of these processes. At zero temperature, q is some 0 confinement scale below which perturbative QCD is no longer applicable. At finite temperature, q is very likely to be replaced by screening masses. Since the largest 0 momentum scale in a system at finite temperature is q2 T2, the leading correction ∼ from 2 2+n processes must be, ↔ Γ Γ [α ln(1/α )]n. (30) n 0 s s ∼ Such corrections thereforearehighorders inα ln(1/α ) andarenegligible intheweek s s couple limit. For temperatures not far above the QCD phase transition temperature T 200 MeV, the strong coupling constant is not very small. The above contribu- c ∼ tions might not be negligible. However, for an order-of-magnitude estimate, we can neglect these higher order contributions. If one considers the chemical equilibration of a kinetically thermalized system as in Ref. [6], n m multiplication processes ↔ become very important. The leading contribution to the chemical equilibration, in this case, comes from 2 3 processes. ↔ 3 Approach to thermal equilibrium Letusconsider theearlystageofavery heavyioncollisionwheretransverseexpansion is not important yet. We then can treat the system as a one-dimensional system. We assume along with Bjorken [30] a scaling flow velocity x u = µ = (coshη,sinhη,0 ) (31) µ ⊥ τ in the longitudinal direction, where 1 t+z τ = √t2 z2 , η = ln (32) − 2 t z (cid:18) − (cid:19) are the proper time and spatial rapidity, respectively. In terms of these new vari- ables, the Boltzmann equation in the relaxation time approximation for a system near thermal equilibrium becomes ∂f tanhξ∂f f feq = − , (33) ∂τ − τ ∂ξ − θ where ξ = η y and y is the rapidity of a particle, − 1 E +p z y = ln . (34) 2 E pz! − Since p u = p coshξ, we can see that the solution to the above Boltzmann equation T · is a function of ξ, τ and p , and is Lorentz invariant under longitudinal boost. One T may then, as done in [25], solve the Boltzmann equation in the central slice only, i.e., η = 0. To simplify the problem, we assume the system is already in chemical equilibrium so that the gluon chemical potential vanishes. This is not always true in a realistic situation as shown by numerical simulations [5] of initial parton production at around RHIC energies. However, at LHC energies, the small-x behavior of the partondistributions as measured by recent HERA experiments [31] gives much higher densities of initially produced partons very close to chemical equilibrium [32]. For a given momentum and time dependence of the relaxation time, one can find the solution to Eq. (33) in an integral form, τ sinhξ f(p ,ξ,τ) = eχ0−χf (p ,sinh−1( )) T 0 T τ 0 χ τ sinhξ + dχ′eχ−χ′feq(p ,sinh−1( ),T′) . (35) T τ′ Zχ0