Universal behaviour, transients and attractors in supersymmetric Yang-Mills plasma

Numerical simulations of expanding plasma based on the AdS/CFT correspondence as well as kinetic theory and hydrodynamic models strongly suggest that some observables exhibit universal behaviour even when the system is not close to local equilibrium. This leading behaviour is expected to be corrected by transient, exponentially decaying contributions which carry information about the initial state. Focusing on late times, when the system is already in the hydrodynamic regime, we analyse numerical solutions describing expanding plasma of strongly coupled N=4 supersymmetric Yang-Mills theory and identify these transient effects, matching them in a quantitative way to leading trans-series corrections corresponding to least-damped quasinormal modes of AdS black branes. In the process we offer additional evidence supporting the recent identification of the Borel sum of the hydrodynamic gradient expansion with the far-from-equilibrium attractor in this system.


Introduction
Given the critical role of hydrodynamics in the physical picture behind the evolution of quark-gluon plasma created in heavy ion collisions, it is very important to have a reliable understanding of how relativistic systems in highly nonequilibrium initial states tend toward local equilibrium. It is useful to view this process as consisting of two stages. The first is characterised by the quasi-exponential decay of transient effective degrees of freedom -the nonhydrodynamic modes, which in the case of N = 4 supersymmetric Yang-Mills theory (SYM) plasma are in one-to-one correspondence with the quasinormal modes of AdS black branes [1,2]. The second stage is dominated by slowlyevolving hydrodynamic modes, whose time evolution is determined by conservation laws. This pattern holds both in microscopic theories and at the level of hydrodynamic models, such as the Müller-Israel-Stewart (MIS) theory [3,4] and its variants. This qualitative picture can be made quantitative and appreciated most clearly when examining certain special observables which enjoy universal behaviour at late times: they evolve (up to exponentially small corrections) in a way determined solely by some microscopic model and independently of the initial conditions. Such quantities, to which one may refer as universal observables, exhibit distinct attractor behaviour: for a wide range of initial conditions numerical solutions tend to a distinguished attractor solution in an approximately exponential fashion [5]. This universal behaviour often sets in while the system is still highly non-isotropic, and at later times coincides with the prediction of Navier-Stokes hydrodynamics. It is natural to interpret hydrodynamization in the sense of Refs. [6,7] as a manifestation of this phenomenon.
This set of ideas leads to the notion that the decay of transient modes describes the way the system approaches a far-from-equilibrium attractor [5,8,14,16] which can be thought of as "hydrodynamics beyond the gradient expansion" in the sense of Refs. [12] and [5]. The simplest and most widely studied example appears in the context of conformal Bjorken flow, where the pressure anisotropy turns out to have such universal behaviour when expressed as a function of the dimensionless variable w ≡ τ T (which is the proper time in units of the shear-stress relaxation time, up to a constant factor). This behaviour is captured by the gradient expansion of the pressure anisotropy, which in this particular case takes the form of a series in 1/w. More general expressions for such universal observables -valid also without assuming conformal symmetry or boost invariance -were recently described by Romatschke [9].
An important point is that the gradient expansion itself does not depend on the initial state of the system. However, this asymptotic solution receives exponentially-suppressed corrections which on the one hand reflect the spectrum of nonhydrodynamic modes, and on the other carry information about initial conditions. Specifically, at the linearized level each black-brane quasinormal mode introduces a transseries sector, which is an infinite series of exponentially damped corrections, with the damping rate determined by the quasinormal mode frequency. Each such sector enters with an amplitude dependent on the initial state. The entire solution takes the form of a transseries [5,10] whose leading element is the hydrodynamic gradient expansion.
As first shown in Ref. [24], the hydrodynamic gradient expansion in N = 4 SYM plasma is divergent, and the precise form of the large order behaviour of the expansion coefficients contains information about the nonhydrodynamic modes of the system. In fact, each transseries sector contains an asymptotic series in 1/w. The transseries structure, together with resurgence relations connecting expansion coefficients in different sectors, gives a consistent solution when all the divergent power series are properly summed (e.g. using Borel techniques). The Borel summation of the hydrodynamic gradient series itself gives an approximation of the far-from-equilibrium attractor [13], while the non-trivial transseries sectors describe the dissipation of initial state information as the attractor is approached. Indeed, the transseries structure -which was developed in the study of asymptotic series (see e.g. [?] and references therein) and has a wide array of applications -is perfectly suited to capture this phenomenon.
It was shown in Ref. [10] that the complex-conjugate pair of least-damped quasinormal modes implies a specific form of the leading transseries correction to the hydrodynamic solution. Thus, in the case of N = 4 SYM the dominant terms are known precisely, apart from two amplitudes which reflect particular initial conditions. The calculations reported here are aimed at comparing the form of this leading correction with numerical solutions of the full time evolution starting from some randomly chosen initial states. Such calculations, based on the AdS/CFT correspondence, were performed in Refs. [6,7]. It is well-known that in the first approximation, at sufficiently large times the behaviour of these solutions approaches the prediction of the leading order of the gradient expansion. The issue addressed here is how the time evolution matches expectations based on a quantitative analysis of late-time asymptotics.
In a typical initial state many of the nonhydrodynamic modes will be excited, leading to the complex patterns seen at early times in the solutions of Refs. [6,7]. As transients decay, the system tracks the attractor more and more closely. Soon beyond the regime where hydrodynamization occurs, one should expect that the least-damped nonhydrodynamic modes will become the dominant deviation from purely hydrodynamic behaviour [11]. These subtle traces of the initial state must be present even deep in the hydrodynamic regime. The aim of the work reported here was to verify that these effects are indeed present. We show how one can isolate these contributions and confirm that they have the expected form, so one can say that late-time behaviour of the pressure anisotropy describes damped oscillations around the hydrodynamic attractor. We demonstrate that to see this effect one really needs more than just the gradient expansion truncated at some low order. Finally we compare the attractor of N = 4 SYM to the attractor of conformal BRSSS hydrodynamics [23] and find that the two are essentially identical once the transport coefficients are matched, which may be a part of the reason why relativistic second order hydrodynamics is so successful in the description of quark-gluon plasma [19].

Universal observables
The notion of a universal observable is implicit in much of the recent work on hydrodynamic attractors [5,8,9,[14][15][16][17][18][19][20][21][22], and is rooted in observations made in Ref. [6]. It is best illustrated in the case of Bjorken flow in conformal BRSSS hydrodynamics [23]. In terms of the longitudinal and transverse pressures P L , P T and equilibrium pressure at the same energy density P = E/3, the pressure anisotropy is defined by This quantity exhibits universal behaviour [7,24] at late times when expressed in terms of the dimensionless "clock variable" This can be inferred directly from the equations of BRSSS hydrodynamics, which imply the following equation for the pressure anisotropy A(w) [5,8]: where C η ≡ η/s is the ratio of the shear viscosity to entropy density, and C τπ , C λ 1 are dimensionless second-order transport coefficients. At late times, which corresponds to large values of w, this equation possesses an asymptotic solution of the form This solution contains no trace of initial conditions, and thus it is universal in the sense that all solutions tend to it regardless of the initial state. One can go beyond the asymptotic solution (4) by including exponential corrections -this is the subject of the following section. The leading term describing the way in which the pressure anisotropy approaches zero at late times depends on the shear viscosity to entropy ratio, so it will depend on the physical system under consideration. However, if one rescales the clock variable by introducing [15] w ≡ w 4πC η (5) then the leading asymptotic behaviour of the pressure anisotropy is universal not only in the sense of being independent of initial conditions, but is also independent of any transport coefficients. The choice of scaling in Eq. (5) is such that the leading term in Eq. (6) takes the form of Eq. (4) when the shear viscosity to entropy density ratio assumes the value C η = 1/4π. Although the discussion in this section was framed in the context of BRSSS fluid dynamics, it is in fact far more general. The asymptotic behaviour of universal observables expresses the fact that the predictions of any sensible theory of hydrodynamics approach those of Navier-Stokes theory at late times. While this is obvious by construction, it is only apparent if one examines a universal observable such as A(w). In the evolution of non-universal quantities such as the temperature, this basic fact will be obscured by the dependence on initial conditions.

The form of transients
The universal behaviour of A(w) makes it possible to study transient effects -the decay of nonhydrodynamic modes -in an unambiguous way. Since the asymptotic expansion given in Eq. (4) is independent of initial conditions, all solutions are guaranteed to behave in accordance with it up to exponential corrections. These effects, at least in hydrodynamic theories and in strongly coupled N = 4 SYM plasma, are captured by the transseries representation [5,10]. The need to include transseries corrections is also a reflection of the fact that the gradient expansion in Eq. (4) is divergent both at the microscopic level [15,24] (see also Ref. [25]) and in hydrodynamics [5,10].
The pressure anisotropy can be written as a sum of two contributions where A H represents the purely hydrodynamic, universal part of the solution and δA is a correction which will depend on the initial conditions, and will thus be different for different solutions. Due to the divergence of the gradient expansion A H has to be understood either in the sense of a truncation of the gradient series at some low order, or as the result of a applying a procedure such as Borel summation. The precise form of the correction depends on the spectrum of nonhydrodynamic modes. Roughly speaking, at the linearized level each nonhydrodynamic mode makes a contribution of the form δA(w) = σw β e −Aw Φ(w) (8) where Φ is an infinite series in powers of 1/w whose coefficients are determined in terms of the parameters of the theory, as are the constants A and β. The amplitude σ however depends on the initial conditions. In general all these parameters are complex. In the case of N = 4 SYM each AdS black-brane quasinormal mode makes a contribution of the form (8), where the parameter A is proportional to the complex QNM frequency. The precise form of the correction due to the complex-conjugate pair of least-damped nonhydrodynamic modes can be inferred from Refs. [10,11,24,26] and reads In this equation, Φ ± (w) denote infinite series of the form similar to Eq. (4) appearing in the first transseries sectors, but with a constant leading term: The quantities σ ± play the role of integration constants. The remaining coefficients a (±) n are independent of initial conditions, but their specific values are not relevant for the present study and we will approximate the entire series in Eq. (10) by constant amplitudes, i.e. by the leading contributions given by the integration constants σ ± . In what follows we will fit these amplitudes to numerical data from AdS/CFT simulations.
The values of the remaining parameters appearing in Eq. (9) are known. The least-damped AdS black brane QNM frequencies [27] (corresponding to an operator of conformal weight ∆ = 4, apart from a factor of π) are given by: Ω R ≈ 9.800, Ω I ≈ 8.629 (11) and, furthermore β R ≈ 0.6866, β I ≈ 0.7798 (12) which can be extracted from Ref. [24]. There are two types of correction to Eq. (9): nonlinear terms, which are damped by powers of the exponential already appearing there and further contributions of a similar form, but with parameter values corresponding to more strongly-damped QNM of the AdS black brane. These contributions are subleading and will be ignored in the following.

Matching numerical evolution to the leading QNM
We now turn to examining numerical solutions of Bjorken flow with the aim of confronting their late time behaviour with the expectations discussed in the previous section. Suitable solutions can be obtained using the approach of Ref. [7] (using methods developed earlier in Refs. [6,28]). The focus of that study was on hydrodynamization, so randomly generated initial conditions were evolved only until hydrodynamic behaviour was identified -this occurred overwhelmingly for w < 1. Note that one could contemplate interpreting solutions at times preceding hydrodynamization in terms of compositions of quasinormal modes, in the spirit of what was done in Refs. [29,30] for the isotropisation problem. In the present study however, we are interested in transient effects after hydrodynamization, which requires evolving the asymptotically-AdS geometry to much later times, where only the least damped pair of quasinormal modes is relevant. This is based on the observation that while at early times all the quasinormal modes can a priori play an important role in a given solution, some time after hydrodynamization occurs only the least-damped QNM should dominate deviations from the asymptotic form given in Eq. (4).
To isolate this effect explicitly we will first make use of the observation also used in Ref. [31], which is that the universal part A H will cancel in the difference of any pair of solutions A 1 (w), A 2 (w), so that where the coefficients C,C are given by the differences of amplitudes σ ± of the QNM contributions to each of the solutions A 1 (w), A 2 (w). Since apart from these amplitudes all the parameters appearing in Eq. (13) are known, one can fit the constants C,C and check if one can reproduce the behaviour of differences of numerical solutions. As seen from Fig. 1, at early times this fails, but already at w > 3/2 it works very well. This is completely in line with expectations. Only one pair of solutions was used to generate Fig. 1, but we considered a number of such pairs and confirmed that they all lead to the same conclusion.
There is an alternative to the approach described above. Instead of considering the difference of two numerical solutions one can consider the difference between a numerical solution and some approximate representation of the universal hydrodynamic contribution A H . As mentioned earlier, there are two obvious choices here: a truncation of the series in Eq. (4), or its Borel sum. If either of these choices provides an accurate representation of A H we should find that subtracting it from any numerical solution will again leave an exponentially-damped correction of the form (13).
If one considers, instead of Eq. (13), the difference between a numerical solution and a truncation of the gradient series, one finds a situation such as that depicted in left plot of Fig. 2, which results in the case of truncating the gradient expansion at first-order. Clearly, such an approximation is inadequate for the present purpose. It turns out that going to second order does not improve the situation.
The other option is to use the Borel sum of the gradient expansion of N = 4 SYM which was recently computed in Ref. [13]. It was also shown there that the Borel sum acts as a far-fromequilibrium attractor for numerical solutions such as those discussed in the previous section. This fact can be taken as supporting evidence for the idea that this sum represents "hydrodynamics beyond the gradient expansion" and gives motivation for considering it as a useful estimate of A H . Indeed, as seen from the right-hand plot in Fig. 2, subtracting the Borel sum is, at least qualitatively as effective as subtracting a full numerical solution.
Note that if we were to study the difference between some numerical solution and the attractor determined as a solution of the AdS equations of motion (such as was done in Ref. [16]) this would just be a special case of the argument given earlier in this section, since the attractor is (after all) also a solution -albeit a rather special one. However, if we instead take the Borel sum of the gradient expansion [13] this can be interpreted as a test of whether such a procedure yields a result which is exponentially close to a full numerical solution.

Conclusions and outlook
The existence of universal observables which exhibit attractor behaviour and the availability of very precise numerical simulations of Bjorken flow [7] makes it possible study transient effects in the late time behaviour of numerical solutions of the full nonlinear evolution equations based on the AdS/CFT correspondence. The fact that these effects can be explicitly detected in the numerical simulations supports the general picture of hydrodynamization developed in a number of recent works [?, 5, 7, 8, 14-22, 24, 26, 32]. It is gratifying that these effects can be found in the precise form expected on the basis of the identification of transient, nonhydrodynamic modes of the expanding plasma with the quasinormal modes of AdS black branes and the realization that the hydrodynamic gradient expansion is the leading element of a transseries.
We have also demonstrated that the Borel sum of the gradient series, which was calculated and shown to act as an attractor in Ref. [13], is a useful proxy for the notion of "hydrodynamics beyond the gradient expansion" [5,12]. In view of this it is interesting to ask to what extent the attractor of N = 4 SYM coincides with the attractor of MIS theory, which is most commonly used to model the hydrodynamic stage of evolution of quark-gluon plasma created in heavy-ion collisions. The result of such a comparison is presented in Fig. 3, where the BRSSS variant of MIS theory is used, with the shear viscosity, relaxation time and λ 1 transport coefficients fitted to their N = 4 SYM values known from fluid-gravity duality [33]: It is striking that the two attractors coincide almost as soon as the Borel summation can be considered reliable (and coincides with the result of a calculation reported in Ref. [16] which was based on a different approach). This is in contrast to the attractor of kinetic theory in the relaxation time approximation (RTA) [16,34,35], which is not very well reproduced by BRSSS hydrodynamics. Indeed, in a recent study [17] the attractor of RTA kinetic theory was compared to the MIS attractor, as well as to the attractor of anisotropic hydrodynamics (see e.g. Ref. [36,37]). In that case the MIS attractor was found to converge to the exact kinetic theory result rather late, with anisotropic hydrodynamics reproducing it well already at very early times. Similar results were found in Ref. [18], where also the DNMR hydrodynamic theory [?] was considered. In this context note that the recent study of Ref. [38] suggests that despite similarities in hydrodynamization [15] the approaches based on kinetic theory and on the AdS/CFT correspondence are quite far apart by some metrics. It is an open question at the moment which of these two paradigms captures essential features of quantum chromodynamics more closely.