Plasmon mass scale in two dimensional classical nonequilibrium gauge theory

We study the plasmon mass scale in weakly coupled strongly interacting nonabelian gauge theory in a two dimensional configuration that mimics the boost invariant initial color fields in a heavy ion collision. We numerically measure the plasmon mass scale using three different methods: a Hard Thermal Loop (HTL) expression involving the quasiparticle spectrum constructed from Coulomb gauge field correlators, an effective dispersion relation and the measurement of oscillations between electric and magnetic energies after introducing a spatially uniform perturbation to the electric field. We find that the hard thermal loop expression and the uniform electric field measurement are in rough agreement. The effective dispersion relation agrees with other methods within a factor of two. We also study the dependence on time and occupation number, observing similar trends as in three spatial dimensions, where a power law dependence sets in after an occupation number dependent transient time. We observe a decrease of the plasmon mass squared as $t^{\frac{-1}{3}}$ at late times.


I. INTRODUCTION
The color glass condensate (CGC) [1,2] is an effective theory of quantum chromodynamics (QCD) at high energy. One of the remarkable predictions of CGC is that two colliding CGC sheets create gluon states with nonperturbatively high occupation numbers (1/g 2 ) [3], which can be described using classical fields [4].
The question of how this strongly interacting matter eventually isotropizes and thermalizes has been a longstanding question in the theory of ultrarelativistic heavy-ion collisions. According to our current understanding, the classical picture is valid for a short time after the initial collision, until the occupation numbers fall below unity. The out-of-equilibrium matter admits a kinetic theory description when the occupation numbers are ≪ 1/g 2 . Fortunately, the kinetic theory description also has an overlapping range of validity with classical simulations [5][6][7], meaning that the two can be smoothly matched. In recent simulations the equilibration process has been studied by matching kinetic theory and relativistic hydrodynamics with promising results [8].
The aim of this paper is to compare different methods to estimate the plasmon mass in classical Yang-Mills systems in the case where a three-dimensional theory exists in a two-dimensional configuration. In practice we implement this by using a three-dimensional calculation on a lattice with only one point in the z direction. This configuration mimics the very anisotropic, nearly boost-invariant, field configuration predicted by the CGC for the initial stage of a heavy-ion collision. At this stage we neglect the longitudinal expansion of the system in a heavy-ion collision and work on a fixed size lattice. Instead of calculating in the asymptotically large time regime, where a clear separation between the hard and Debye scales has developed, we also want to address relatively early times where it is less obvious that such a scale separation exists. At this regime of earlier times one can also study the dependence on the parameters of the initial condition, the occupation number, and the hard scale, separately.
To extract the plasmon mass we will systematically compare three methods, which we already compared in three dimensions [41]. The first method will be to use a formula one can derive in hard thermal loop (HTL) perturbation theory. Even though the HTL-like scale separation is not guaranteed by the weak coupling in the classical theory, the plasmon mass scale nevertheless exists [39,40,[42][43][44]. The second method involves perturbing the system with a spatially uniform electric field (UE) [45] and measuring the response to this zero momentum perturbation. The third method involves the effective dispersion relation (DR) [42], which we can extract from the Coulomb gauge correlators of the fields.
The HTL formula relies on the quasiparticle spectrum, which is also typically extracted from Coulomb gauge correlators of the classical fields. However, we will notice that at high occupation numbers in two spatial dimensions the gauge fixing has a deforming effect on the observed quasiparticle spectrum. Thus, we will argue for the need to use gauge-invariant observables for measuring the typical momentum scale and occupation number.
We will first briefly introduce the numerical methods and initial conditions in Sec. II. In Sec. III we will introduce the three methods we use to extract the plasmon mass. Then we move on to dependencies on the lattice cutoffs, time, and occupation number in Sec. IV. Finally, we conclude in Sec. V.

A. Equations of motion in the temporal gauge
We have done all the numerical simulations with the SU(2) gauge group for numerical convenience. Several studies have compared SU(2) and SU(3) gauge groups, demonstrating qualitatively similar results [46][47][48]. We use the standard pure gauge Wilson action on a threedimensional lattice N is the number of colors and U x;i are the link matrices defined as The plaquette variables appearing in (1) are defined as □ i;j x ≡ U x;i U xþi;j U † xþj;i U † x;j . The spatial lattice spacing is a s , and the time step is denoted by a t . We use standard normalization for the SU(2) generators, Trðt a t b Þ ¼ 1 2 δ ab . The relationships between the lattice variables and the actual fields are given by the following equations: We refer to the lattice spacing in the μ direction with a μ . The equations of motion of the electric field are obtained by varying the action (1) with respect to the spatial links where □ j;k x ¼ U x;j U † xþj−k;k U † x−k;j U x−k;k . When the electric field on the next time step is known, we can easily construct the temporal plaquette [when using the SU(2) symmetry group] by decomposing the temporal plaquette into two parts: In the temporal gauge the temporal plaquette simplifies to a product of link matrices at two different time steps, making it easy to solve for the link at the next time step. Color charge conservation is encoded in Gauss's law, which is preserved by the discretization algorithm.

B. Quasiparticle distribution
We extract the quasiparticle spectrum by eliminating the residual gauge freedom with the Coulomb gauge condition. The gauge fixing is done by a Fourier accelerated algorithm [49]. However, even with gauge fixing, there is no unique way to determine a quasiparticle distribution from a given classical field configuration (see also the discussion in Refs. [45,50]). If our system can be described by weakly interacting quasiparticles, the energy density of the system should be obtained as Here the factor 2ðN c 2 − 1Þ accounts for the number of color and transverse polarization states in the system. The number of physical polarization states the plasmons have is 3. However, the longitudinal mode is only present for modes close to the Debye scale and is not expected to contribute significantly to the total energy density in (10). Thus a factor 3 would lead to a significant underestimation of the occupation number of hard modes. The total energy of the system is given by the Hamiltonian We now keep only the terms which are quadratic in the gluon field and equate the energy given by the Hamiltonian with the one given by the quasiparticle spectrum. Solving for the quasiparticle spectrum gives Here jE C ðkÞj is the Coulomb gauge electric field and jA C ðkÞj is the gauge field in Coulomb gauge. This procedure also removes the magnetic part of the longitudinal polarization state. The energy of a mode with momentum k (the dispersion relation) is given by ωðkÞ.
When extracting the quasiparticle spectrum we will assume a massless dispersion relation ωðkÞ ¼ k. It is not immediately obvious what would be the correct procedure to selfconsistently include a plasmon mass in the dispersion relation used here. In any case the effect of such a correction would be of the same order as the higher order terms in the gauge potential that we are already neglecting in Eq. (12). Alternatively, we can also extract the quasiparticle spectrum using only the gauge fields, or only electric fields. These two should be equivalent above the Debye scale after the system has been evolved in time for a few Qt. The alternative definitions for the occupation number are for the electric estimator and for the magnetic estimator. We can also use a combination of electric and magnetic fields, We compare these different definitions for the occupation number in Fig. 1. We observe that the electric occupation number diverges in the infrared, whereas the magnetic estimator is IR finite. The combination of electric and magnetic fields (f EA ) also behaves reasonably well in the IR. This behavior is to be expected since below the Debye scale we would not expect the dispersion relation to be massless. We will return to this question later in the context of using these distributions. Figure 1 also demonstrates the transverse electric occupation number. We clearly see that it is very close to the total electric occupation number, especially at higher momenta, indicating that the quasiparticle spectrum is dominated by transverse plasmons.

C. Initial conditions
The initial gauge fields are sampled from the distribution Averaged over 20 configurations.
Since our system is two dimensional, the z components of the momenta are in fact always zero, and the exponential Here Q is the dominant momentum scale, V is the lattice volume, and n 0 is a parameter describing the typical occupation number of the system. This initial condition is a momentum distribution clearly peaked around Q, and it behaves very well in the ultraviolet and infrared regions. It also trivially satisfies Gauss's law since it contains only magnetic energy. In the context of the early stages of heavy-ion collisions, we should consider Q as analogous to the saturation scale Q s [51]. For the classical approximation to be valid, the occupation number f ∼ n 0 /g 2 should be greater than of the order of 1, i.e., n 0 ≫ g 2 .
Our choice of initial fields, using (12) as the quasiparticle distribution, results in the following form, 1 where a s is the length of the system in the z direction.
We are using the same initial condition as in Refs. [52,53]. This distribution is strongly cut off in the UV and is, in this regard, similar to initial conditions used in Refs. [54,55]. For more realistic initial conditions we refer the reader to, for example, Refs. [56,57]. The precise functional form of the initial condition does not very strongly affect the late-time behavior of the system (unless one changes the occupation number by a very large margin) because overoccupied classical Yang-Mills systems eventually evolve into a well-known scaling solution in a time of the order of a few Qt [45,54,58]. Similar findings have also been made in scalar field theories [55,59,60].

D. Measured observables
Our initial condition (16) is constructed with the gauge fields in momentum space. These are then Fourier transformed back to coordinate space and exponentiated to form the link matrices that are the actual variables in the calculation. In order to measure the quasiparticle spectrum (12), one then fixes the Coulomb gauge and calculates the anti-Hermitian traceless part of the link matrix to get the gauge potential A i appearing in the expression (12). This process is very nonlinear, and in the high density regime it is not obvious that one recovers the same quasiparticle distribution that one put in as an initial condition. In order to control the limits in this process, it is useful to compare the measured quasiparticle distribution at t ¼ 0 to the one used in the initial condition. This comparison is shown in Figs. 2 and 3. Figure 2 shows the dependence of the gaugefixed quasiparticle spectrum on the initial momentum scale. We observe that the gauge-fixed spectrum deviates the most from the analytical initial condition when Qa s is small. Figure 3 shows the gauge-fixed spectrum for different occupation numbers. Here we find a similar effect-at higher occupancy the gauge fixing significantly drives the spectrum away from the analytical initial condition by redistributing the energy to higher momentum modes and decreasing the occupation number correspondingly. This effect seems to be more dramatic in two dimensions than in three dimensions. (In three dimensions we did not observe this effect with n 0 ¼ 1, but we presume that this effect can be observed also in three dimensions, provided that one goes to high enough occupation numbers.) In order to reduce the gauge fixing effects, we project out the longitudinal components of the initial gauge field using the standard transverse projection operator Here p i and p j correspond to the (complex) eigenvalues of the discretized derivative operator. However, the anti-Hermitian parts of the links constructed by exponentiation from the transverse fields do not necessarily satisfy the Coulomb gauge condition to the desired accuracy. This problem is especially severe when the occupation numbers are high. The effect is also visible in Fig. 3. The curve with n 0 ¼ 0.5 overlaps with the analytical initial condition because here the transverse projection done on the initial gauge field also keeps the lattice Coulomb gauge violation so small that no additional gauge fixing is needed. However, when one goes to higher occupation numbers, additional gauge fixing becomes necessary and we start to observe deviations from the analytical distribution function. Due to this fact we want to measure the momentum scale and occupation number in a gauge-invariant manner. We estimate the typical momentum of the chromomagnetic field squared as done in [36,45], For our initial condition we can estimate this perturbatively, We define the effective momentum scale in such a way that it matches the initial Q in the dilute limit In order to estimate the occupation number, we first compute the initial energy density in terms of Q and n 0 using Eq. (10), This means that we can use the gauge-invariant momentum scale Q eff and the two-dimensional energy density ϵ 2d ≡ a s ϵ to define a gauge-invariant measure of the typical occupation number of gluons as The normalization of Q eff and n eff 0 has now been chosen in such a way that in the dilute limit and at t ¼ 0, they agree with the input parameters Q and n 0 . Away from the dilute limit we want to perform simulations by varying the lattice parameters in a way that maintains, as much as possible, fixed values of the gauge-invariant scales Q eff and n eff 0 . In order to do this, we must map the relation between the input parameters and the gauge-invariant scales. To this end we have performed a series of measurements at the initialization time. Figure 4 shows the dependence of the gauge-invariant scale Q eff on Q for various occupation numbers n 0 . The gauge-invariant scale is linear in the initial scale, but an increase in the initial occupation number n 0 also results in an increased gauge-invariant scale. This reinforces our interpretation of the phenomenon seen in the gauge-fixed spectra in Fig. 3, where we observed that, for a fixed input parameter Q, the resulting effective momentum scale increases when one increases the occupation number. It seems that the system resists attempts to increase the gluon density by increasing the amplitude n 0 , instead transferring the additional energy into modes with higher momenta. Figure 5 shows the connection between the gauge-invariant occupation number and the occupation number n 0 as a function of Q. The main conclusion of this figure is that at low Qa s it is impossible to go to high occupation numbers, but instead n eff 0 saturates to a certain value. At higher Q we can go to higher effective occupation numbers, but we are also further away from the continuum limit. It is also useful to compare Figs. 5 and 2. When keeping n 0 fixed in Fig. 5 and going to higher Q, we observe that also the effective occupation number also increases. Thus, the behavior of the gauge-fixed spectrum FIG. 4. Effective momentum scale as a function of the scale set by the initial condition for different n 0 . We observe a linear dependence, but increasing n 0 also increases the observed Q eff .
in Fig. 2 is in line with the behavior we find in gaugeinvariant occupation number in Fig. 5. To summarize, with the help of Figs. 4 and 5 we can establish a link between the initial simulation parameters and the initial measured parameters. However, the reader should bear in mind that, for example, taking the continuum limit (Q eff a s → 0 for fixed n eff 0 ) in the high occupation number regime is actually impossible due to the saturation in n eff 0 .
Both of these effective scales Q eff ðtÞ and n eff 0 ðtÞ are functions of time. The time dependence of the momentum scale is shown in Fig. 6. In practice we find that these scales evolve in time as Q eff ðtÞ ∼ t 1/7 − t 1/6 , and consequently n eff 0 ðtÞ ∼ t −3/7 − t −½ in two spatial dimensions. From now on we will be using the notation Q eff ðt ¼ 0Þ ¼ Q eff and n eff 0 ðt ¼ 0Þ ¼ n eff 0 .

III. METHODS FOR EXTRACTING THE PLASMON MASS
We use the same methods to extract the plasmon mass as in our previous work [41] and refer the reader there for more details.

A. Uniform electric field
In the three-dimensional case, it turned out [41] that the best method to extract the plasmon mass is probably the uniform electric field method [45]. In this measurement one introduces, by hand, a uniform chromoelectric field, corresponding to a perturbation with zero momentum. The response to this perturbation is measured in the total energy of the system, and the plasmon mass can be read off from the frequency of the oscillation between electric and magnetic energies. This procedure explicitly breaks Gauss's law, and one has to restore it by hand. The restoration is done using the algorithm described in Ref. [61].
We use two methods to measure the oscillation frequency. The first method is to fit a damped oscillation to the signal, as we have done in Fig. 7. Note that we subtract the time average of the energy in order to move the signal to oscillate around zero. In practice the oscillation gets damped quite quickly, and the whole time interval until the end of the simulation is dominated by noise. To treat this situation we perform two successive fits. First, by fitting to all the available data, we obtain a first estimate of the oscillation period. We then perform a second fit to a time interval of two periods of the first oscillation. The second fit gives the actual estimate for the frequency. The electric energy is proportional to the electric field squared, and thus two oscillation periods correspond to one oscillation in the electric field itself. In this way we could, in principle, also extract the damping rate. In practice we have found that this method systematically overestimates the damping rate, and the fitting procedure itself does not work as well in two dimensions as in three dimensions.
The second variant of the UE method is to compute the autocorrelation function of the electric energy and look at the separation of the maxima. We define the correlation function as where a and v are sequences. If the span of our data set in time is Δt, the correlation is computed from k ¼ −Δt 2 to k ¼ Δt 2 . The sequences are padded with zeros whenever When averaged over long times, it seems that the effective momentum scale approximately follows the t 1/6 power law. However, when looking at only the asymptotic regime, i.e., tQ eff > 500, it seems that Q eff grows faster than t 1/7 , but also definitely slower than t 1/6 . The results have been averaged over five runs. necessary to keep the sum well defined. When the correlation function is defined in this way, we expect the highest peak to appear, when k is half of the length of the sequences. When computed in this way, the autocorrelation function is also symmetric. An example of the autocorrelation function for the electric energy is shown in Fig. 8. Here we have subtracted the time average of the energies. At Q eff Δt ¼ 50 we have the correlation of the signal with itself without any lag. The period of the oscillation can (remembering that the signal is the square of the electric field) be extracted by looking at the distance between the peak at Q eff Δt ¼ 50 and the first maximum.
Note that in Fig. 7, the actual frequency of the oscillation seems to be located between the two methods. Mathematically we would expect the damped oscillation fit to give better results than the autocorrelation method. The reason for this is that damping also shifts the location of the maxima and minima from the values given by the oscillation frequency. The autocorrelation method, in contrast, really looks for the peaks and dips of the oscillation.

B. Dispersion relation
We extract an effective dispersion relation from the electric field and its time derivative as allowing us to separately study the transverse (T) and longitudinal (L) polarization. As we did for three dimensions, one can also compare the electric field and the gauge potential: as was done in [42] for (2 þ 1)-dimensional gauge theory.
In the case of three dimensions we found [41] that Eq. (26) significantly underestimates the plasmon mass. Typical examples of the dispersion relations (25) and (26) are shown in Fig. 9. We also see here that the method using the fields Eq. (26) gives a lower value for the mass gap. The difference is, however, nowhere near as drastic as in three dimensions, and actually the values given by (26) are closer to the values given by other methods in two dimensions. Perhaps more importantly, the evidence for a mass gap with the estimator (26) vanishes at larger momenta, while for Eq. (25) the value remains consistent with the gap at zero momentum. Also, since in Coulomb gauge the gauge potential A i has no longitudinal component, but the magnetic part of the longitudinal plasmon is hidden in the nonlinear terms, we cannot study longitudinal modes separately with (26). Thus, while (26) gives much more reasonable results in two dimensions than it does in three, we still consider the estimate (25) using the time derivative of the electric field a better one and will use it in the following. The extraction of the plasmon mass is done by using a linear fit of the form ω 2 ¼ ω 2 pl þ ak 2 (with two free parameters ω 2 pl and a) to Eq. (25) The maximum momentum FIG. 8. Autocorrelation function of the electric energy. Due to our definition of the correlation function, the peak at Q eff Δt ¼ 50 corresponds to the autocorrelation of the electric field without any lag. The distance between this peak and the first maxima gives the period of the oscillation divided by 2. FIG. 7. Oscillation between the electric (E E ) and magnetic (E B ) energies after the addition of the homogeneous chromoelectric field. Both fields have been shifted to oscillate around 0 by subtracting the time average after the addition of the homogeneous chromoelectric field. The curve labeled as "fit" shows the result of the damped oscillation fit. The curves "Fit freq." and "Autocor. freq." show the oscillations with the frequencies extracted from the fit and autocorrelation function without any damping. The fitting region here is constrained to the part of the oscillation visible in this plot. The uniform electric field was introduced at Q eff t ¼ 160.
we use in the fit is k/Q eff < 0.25. This cutoff has been chosen by experimenting with the fitting procedure. Choosing a significantly larger cutoff results in the fit overshooting the constant parameter. However, one still wants to have plenty of statistics for the fit, which is why one should not choose a cutoff that is too small.

C. HTL resummed approximation
If the HTL kind of separation of scales is valid, the plasmon mass is given by the integral On the lattice the integral is discretized by the standard replacement where k runs over the modes available on the lattice. This method is also widely used in the literature-see, e.g., Refs. [19,43,58]. While estimating the plasmon mass scale using Eq. (27), one can use different definitions for the particle distribution, as discussed in Sec. II B. As we saw from Fig. 1, the occupation number given by the estimators f E or f AþE is IR divergent. Because of the different phase space, this effect is much stronger in two dimensions than it is in three. In fact, an infrared convergent value for the integral (27) with the estimators f AþE or f E would require the EE correlator to approach zero for k → 0. This would be a large suppression compared to the thermal value jEðkÞj 2 ∼ T, and would, in a contradictory fashion, correspond to a vanishing mass gap in the dispersion relation (26). Thus, we use only the estimators f A and f EA that yield a finite value to compute the plasmon mass using Eq. (27). We refer to these as HTL-A and HTL-EA.

IV. DEPENDENCE ON LATTICE CUTOFFS, TIME, AND OCCUPATION NUMBER
Next we study how our results depend on lattice cutoffs. The ultraviolet cutoff is set by Q eff a s , and the infrared cutoff is controlled by Q eff L. In practice when we study the dependencies we vary one of these cutoffs and keep the other fixed. The fact that we are using the gauge-invariant observables Q eff and n eff 0 makes the choice of initial parameters Q and n 0 rather complicated. The necessity to invert the relation shown in Figs. 4 and 5 to choose Q and n 0 corresponding to fixed Q eff and n eff 0 introduces some additional uncertainty into these estimates. The error bars shown in the figures in this section are statistical errors computed as the standard error of the mean. It turns out that these errors are insignificant for HTL and UE methods, but for the DR method the statistical errors are larger. The infrared cutoff dependence is shown in Fig. 10. All observables seem to be well behaved in the infinite volume limit. Our results on the ultraviolet cutoff dependence are shown in Fig. 11. A remarkable feature is that values given by all methods seem to increase when we approach the continuum limit, in contrast to the behavior observed in the three-dimensional system [41]. However, the results do not seem to be UV divergent.
The time dependence of the plasmon mass using various methods at early times is shown in Fig. 12. We find that the time evolution computed with the HTL-A method agrees reasonably well with a t −1/3 power law. However, it seems that the UE method decreases slightly faster than the HTL-A method and the power law. We also find that the difference between DR and the other methods persists even at late times. In order to study the asymptotic power-law behavior on longer time scales, we do a study using only the HTL-A method (using the UE method would require prohibitively many separate runs). The results are shown in Fig. 13. We find that at late times the HTL-A method is in agreement with the proposed power law.
The occupation number dependence is depicted in Fig. 14. The dependence seems qualitatively similar to FIG. 12. Time dependence of the plasmon mass scale with various methods, with labels as in Fig. 10. We also show a curve corresponding to a t −1/3 power law. We find that the HTL-A method almost agrees with this power law. It seems that the value given by the UE measurement decreases slightly faster than that of HTL methods and the power law. Results are averages over five runs. We also show the t −1/3 power law here for comparison. It seems that at late times the HTL-A method agrees with the power law. The results are averaged over five runs.
FIG. 11. The UV-cutoff dependence of the plasmon mass, with labels as in Fig. 10. Here QL ranges from 246 to 270 and n eff 0 ranges from 0.147 to 0.162. We observe that the results given by all methods increase when we approach the continuum limit. The lattice sizes we used were 376 2 ; 514 2 ; 646 2 ; 763 2 ; 1008 2 ; 1326 2 ; 1575 2 ; 1800 2 with averages taken over 10, 9, 8, 7, 6, 5, 4, 3 simulations. FIG. 14. Dependence of plasmon mass (scaled by the occupation number n eff 0 ) on the occupation number n eff 0 for the different methods of evaluating the plasmon mass scale at fixed time, with labels as in Fig. 10. We observe a similar trend as we did in three dimensions: the plasmon mass scale at a fixed time decreases as a function of the occupation number. Here Q eff L ¼ 432 − 480, Q eff ¼ 0.240 − 0.267, and ðQ eff tÞ UE ¼ 205 − 227. The bump at n 0 ≈ 0.24 is caused by a deviation in these parameters. These results are averaged over five configurations. the three-dimensional case-the plasmon mass scale (normalized by the occupation number) falls as a function of increasing n eff 0 . The bump in the figure roughly at n eff 0 ¼ 0.24 is caused by a deviation in the input parameters. It also seems that the differences between the methods are independent of occupation number. This is suggestive of a similar interpretation as the three-dimensional results [41], namely, that larger occupation numbers lead to a more rapid start of the large time scaling regime where the plasmon mass begins to fall with time.

V. CONCLUSIONS AND OUTLOOK
In this paper we have first argued for a need for gaugeinvariant observables to measure the occupation number and momentum scale. This need arose due to gauge fixing effects, which turned out to deform the spectrum more than one might have expected based on our previous threedimensional simulations.
According to our observations the UE and HTL-A methods are in rough agreement, and they are almost equally good in measuring the plasmon mass scale. However, the HTL-A method is computationally much cheaper and easier to implement, so we recommend that the reader use it when measuring the plasmon mass scale in two-dimensional simulations. However, one should bear in mind that in two spatial dimensions, the dependence on the precise definition of the quasiparticle distribution is much larger than in three because the HTL integral is more dominated by the infrared. The DR method works in a similar fashion as in the three-dimensional case-the agreement with the other methods is within a factor of 2.
We have also studied the cutoff dependencies of our results. We find no cutoff dependencies when we vary the infrared cutoff. However the choice of the ultraviolet cutoff does influence the results we obtain. It turns out that the values given by all methods increase when we get closer to the continuum. This behavior is different from the three-dimensional case, where we observed that in the continuum limit the agreement between the UE and HTL results seemed to improve. We find no evidence of this in two spatial dimensions.
The dependence on the occupation number is similar to the behavior in three dimensions, and the differences between the methods persist regardless of occupation number. We observe that the plasmon mass scale squared seems to decrease like t −1/3 at late times when measured with the HTL-A method. We also observe a reasonably good agreement at earlier times. However, at earlier times the UE measurement seems to decrease slightly faster than this power law.
It would be interesting to study whether two-dimensional classical Yang-Mills system can be understood in a kinetic theory framework. However, in two dimensions the contribution of the modes below the Debye scale becomes a lot more important than in three dimensions, which we can see from the deviations between the different HTL estimates.
We are planning to address the dispersion relation and spectral function in classical Yang-Mills systems in a future publication using our recently developed techniques to perform simulations with fluctuations on top of the classical background [62]. In this way we could also address the existence of quasiparticles in classical Yang-Mills theory, which would shed light on the kinetic theory description of classical Yang-Mills systems.