Simulating plasma instabilities in SU(3) gauge theory

We compute nonequilibrium dynamics of plasma instabilities in classical-statistical lattice gauge theory in 3+1 dimensions. The simulations are done for the first time for the SU(3) gauge group relevant for quantum chromodynamics. We find a qualitatively similar behavior as compared to earlier investigations in SU(2) gauge theory. The characteristic growth rates are about 25 % lower for given energy density, such that the isotropization process is slower. Measured in units of the characteristic screening mass, the primary growth rate is independent of the number of colors.


Introduction
The theoretical understanding of the apparent fast thermalization in collision experiments of heavy nuclei at the Relativistic Heavy Ion Collider provides a challenge for theory. It was noticed that plasma instabilities in the locally anisotropic medium might play an important role for a rapid isotropization of the equation of state, which is relevant to explain the observed hydrodynamic behavior [1,2,3,4,5]. Extensive studies have been carried out for SU (2) pure gauge theory using the hard-loop effective theory of soft excitations, which is based on collisionless kinetic theory for hard particles coupled to a soft classical field. 1 This approach neglects quantum corrections and may also be considered as an approximation of the classical-statistical field theory limit of the respective quantum gauge theory [3,4]. Classical-statistical lattice gauge theory provides a quantitative description in the presence of sufficiently large energy density or occupation numbers per mode. The simulations are done by numerical integration of the classical lattice equations of motion and Monte Carlo sampling of initial conditions.
In this work we present for the first time classical-statistical lattice gauge theory simulations for the SU(3) gauge group in 3+1 dimensions, relevant for quantum chromodynamics (QCD). Simulations using the hard-loop approximation for SU (3) gauge group in 1+1 dimensions have been done in Ref. [5]. Though the dynamics in one and three spatial dimensions are known to be very different, it is not expected that SU(3) results will qualitatively change as compared to those previously obtained for SU (2), however a quantitative estimate of the involved time scales seems imperative. This work is a follow-up to our previous studies concerning the numerically less demanding gauge group SU(2) [4]. In this Letter we only describe the relevant changes and refer to that reference for further computational details. We find a qualitatively similar behavior as compared to the earlier investigations in SU(2) gauge theory, but for given initial energy density the characteristic growth rates are about 25% lower such that the isotropization process is slower.
The paper is organized as follows. In Section 2, we outline the setup of our calculations, and the algorithm for solving the equations of motions for the SU(3) group. In Section 3 we present the results of the calculations and explain the dependence of growth rates on the gauge group in terms of diagrammatics. We conclude in Section 4.

Classical-statistical gauge field theory on a lattice
Following Ref. [4] we use the Wilsonian lattice action on a 3+1 dimensional Minkowskian lattice, Tr Tr written in terms of plaquette variables where U † x,νµ = U x,µν . Here U x,µ denotes the link variable which is the parallel transporter associated with the link from the neighboring lattice point x +μ to the point x ≡ (t, x) in the direction of the lattice axis µ = 0, 1, 2, 3. The lattice parameters are defined as where γ ≡ a s /a t is the ratio of the spatial and temporal lattice spacings, and we will consider g 0 = g s = g as the coupling constant of the lattice theory.
Varying the action (1) w. r. t. the spatial link variables U x, j yields the leapfrogtype equations of motion where λ b are the generators of SU(N) satisfying Tr(λ a λ b ) = 1 2 δ ab . In the case of SU(2) these can be related to the Pauli matrices σ a through λ a = 1 2 σ a . For SU(3), the generators are given by the Gell-Mann matrices divided by two. The electric fields in Eq. (4) are Varying the action (1) w. r. t. to a temporal link gives the Gauss constraint Using the gauge freedom the dynamics is computed in temporal axial gauge, which is defined by the gauge condition U x,0 = 1. This choice still leaves the possibility of a time-independent gauge fixing which we do not specify. The algorithm for the numerical simulation may be summarized as follows: 2 , j is determined at the next time step from the defining equation of the temporal plaquette in temporal axial gauge: The main difference between the algorithms for SU(2) and SU(3) occurs in step 2. For U ∈ SU(2) and given ImTr(Uσ a ) = 2b a , with a = 1, 2, 3, one can represent U as U = √ 1 − b a b a 1 + ib a σ a for the solution which is close to the unit matrix. For SU(3) (and SU(N) with N ≥ 3 in general) there is no simple corresponding procedure in step 2 and we solve for the temporal plaquette numerically. For a general complex 3 × 3 matrix, which we use in the program, there are 18 real variables (3 × 3 complex numbers) and 28 real, but not independent, equations (8 from (5), 2 from detU x,0 j = 1, and 18 from U x,0 j U † x,0 j = 1). Since one only needs 18 equations, we use part of the unitarity equation: components 12, 13, and 23, and the real part of components 11 and 22 because the imaginary parts do not give any constraints on the variables. These equations are solved with the multidimensional Newton method [8]. The fastest simulation speed was achieved by using the unit matrix as a starting point of the Newton method, which typically leads to convergence in about 3-4 steps. Since for every new link variable one has to iterate the Newton method in 18 dimensions a few times, the SU(3) simulation is considerably slower than the corresponding SU(2) simulation.
The initial Gaussian probability functional is chosen as in Ref. [4] such that where the gauge fields are calculated from the link variables using We consider E a j = −Ȧ a j = 0 at t = 0 fulfilling the Gauss constraint (6). We typically choose ∆ z ≪ ∆, and the distribution is practically δ(p z )-like on the lattice. Here ∆ determines the typical transverse momentum of the gluons and may be associated with the saturation scale Q s at time Q −1 s in the saturation scenario [7]. We solve the equations of motion with the aforementioned algorithm for a set of initial configurations sampled according to Eq. (7) and compute expectation values as averages over the results from the individual runs.
The local energy density in lattice units is determined by the action aŝ and we denote its average value byε ≡ ε(t, x) . We choose the factorÃ appearing in Eq. (7) such thatε = 0.05. For the conversion to physical units we fix the lattice spacing a s from the relation between the physical average energy density and its lattice analogue according to ε =ε · a −4 s . If g is taken to be different from one then the lattice spacing is altered by a factor 1/ √ g , which follows from Eq. (9). The values for a s will not be altered significantly as long as g ∼ O(1). For later reference, we note that the relation betweenε, N andÃ is approximately given byε for our initial conditions, where the magnetic fields are calculated from the spatial plaquettes according to similarly to the electric fields (5).   The plotted low-momentum modes clearly show exponential growth starting at the very beginning of the simulation. In contrast to these "primary" instabilities occurring at low momenta, one observes from Fig. 1 that gauge field modes at sufficiently high momenta do not grow initially. The higher wave number modes typically exhibit exponential growth at a "secondary" stage that sets in later, but with a significantly larger growth rate. The secondaries arise from fluctuation effects induced by the growth in the lower momentum modes, which can be explained by taking into account (2PI) resummed loop diagrams beyond the hard-loop approximation as discussed in detail for SU(2) gauge theory in Ref. [4]. In that case it has been shown that the exponential growth saturates when all loop diagrams become of order one, which leads to a subsequent slow evolution towards a power-law regime. 3 Since the fastest growing mode in SU(2) gauge theory grows earlier to sizes where loop effects for higher modes become important, one observes from Fig. 1 that the secondary growths also start earlier compared to the SU(3) case. Fig. 2 displays the momentum dependence of the growth rates for |A(t, p)| 2 Figure 3: The tadpole diagram whose non-vacuum part determines m T as described in the main text.

Results
obtained from a fit to an exponential, which is done separately for the primary and secondary growth rates. One observes that while the primary rates are approximately 25 % bigger for SU(2), the secondary rates can be similar within the given errors. The dependence of the primary growth rates on the gauge group can be understood from a diagrammatic analysis as follows. As outlined in Section 2, the overall gauge field amplitudeÃ defined in (7) is chosen such that the average lattice energy densityε has a prescribed value. For simulations of different gauge groups SU(N 1 ) and SU(N 2 ) at the sameε using initial conditions sampled from (7) with identical ∆ and ∆ z , Eq. (10) implies that the respective overall gauge field amplitudesÃ are related according to Plugging in the relevant numbers we find for SU(2) and SU(3) In Ref. [4] we have analyzed for SU(2) gauge theory relevant diagrammatic contributions to the nonequilibrium evolution at various characteristic time scales, which holds along the same lines for SU (3). As a characteristic self-energy contribution relevant for the primary growth rates at early times, we consider the non-vacuum part m 2 T of the tadpole diagram shown in Fig. 3. For SU(N) gauge theory the parametric N-dependence of that contribution is which implies m 2  In order to obtain an estimate in physical units one may consider an initial energy density of about 5-25 GeV/fm 3 for RHIC experiments, and a projected factor of about two more for LHC energies. The inverse of the maximum primary growth rate for | A(t, p) | 2 then leads to the characteristic time scales The results are rather insensitive to the precise value of the initial energy density because they scale with the fourth root of ε. For comparison, a time scale associated with the largest observed secondary growth rate is about a factor of three shorter than what is given in (17). However, even though secondaries can reach considerably higher growth rates than primaries, they start later. As a consequence, a certain range of higher momentum modes can 'catch up' with initially faster growing infrared modes before the exponential growth stops, as seen in Fig. 1. This leads to a relatively fast effective isotropization of a finite momentum range, while higher momentum modes do not isotropize on a time scale characterized by plasma instabilities.
We study the process of isotropization in terms of the local energy density (9) which is a gauge invariant quantity. For this, we Fourier transformε(t, x) with respect to x and compute the absolute value of the ratiô If at some time t the system reaches an isotropic state the mean absolute value of (19) has to be one for all momenta. The time evolution of (19) is depicted in Fig. 5 for several momenta. Similarly to the SU(2) case one observes that the lowmomentum modes of the energy density isotropize, while the high-momentum modes are still anisotropic after the instability growth ends. Inspecting Fig. 5 one sees that the Fourier coefficients ofε become isotropic up to approximately | p | ≃ 1 − 2 ε 1/4 at the time of saturation. Using the same values for the physical energy density ε as for (17) one finds effective isotropization up to a characteristic momentum of about which approximately agrees with the case of the SU(2) gauge theory.

Conclusions
Nonabelian plasma instabilities in SU(3) gauge theory relevant for QCD exhibit a qualitatively similar behavior as previously observed for the SU(2) group. The main quantitative differences concern the reduction of primary growth rates by about 25 % for given initial energy density. These differences can be related to the parametric dependence on the number of colors of one-loop corrections to the self-energy. As a consequence, we find that the properly rescaled primary growth rates agree even quantitatively rather well for the different gauge groups. Though the nonlinear dynamics underlying the secondary growth rates leads to remaining quantitative differences between SU(3) and SU(2), this fact is less important for phenomenology. The primary growth rates determine the characteristic time scales for isotropization at low momenta, which turns out to be (even) slower than previously suggested by SU(2) results.