Hydrodynamization in kinetic theory: Transient modes and the gradient expansion

We explore the transition to hydrodynamics in a weakly-coupled model of quark-gluon plasma given by kinetic theory in the relaxation time approximation with conformal symmetry. We demonstrate that the gradient expansion in this model has a vanishing radius of convergence due to the presence of a transient (nonhydrodynamic) mode, in a way similar to results obtained earlier in strongly-coupled gauge theories. This suggests that the mechanism by which hydrodynamic behaviour emerges is the same, which we further corroborate by a novel comparison between solutions of different weakly and strongly coupled models. However, in contrast with other known cases, we find that not all the singularities of the analytic continuation of the Borel transform of the gradient expansion correspond to transient excitations of the microscopic system: some of them reflect analytic properties of the kinetic equation when the proper time is continued to complex values.

1. Introduction and summary.-Heavy-ion collisions at RHIC and LHC provide an outstanding opportunity to test our understanding of QCD. Perhaps unsurprisingly, a fully ab initio theoretical description has turned out to be very challenging. This has led to the exploration of models of increasing complexity and often nonoverlapping domains of validity -see, e.g., Refs. [1,2] for a review of such models in the context of the hydrodynamic description in ultrarelativistic heavy-ion collisions.
In this Letter we focus on the poorly understood transient far-from-equilibrium regime, which precedes viscous hydrodynamic evolution of quark-gluon plasma (QGP). There are two approaches to the study of the transition to hydrodynamics (hydrodynamization) in non-Abelian gauge theories like QCD: a weakly-coupled description based on effective kinetic theory (EKT) [3] (see also Refs. [4][5][6][7][8][9]) and a strongly-coupled plasma paradigm based on holography [10] (see Refs. [11][12][13][14] for sample results). They involve very different physical pictures and mathematical frameworks: the first relies on the Boltzmann equation while the second makes use of higherdimensional Einstein equations. Since under experimental conditions the QCD coupling is neither parametrically small nor large, it is crucial to understand which implications of these approaches can be viewed as universal.
Our aim is to shed light on equilibration in weaklycoupled systems by examining large-order behaviour of the hydrodynamic gradient expansion [15] in the framework of the kinetic theory model given by the Boltzmann equation in the relaxation time approximation (RTA) [16]. We also assume conformal symmetry, as its breaking at strong coupling does not significantly alter equilibration processes [17,18], whereas at weak coupling at vanishing quark masses it is a NLO effect (see, e.g., Ref. [5]). The key result of this Letter is demonstrating the vanishing radius of convergence of the gradient expansion in this kinetic theory model (see also Ref. [19] which studied in this context the RTA kinetic theory with constant relaxation time) and understanding some puzzling features revealed by these studies. The reason for the divergence turns out to be the same as in the case of holographic plasma: the presence of fast-decaying (nonhydrodynamic) modes [20][21][22][23] whose relaxation controls the emergence of hydrodynamic behaviour (and in particular, its applicability to the physics of heavy-ion collisions [24][25][26]). This is connected with the existence of attractors which govern the evolution far from equilibrium [21,[27][28][29][30][31][32][33][34][35]. In this context, see also Refs. [36,37] for studies of possible manifestations of analogous fastdecaying modes in trapped Fermi gases close to unitarity, as well as Ref. [38] where various features of transient modes are analyzed as a function of the microscopic interaction strength in a holographic toy model.
Quite remarkably, our analysis of the Borel transform of the gradient expansion in the RTA kinetic theory reveals not only the expected purely decaying mode [22], but also singularities which could naively be interpreted as transient contributions to the energy-momentum tensor exhibiting damped oscillatory behaviour, similar to the findings of Ref. [39]. This would however be surprising, since the mechanism described there does not apply to the RTA theory [40]. In fact, we demonstrate below that these singularities are instead a manifestation of analytic properties of the evolution equations in complexified time. This feature is related to what has been observed in other contexts where large-order behavior of perturbative series expansions is used to draw conclusions about non-perturbative effects (see, e.g., Refs. [41][42][43][44]).
Our conclusions concerning hydrodynamization provide strong motivation for comparing numerical solutions of the RTA evolution equations with the EKT results at intermediate coupling reported in Refs. [8,9] and with arXiv:1609.04803v4 [nucl-th] 10 Jun 2018 the AdS/CFT-based simulations of Ref. [45]. We uncover semi-quantitative agreement and, as a byproduct of this analysis, we present a new and effective way of visualizing (see Fig. 2) the correlation between the hydrodynamization time and the value of the η/s ratio noted in Ref. [46].
2. Kinetic theory.-We address the issues discussed above in the context of Bjorken flow [47], which is conveniently formulated in (proper time)-rapidity coordinates τ -y. They are related to Minkowski lab-frame coordinates t-z by t = τ cosh y and z = τ sinh y where z is the collision axis. Assuming translation symmetry in the transverse plane x T , the on-shell distribution function f depends only on the proper time τ , the modulus of the transverse momentum p ≡ |p T | and the boost-invariant variable u = τ 2 p y .
In the RTA the collision kernel appearing in the Boltzmann equation is linearized around the the equilibrium distribution which for simplicity we take to be Boltzmann .
The RTA Boltzmann equation takes the form To ensure conformal symmetry, we assume that the relaxation time is of the form where γ is dimensionless and is the only parameter of this model. The dependence of temperature on the proper time is determined dynamically by imposing the Landau matching condition [48,49] where is the energy density (per particle species).

Hydrodynamic Gradient
Expansion.-In a conformal theory the eigenvalues of the expectation value of the energy momentum tensor in a boost-invariant state are functions of the proper time τ alone. They are given by the energy density E and the longitudinal and transverse pressures P L and P T : Away from equilibrium P L and P T differ from the equilibrium pressure at the same energy density P ≡ E/3.
It is convenient to study the approach to equilibrium by examining the behaviour of the pressure anisotropy as a function of the dimensionless variable w ≡ T τ . The gradient expansion of A takes the form This follows directly from Eq. (7) and Eq. (6) if we use the fact that in conformal theories near equilibrium E ∼ T 4 and for boost invariant flow To determine the coefficients a n we look for a solution of the Boltzmann equation Eq. (2) in the form Inserting Eq. (9) into the Boltzmann equation Eq. (2) one can algebraically determine the functions h n in terms of the unknown coefficients a n . A key step in doing this is to eliminate proper time derivatives of temperature in favor of A and then using Eq. (8).
The Landau matching condition Eq. (4) implies that at each order n > 0 This condition amounts to a linear, algebraic equation which determines the expansion coefficient a n−1 .
Proceeding this way we have calculated the expansion coefficients analytically up to order 426 (we include the result of this calculation in the supplemental material). As a cross-check, we have also calculated the first 10 terms using two independent, albeit slower, methods from Refs. [50] and [48,49], noting perfect agreement.
The leading expansion coefficients read They can be used to match the transport coefficients of BRSSS hydrodynamics [20] to the RTA model. In particular, one finds At large orders, the expansion coefficients in Eq. (8) exhibit factorial growth. This demonstrates the vanishing radius of convergence of the hydrodynamic series, which parallels similar findings obtained numerically in N = 4 SYM using holography [15,51] as well as in hydrodynamics [21,52,53]. This is to be expected on general grounds, since the RTA theory contains, apart from hydrodynamic excitations, also a short-lived mode [22] whose physics is not captured by the truncated gradient expansion, as shown in Refs. [15,21]. We have checked that the structure of singularities remains stable as the number of terms kept in the series is varied. We have also checked that the residues of the pictured poles lie well above what was set as the numerical accuracy, i.e. that they are not numerical artefacts. The depicted singularities are discretizations of branch cuts with branch points at values of ξ given in Eqs. (16) and (20), as well as at ξ0+ξ± and 2 ξ± (see also Ref. [34] for a related statement based on a smaller data set from the previous version of the present manuscript).
The Borel transform removes the leading order factorial growth of the coefficients The inverse transform is given by the Borel summation formula and is not uniquely defined since the analytic continuation of A B (ξ) necessarily contains singularities that are responsible for the vanishing radius of convergence of the original series. We analytically continue the series from Eq. (13) truncated at 426 terms by means of Padé approximants. Fig. 1 shows the poles of the Padé approximant, which condense in a way known to signify a branch-point singularity [54]. In Ref. [21] (see also Refs. [15,[51][52][53]55]) the ambiguity in the Borel summation associated with singularities of the analytic continuation of the Borel transform has been argued to disappear once the gradient expansion is supplemented with exponentially decaying terms where ξ 0 denotes the beginning of a cut in the complex Borel plane. In the case under consideration the cut closest to the origin starts, up to 5 decimal places, at ξ 0 = 1.5000/γ.
The constants ξ 0 appearing in Eq. (15) can in general be complex (they then come in conjugate pairs) and in the examples analyzed so far in the literature they appear in positive integer multiples. Each such term comes with an infinite gradient expansion of its own and the term with the lowest ξ 0 in a given family comes with an independent complex integration constant. In all known cases [15,21,[51][52][53] (see Ref. [1] for a review), those least damped modes within a given family coincide with singularities of retarded equilibrium 2-point functions of the energy-momentum tensor at vanishing momentum.
In the context of the RTA kinetic theory, the studies of Ref. [22] reveal the presence of a zero-momentum fast-evolving mode in the isotropization of the energymomentum tensor to its equilibrium form: The nontrivial background flow is known to modify the above equation to the form [15,56,57] δ T µν ∼ e −ω0 T (x) uµ dx µ , which for the case of Bjorken flow, neglecting subleading terms at large values of w, reduces to This, together with Eq. (17), reproduces Eq. (16). Let us also note here that the analysis in Ref. [58] reveals that the cut along the real axis seen in Fig. 1 must be, in fact, an infinite collection of independent cuts. They all start at the same branch point, i.e. ξ = ξ 0 , but are characterized by different discontinuities and are interpreted as an infinite set of modes carrying information about the initial distribution function to late times w. Importantly, Fig. 1 contains also a pair of singularities characterized by If one were to apply straightforwardly the lessons from earlier studies of other models of expanding plasmas, these complex values of ξ 0 would be interpreted as oscillatory-type transient contributions to the pressure anisotropy. As shown in the next section, where we argue that these are unphysical, this natural-looking conclusion is premature. This is an important point, since such singularities appear also in other kinetic theory models, such as those with τ rel ∼ T −∆ , for 0 < ∆ < 3. Notably, for ∆ > 2, the unphysical modes are actually closest to the origin, so naively they would correspond to the dominant nonhydrodynamic corrections [58].
5. Analytic properties of RTA kinetic theory-In order to explain the singularities of the Borel transform at ξ ± , see Eq. (20) and Fig. 1, we use the integral equation [49] which follows from the Boltzmann equation [48] and directly determines the local energy density E(τ ) In the expression above, E 0 (τ ) carries information about initial conditions, but will not be relevant in the following analysis. The object of interest is the energy density E(τ ), which appears in as well as in τ rel (τ ) through Eqs. (3) and (4). The function H(q) originates from the second moment of the equilibrium distribution function, Eq. (1), and reads What will be crucial in the following is the analytic structure of H. Since the natural variable in our considerations is w and at late times w ∼ τ 2/3 , we shall write the argument of H as q = (w /w) 3/2 ≡ ζ 3/2 . Among the singularities of H ζ 3/2 , the one of interest is the branch point stemming from the square root in the denominator. It appears as a singularity when the arctan function in the numerator is taken in a nonprincipal branch. As a result, one finds which has branch points at third roots of unity. We will be interested in the ones located at The key observation is that the presence of singularities in the complex ζ plane leads to inequivalent choices of integration contours in Eq. (21), the only physical choice being homologous to the integration along the real axis. If one uses the late time solution for g(τ (w)) obtained from Eq. (11) under the integral in Eq. (21) and considers two inequivalent contours around ζ − (or, similarly, ζ + ), denoted C 1 and C 2 , one finds where the elipsis denotes terms subleading at large w.
Similarly to the analysis around Eq. (14), the branch cuts lead to contributions to δg of the form where we truncated subleading terms at large values of w. Finally, note that δg and δE are related through Eq. (22) with, at late times / large values of w, e τ τ 0 dτ τ rel (τ ) ∼ e 3 2 γ w which ultimately gives Comparing the exponent in the equation above with Eq. (20) one observes a remarkable agreement up to 4 decimal places. As a further way of corroborating this result, one can use the techniques utilized earlier in this context in Refs. [15,21] to match the square root character of the branch cut in Eq. (24) with the leading powerlaw correction in w to contributions from the ξ ± singularities in the Borel plane, with very good agreement. All this taken together gives us confidence that the correct interpretation of the singularities given in Eq. (20) is that they correspond not to physical excitations,ut rather to analytic properties of kinetic theory for complexified values of the w variable. The physical integration contour in Eq. (21) along real values of τ does not pick up these contributions (and, no wonder, they are not seen in solutions of the initial value problem discussed in Ref. [58]). The gradient expansion itself does not preclude unphysical choices of contour, and this is reflected in its largeorder behaviour. Similar phenomena can be seen in the integral equation considered in Ref. [59] or in the ghostinstanton story of Ref. [43]. Their origin goes back to the fundamental point of resurgence: all non-perturbative information is encoded in the large-order behavior of the perturbative sector.
6. Hydrodynamization compared-In previous sections we have demonstrated that hydrodynamization in the RTA takes place through the decay of transient nonhydrodynamical modes in complete analogy to the situation at strong coupling. On the other hand, the RTA shares structural similarities with EKT -while the collision kernels of the QCD effective kinetic theory have much richer structure than the RTA, many qualitative features at the level of kinetic theory coincide. In this section we further strengthen the connection between RTA and EKT by noting that the similarity between these two theories goes beyond abstract structural similarity and that they agree semi-quantitatively when the values of η/s are matched between the two theories.
To demonstrate the extent of quantitative agreement between different models, inspired by Ref. [46] we introduce a new, rescaled variablew ≡ w 4πη/s . This is useful for such comparisons, since the late time behaviour of the The dashed magenta curve represents first order hydrodynamics, the blue line is the holographic result, and the red, dashed-dotted line is from EKT. The green dotted curve stands for a solution of RTA starting from similar initial distribution as EKT and with the same shear viscosity, η/s = 0.624. Despite differences in microscopic dynamics, one sees significant qualitative and some quantitative similarities between different theories.

pressure anisotropy A is given by
The leading behaviour is completely universal and does not depend on the value of η/s. Deviations from the asymptotic form characterize contributions arising beyond first order hydrodynamics and indeed we say that the system has reached the hydrodynamic regime when for a given state the relative difference between A and A H remains smaller than some threshold value. Fig. (2) shows a comparison of the time evolution of the system evolved in the EKT from [9], RTA using the methodology of [48,49], and numerical AdS/CFT calculation of [13,45,60]. For the EKT and RTA simulations we took the initial condition used in [9], whereas for the AdS/CFT simulation we took a typical initial conditions from [45]. We evolved the systems using EKT with λ = 10 corresponding to η/s ≈ 0.642, holography with η/s = 1/4π and RTA with γ fixed to reproduce the value of η/s of either model. Note that every kinetic theory curve in Fig. (2) corresponds to a different initial distribution function. In all these models the evolution is similar, but distinct. Remarkably, in each case -despite vastly differing microphysics -the evolution converges to first order viscous hydrodynamics roughly at samẽ w ≈ 1. It is striking that in all these cases the pressure anisotropy at the time of hydrodynamization is as high as A ≈ 0.6 − 0.8.
The structural and quantitative similarities of EKT and RTA suggests that the gradient expansion in EKT also exhibits zero radius of convergence and that the weak coupling hydrodynamization is driven by the same qualitative process. This is connected with the notion of attractors which have been explored in both RTA and holography [28,29]. To what extent these insights translate to EKT is an important problem.
7. Acknowledgments-We would like to thank the authors of Ref. [39] as well as I. Aniceto