Heavy Operators and Hydrodynamic Tails

The late time physics of interacting QFTs at finite temperature is controlled by hydrodynamics. For CFTs this implies that heavy operators -- which are generically expected to create thermal states -- can be studied semiclassically. We show that hydrodynamics universally fixes the OPE coefficients $C_{HH'L}$, on average, of all neutral light operators with two non-identical heavy ones, as a function of the scaling dimension and spin of the operators. These methods can be straightforwardly extended to CFTs with global symmetries, and generalize recent EFT results on large charge operators away from the case of minimal dimension at fixed charge. We also revisit certain aspects of late time thermal correlators in QFT and other diffusive systems.


Introduction
The analytic conformal bootstrap has uncovered universal features in sparse corners of the spectrum of conformal field theories (CFTs), at large spin [1,2] or large charge [3]. The 'middle' of the spectrum is instead exponentially dense, but reveals universal properties as well [4,5]. Some of these advances were guided by the existence of semiclassical descriptions, such as weakly interacting probe particles in AdS [1] for large spin states, or a superfluid effective field theory (EFT) for large charge states of certain CFTs [6,7,8,9]. The middle of the spectrum also enjoys a natural semiclassical description: thermodynamics [10], and more generally hydrodynamics. The subject of this paper is to study the consequences of this description.
Hydrodynamics is expected to emerge as the late time dynamics of any non-integrable quantum field theory (QFT) at finite temperature. The first theoretically controlled demonstration of this phenomenon is possibly Landau's two-fluid model [11]; for weakly coupled QFTs the emergence of hydrodynamics is now well understood within the framework of Boltzmann kinetic theory [12,13]. The fluid-gravity correspondence is a more recent example [14,15,16], for strongly coupled holographic theories. Although an analogous proof in generic CFTs may be too formidable a task for the conformal bootstrap, analytic methods may be able to place constraints on hydrodynamics, such as bounds on transport parameters [17].
The approach followed here is instead to work from the bottom-up, with the hope to guide future efforts from the analytic or numerical bootstrap. Hydrodynamics tightly constrains the thermal correlator of any light neutral operator (e.g. any Z 2 -even light operator in the 3d Ising model) at late times. This regime is difficult to address with conventional CFT methods because large Lorentzian times t β are far outside of the radius of convergence of the operator product expansion (OPE) [18]. In the microcanonical ensemble, hydrodynamics controls heavy-light four-point functions HLLH far from the LL OPE limit. Assuming typicality of heavy operators, hydrodynamic predictions can be recast as expressions for off-diagonal heavy-heavy-light OPE coefficients C HH L . Our results, summarized below, should hold in any non-integrable unitary CFT in three or more spacetime dimensions.

Summary of results
We consider thermalizing (or chaotic) CFTs in d + 1 spacetime dimensions. Operators that do not carry any internal quantum numbers acquire thermal expectation values: for example a neutral dimension ∆ O scalar satisfies where β is the inverse temperature, and b O a coefficient that is generically O (1). As argued in Ref. [10], consistency with the microcanonical ensemble implies that diagonal heavyheavy-light OPE coefficients are on average controlled by the thermal expectation value (1.1). Assuming typicality of heavy eigenstates allows one to drop the averages and leads to the prediction [10] (dropping numerical factors) (1.3) We show under the same assumptions that lead to (1.2) that this hydrodynamic correlator implies , (1.4) for the OPE coefficient of the light operator O with heavy operators of dimension ∆, ∆ and spin J, J . Off-diagonal OPE coefficients are exponentially suppressed in the entropy , as expected on general grounds [4]. We have dropped subexponential dependence on ∆, but instead emphasize the singular dependence on ∆ − ∆ and J − J featuring the hydrodynamic sound pole. This result holds for heavy operators satisfying . (1.5) The difference in spin must satisfy the same upper bound J − J (∆/b T ) 1/(d+1) . This upper bound comes from the UV cutoff of hydrodynamics, which only describes dynamics at times larger than the thermalization time t τ th . The lower bound comes from IR effects which resolve the singularity in (1.4). In (1.5) we have assumed τ th ∼ β; weakly coupled CFTs have τ th β and the window (1.5) is parametrically smaller.
Hydrodynamics pervades late time correlators, and not just those of scalar operators. In a thermal state, neutral operators of any integer spin can decay into composite hydrodynamics operators -this is illustrated in Fig. 1. Consider an operator of spin . Its component with an even number¯ of spatial indices with 2 ≤¯ ≤ has the same quantum numbers as composite hydrodynamic fields involving the stress tensor T µν This equation is not meant as a microscopic operator equation in the CFT, but rather as an operator equation in the low-energy (dissipative) effective theory around the thermal state.
The first term shows that the operator overlaps linearly with hydrodynamic excitations. Its two-point function will therefore contain hydrodynamic poles, leading to OPE coefficients similar to (1.4). If we consider this operator at vanishing wavevector k = 0, then the leading term drops because it is a total derivative and the operator no longer overlaps linearly with hydrodynamic modes. However, it can still decay into the second composite operator which leads to a hydrodynamic loop contribution to its correlator (1.7) Although this universal late-behavior for thermal correlators of generic operators in QFTs can be straightforwardly derived using the time-honored framework of fluctuating hydrodynamics, it has to our knowledge not appeared previously in the literature.
We further derive generalizations of Eqs. (1.4) and (1.8); these results apply to any non-integrable CFT in spatial dimensions d ≥ 2, without additional continuous global symmetries. Continuous global symmetries G can be incorporated straightforwardly: they lead to additional hydrodynamic modes which can give further contributions to OPE coefficients. We illustrate this with the case G = U (1). OPE coefficients involving charged heavy operators are similar to (1.4) and (1.8), with some differences for odd-spin light operators which receive larger hydrodynamic contributions because of the new slow density.
The U (1) symmetry can be spontaneously broken in the state created by the heavy operator of large charge. In this case, the hydrodynamic description includes a Goldstone phase. This allows us to connect to the large charge program [6,7,8,9,3,22], which can be thought of as a special case where a hydrodynamic (or semiclassical) description survives the T → 0 limit thanks to the spontaneous breaking of the U (1) symmetry. The various possible phases created by heavy operators are shown in Fig. 2.
The rest of this paper is organized as follows: Fluctuating hydrodynamics is reviewed in Sec. 2, and applied to relativistic QFTs. A few novel results are also obtained there,  Existing analytic methods to study various regions of the spectrum include the light-cone bootstrap [1,2], Tauberian theorems [4,5], and the large charge limit [6]. βµ < ∞ to connect the hydrodynamic description, and speculate on signatures of thermal phase transitions in the spectrum of heavy operators.

Hydrodynamics in QFT
Hydrodynamics governs the late time dynamics of non-integrable QFTs at finite temperature.
The simplicity of the hydrodynamic description arises from the fact that most excitations are short-lived at finite temperature, with lifetimes of order of the thermalization time τ th .
This allows for an effective description of the system for times in terms long wavelength fluctuations of the variables characterizing thermal equilibrium, namely temperature and velocity β(x), u µ (x), or their associated densities T 00 (x), T 0i (x).
Additional continuous global symmetries would lead to more conserved quantities. These modes are parametrically long lived because their lifetime grows with their wavelength 1/k.
We define the thermalization length th as the length scale where hydrodynamic modes are no longer parametrically longer-lived than τ th . We will then focus on modes satisfying These time and length scales are parametrically long when the microscopics is weakly coupled, for example th ∼ τ th ∼ β g 4 in (3+1)d gauge theories with coupling g 1 [13] . For strongly interacting QFTs (with speed of sound ∼ 1) one expects th ∼ τ th β, see e.g. [23].
We briefly outline the construction of hydrodynamics for relativistic QFTs, see [24] for a self-contained introduction. Correlation functions for the conserved densities are obtained by solving continuity relations These equations also involve the currents T ij . They can be closed by writing constitutive relations for the currents in a gradient expansion -in the Landau frame one has where P is the pressure, the energy density, ζ, η the bulk and shear viscosities, and the velocity satisfies u µ u µ = −1. We defined the projector ∆ µν ≡ η µν + u µ u ν . The ellipses denote terms that are higher order in derivatives.
Hydrodynamic correlation functions can be found by expanding fields around equilibrium. These correlation functions are therefore obtained after two expansions: a derivative expansion, apparent in (2.4), and an expansion in fields that we will perform below.
where the entropy density is given by s = β( + P ) and the speed of sound c 2 s = ∂P ∂ . This leads to the retarded Green's function where · · · denotes terms that are analytic or subleading when ωτ th , k th 1. The long lived densities T 00 , T 0i carry a sound mode with attenuation rate Γ s = β · ζ + 2(d−1) d η /s, and a diffusive mode with diffusion constant D = β · η/s. Other two-point functions can be obtained from the fluctuation-dissipation theorem: the Wightman Green's function for 1) implies that we are working at small frequencies βω 1). Its Fourier transform will be used below: (2.7) † Neither fixed point describes CFTs in d = 1, where the enhanced symmetries completely fix thermal physics in the thermodynamic limit R/β 1.
For the present purposes it will be useful to understand the constitutive relation (2.4) as an operator equation. Namely, using (2.5) we can write the traceless spatial part as The operator on the left is studied in the IR by expanding it in terms of composites of IR operators T 00 , T 0i with the same quantum numbers (here the quantum number being matched is spin under spatial rotations SO(d)). Correlation functions of both operators will match in the IR. This is routinely done in EFTs, e.g. in chiral perturbation theory where UV operators are represented in the IR in terms of pion degrees of freedom. A similar strategy was followed in [8] where operators with small global charge were represented in terms of operators in the superfluid effective field theory. Although this distinction of UV and IR operators may seem awkward for components of the stress-tensor, we will see that it becomes a useful concept when studying other operators.
In the case at hand, the linear overlap of T ij with IR degrees of freedom implies that the two-point function of T ij will contain the hydrodynamic poles in (2.6) (as can be checked explicitly, see e.g. appendix A in Ref. [23]). At k = 0, the linear term vanishes, but T ij can still decay into a composite of hydrodynamic operators via the second term in (2.8).
It was found [28] (see [29] for a more recent relativistic exposition) that this term leads to 'long-time tails' in the two-point function (2.9) where A ijkl = δ ik δ jl + δ il δ jk − 2 d δ ij δ kl and the integral was computed using (2.7), dropping terms that decay exponentially fast in time. In the first step, we assumed the theory was Gaussian in the hydrodynamic variables, in which case the symmetric Green's functions factorize [30]. This is of course not the case; the same term in (2.8) that leads to longtime tails is responsible for hydrodynamic interactions (classically, these non-linearities are responsible for turbulence in the Navier-Stokes equations). The framework of fluctuating hydrodynamics addresses these interactions. Although hydrodynamics has been understood as a field theory since the work of Euler, the formulation of dissipative hydrodynamics as an EFT is somewhat more recent [31,32,25] and was motivated by the observation of long-time tails in numerics [33], which are now understood as hydrodynamic loops as in Eq. (2.9).
Recent developments in dissipative EFTs for hydrodynamics include [34,35,36,37] (see [38] for a review, and e.g. [39,40,41] for alternative approaches). These constructions allow for a systematic treatment of interactions to arbitrary order in perturbations. Here, we will be working in dimensions where interactions are irrelevant, and will only be interested in the leading hydrodynamic contribution to correlation functions at late times. In this sense we are justified in approximating the action as Gaussian in evaluating (2.9) and in the following.
Systematically accounting for corrections to our results would require knowing the structure of interactions in the effective field theory -this was done for simple diffusion in [42].

Late time correlators from hydrodynamics
How do the thermal correlators of other simple operators behave at late times? The central assumption of thermalization and hydrodynamics is that after short time transients, the only long-lived dynamical degrees of freedom are the densities (2.5). Hence any simple operator will be carried by these densities at late times. For example, any neutral spin-2 operator O µν will have a constitutive relation similar to (2.4) -the stress tensor is only distinguished by the coefficients in its constitutive relation which are fixed in terms of thermodynamic and transport parameters. More generally, consider a traceless symmetric tensor O µ 1 ···µ with even spin (odd spin is mostly similar and is treated in appendix A.2). Its constitutive relation has the schematic form where all terms should be understood to be symmetrized, with traces removed. For some terms there are several possible choices for how the derivatives are distributed -we will be more precise below after determining which terms are most important. The strategy is simply to write all possible composite hydrodynamic operators with the right quantum numbers, in a derivative expansion -we therefore do not explicitly include terms like The λ 0 term in (2.10) was considered in a CFT context in [18] -it is special in that it leads to a non-vanishing equilibrium expectation value O β = 0. However, this term will not always give the leading hydrodynamic contribution to the late time correlators of O, as we show below. In particular this term is forbidden by CPT for odd spin , but odd spin operators still have hydrodynamic tails.
Let us first consider the components of O with zero or one spatial index. Linearizing the constitutive relation (2.10) using Eq. (2.5) shows that these components overlap linearly with hydrodynamic modes: the leading terms are (2.12b) Using (2.6), one finds correlation functions that involve the hydrodynamic poles where we defined the dimensionless entropy density s o ≡ sβ d , and · · · are corrections that are subleading when ωτ th , k th 1.
Now consider correlators involving¯ spatial components of the operator O, with 1 <¯ ≤ .
The constitutive relation (2.10) can again be turned into an operator equation using (2.5)the part that is traceless symmetric in spatial indices is where again all terms should be understood to be symmetrized, with traces removed. There are still several possibilities for how the derivatives act, e.g. in the λ¯ −2 term -this will be specified shortly. We focus on the traceless symmetric part of O i 1 ···i¯ 0···0 , because its traces are related to time components of the operator (e.g.
in turn satisfy similar constitutive relations with fewer indices. This operator matching equation is illustrated in Fig. 1.
We could now proceed by studying the contribution of every operator in (2.14) to the correlator OO . However a simple scaling argument can be used to determine which term in (2.14) is the most relevant: note that (2.7) implies that the densities scale as T 00 ∼ T 0i ∼ k d/2 .
For dimensions d > 2, it is therefore more advantageous to use gradients to build spin. The most relevant operator is the total derivative term λ¯ −1 . We must also keep the term λ¯ ; although it is suppressed when ω ∼ k it can give an enhanced contribution when ω βk 2 , as was shown for¯ = 1 in (2.12b) and (2.13b). Finally, since both of these terms vanish at k = 0, it is also important to keep the most relevant operator that is not a total derivative -when¯ is even this is λ¯ −2 in (2.14) (when¯ is odd, the λ¯ −2 term is a total derivative -this case is treated below). The terms in the constitutive relation (2.14) that give the leading contribution in the hydrodynamic regime ωτ th , k th 1 are therefore λ¯ −2 , λ¯ −1 and λ¯ .
Which term dominates depends on how ω compares to the scales c s k and Dk 2 ∼ Γ s k 2 ; their contributions to the correlator take the form † Here we let O (¯ , ) ≡ O i 1 ···i¯ 0···0 denote components of a spin operator with¯ spatial indices, omitting the corresponding tensor structures; these are treated more carefully in appendix A, see Eq. (A.9). The first line follows from the linear overlaps with the hydrodynamic modes as in (2.13). The second line dominates for k → 0 and comes from a long-time tail contribution to the two-point function from a hydrodynamic loop, as we now explain. The hydrodynamic loop computation is similar to (2.9), with extra gradients acting on the internal legs. Since λ¯ −2 term in (2.14) scales as k d+¯ −2 , one expects a contribution to the two-point function +¯ −2 (note that one must scale ω ∼ k 2 ). The numerical prefactor can be found by performing the loop integral (see appendix A.1 for more details and the tensor structure): where the numerical coefficients a 1 and a 2 are given in (A.8), and were dropped in (2.15). We see that operators with¯ ≥ 2 spatial indices universally decay as 1/t d 2 +¯ −2 in thermalizing QFTs -although this is a straightforward extension of the well-known stress tensor long-time tails (2.9) to operators with higher spin, this result has to our knowledge not appeared † These hydrodynamic contributions also imply that correlators O (¯ , ) O (¯ , ) (x, t) of neutral operators always decay polynomially in real time. Exponential decay of correlators is therefore not a good criterion for thermalization. I thank Erez Berg for discussions on this point. where we have also indicated the subleading corrections O( D 2 k 4 ω 2 ) for small k = 0 (they are computed explicitly in a special case in appendix A.3, where the analytic structure is also discussed). When the number of spatial derivatives¯ is odd, the λ¯ −2 term in the constitutive relation (2.14) is a total derivative, and there is competition between less relevant terms.
Their contribution to the late time correlator can be computed as in the even¯ case -for See appendix A.2 for more details. The first line in (2.15) is then unchanged for¯ odd, but the second line will be given by the Fourier transform of (2.17) instead of (2.16).
In theories with a large number of degrees of freedom such as holographic theories, the suppression in (2.16) by the dimensionless entropy density s o ≡ sβ d ∼ N 2 1 implies that these hydrodynamic tails will only overcome short-time transients ∼ e −t/τ th at times t τ th log s o (the late time limit of correlation functions therefore does not commute with the N → ∞ limit). The stress tensor tails (2.9) were captured in a holographic model in Ref. [43] by computing a graviton loop in the bulk. However certain tails in the holographic correlators of higher-spin operators (2.16) are reproduced more simply, and are direct consequences of large N factorization: consider a holographic model with a single trace scalar φ. In the absence of a φ → −φ symmetry, the scalar will have a thermal expectation value This is achieved in bottom-up holographic models by including a coupling in the bulk between the scalar and the Weyl tensor [44] (see also [45,46] Although s o acts as a loop counting parameter in fluctuating hydrodynamics, we emphasize that the perturbative expansion is controlled even when s o ∼ 1 because hydrodynamic interactions are irrelevant † . In this paper, we do not assume that s o is large.
We focused above on diagonal two-point functions; extending these results to off-diagonal These methods can also be † For example, the quark-gluon plasma has so ∼ 10 [47,48,49].
with t ij ≡ |t i − t j |, and where 'sym.' means symmetrizing † the times t 1 , · · · , t n . When n is odd, the contribution from the sound pole vanishes because the integrand cos n (c s k|t|) oscillates around zero; when n is even the correlator receives an extra contribution from the sound attenuation rate as in (2.16).

The critical dimension d = 2
The results in the previous section apply to any QFT in spatial dimensions d > 2. For d = 2, hydrodynamic interactions are only marginally irrelevant. One manifestation of this is that all terms in the first line of (2.14) have the same scaling. This implies that many terms contribute to the correlator, which however still scales like (2. 16).
An additional subtlety is that the transport parameters D, Γ s now run. For simplicity, let us assume the bulk viscosity ζ = 0 (as is the case for CFTs), so that Γ s = D = βη/s.
The β-function for D is negative [25,24], so that it flows to infinity in the IR ‡ . Indeed, the † In the approach presented here, the correlators are necessarily symmetrized. Correlators with arbitrary time orderings can however still be computed from hydrodynamics using the effective action, see [50,51,52]. ‡ This implies that canonically normalized interactions ∼ 1/D are marginally irrelevant and the theory is 'free' in the IR, in the sense that it is described by regular tree level hydrodynamics like in higher dimensions.
tree-level and one-loop contribution to the Green's function can be found from (2.6) and (2.9) to be (2.20) Interpreting the quantity in brackets D + ∂D ∂ log ω log ω as a running of the diffusion constant one finds [25] where D Λ is the diffusion constant at the scale Λ. In the deep IR It is a striking feature of (2+1)d hydrodynamics that dissipation does not introduce new parameters at the latest times -transport parameters are fixed in terms of the thermody- -however since many other terms in (2.14) contribute to the same order in ω and k when d = 2, we will not attempt to obtain the exact correlator. These logarithmic corrections are negligible for many practical purposes, but have been observed in classical simulations, see e.g. [55,56]. We will mostly ignore logarithmic corrections in applications to CFT data in Sec. 3.

Real time correlators and diffuson cascade
The hydrodynamic correlators in frequency space G(ω, k) obtained above are the ingredients needed for the CFT applications in Sec. 3; the reader interested in these results may therefore directly skip ahead to that section. In this section we take a slight digression to discuss finite temperature QFT correlators in real time. At finite wavevector k, the linearized hydrodynamic correlators (2.7) decay exponentially in time ∼ e −Dk 2 |t| . We will see in this section that even this standard result is drastically affected by hydrodynamic fluctuations, † Dissipation is also tied to thermodynamics in (1+1)d when hydrodynamic fluctuations are relevant [25,26,27], as was recently emphasized in Ref. [53]. ‡ For CFTs, so = bT in the notation of [18]. A free massless scalar has so = 3ζ(3) 2π , so that e −8πso ≈ 5×10 −7 . For the (2 + 1)d Ising model so ≈ 0.459 [54,18] which in a sense are dangerously irrelevant: although they only give small corrections to G(ω, k), they entirely control the leading behavior of G(t, k) at late times.
Let us therefore study the real time thermal correlation function of an operator with¯ spatial indices (2.23) We will take¯ ≥ 2 even for simplicity, but similar results hold for any¯ . Based on the previous section, one expects the polynomial decay of correlation functions (2.16) to still hold for times smaller than the diffusion time 1/Dk 2 of the mode (in this section, we assume for simplicity that D ∼ Γ s ∼ 2 th /τ th ). The small wavevector of the operator in units of the UV cutoff of hydrodynamics k th 1 will allow for a parametric separation between various regimes of the correlator at late times. To find the cross-over time more precisely, we can compare the contributions from the two first terms in Fig. 1 to the two-point function -one finds that the polynomial decay (2.16) holds in the window regime I: . At slightly later times, the correlator is controlled by the linear overlap with the hydrodynamic mode and has the form regime II: The power-law decay therefore plateaus to a constant before starting to decay exponentially around the diffusion time 1/Dk 2 . So far the discussion here mirrors the one in frequency space, see Eq. (2.15). However the result above eventually breaks down at late times. Indeed, consider again the second term in Fig. 1, where the operator decays into two hydrodynamic excitations. Since it is less relevant than the first term, its contribution to the correlator will be more suppressed by 1/t (or k); however its exponential factor is larger: The exponent corresponds to the energy threshold for production of two diffusive fluctuations, which is half that of a single diffusive mode [42]. More generally, the operator O(k, t) can decay into n diffusive modes, distributing its momentum such that each mode carries k = k/n so that the exponential factor becomes e −Dk 2 t n = e − 1 n Dk 2 t . This is the manifestation in real time t of the n-diffuson branch cut, with branch point ω n-diff = − i n Dk 2 . There are similar ω . . . Figure 4: Analytic structure of hydrodynamic correlation functions G R (ω, k). The circles denote hydrodynamics poles ω diff = −iDk 2 and ω sound = ±c s k − i 2 Γ s k 2 , and the crosses denote branch points ω n-diff = − i n Dk 2 and ω n-sound = ±c s k − i 2n Γ s k 2 located at the threshold for production of n hydrodynamic excitations. branch points at the threshold for production of n sound modes ω n-sound = ±c s k − i 2n Γ s k 2 (see appendix A.3). The analytic structure of G(ω, k) is shown in Fig. 4.
For n sufficiently large, the n-diffuson contribution to the correlator has the form † The perturbative expansion in τ th t is asymptotic and this result should therefore only be trusted for n . The largest contribution can be determined by extremizing (2.28) over n. Approximating log n! ∼ n log n and assuming n t τ th d/2 one finds that the largest contribution at time t comes from decay into diffusons. Plugging back into (2.28) produces the correlator regime III: with α ∼ 2d log t τ th ∼ 1. In this regime, operators decay into more and more diffusive excitations, which leads to a stretched exponential decay of correlators.
Although we have focused on the decay of operators with¯ spatial indices in QFTs at finite temperature, similar results apply for two-point functions of generic neutral operators in any diffusive system, including non-integrable spin chains and random unitary circuits with † This expression applies when n is larger than the spin . When n , decaying into one more T0i costs k d/2 but saves a derivative k, so that the suppression is only (k th ) n(d−2) instead of (k th ) nd in (2.28).
conservation laws. In particular Eqs. (2.26) and (2.30) would apply there, after removing the spatial spin dependence k 2¯ −2 → 1. Certain signatures of diffusive tails have been observed numerically in these systems [57,58], and finite k correlators have been studied e.g. in [59], but to our knowledge this diffuson cascade (2.30) and cross-over from e −Dk 2 t to e − √ Dk 2 t has yet to be observed. One issue is that of finite system size, which we discuss below.
In the thermodynamic limit, correlators will decay as (2.30) indefinitely. Indeed, the effective coupling (t/τ th ) d/2 is always smaller than 1/n(t), with n(t) given by (2.29), so that we do not expect a breakdown of the perturbative expansion in hydrodynamic fluctuations, and Eq. (2.30) can be trusted to arbitrarily late times despite the dangerous irrelevance of these fluctuations † .
In a finite volume L d there is a minimal wavelength that the diffusive fluctuations in the loop can carry: k min = 2π L . The correlator will then be controlled by decay of the operator into n max ∼ k/k min diffusive modes with momentum k min , so that at times later than the Thouless time L 2 /D the correlation function has the form regime IV: where s is the entropy density. When t reaches the upper limit, the Green's function is exponentially small ∼ e −S . At this point we expect the correlator to be described by random matrix theory (RMT) [60,61], exhibiting a ramp that levels off to a plateau ‡ , upon averaging (over a few operators with the same quantum numbers, for example). The exponentially small value of the Green's function  [61] is smaller than the time scale above by a factor k/s k th 1. enter in the determination of the thermalization time τ th , and transport parameters such as D. For weakly coupled theories, the early time behavior t τ th can be studied using direct finite temperature perturbation theory or kinetic theory [13] (which can also capture chaos [62]). However it is difficult to observe the regimes I, III and IV directly in a weakly coupled approach, as hydrodynamic fluctuations are not captured by the linearized approximation to the Boltzmann kinetic equation.
We close with a comment on the convergence of the perturbative expansion. The convergence of the hydrodynamic gradient expansion in large N systems (where hydrodynamic interactions can be ignored if one takes the N → ∞ limit first) has been discussed e.g. in [63,64]. Away from the large N limit, one expects that loop effects cause the gradient expansion to be asymptotic, as usual in effective field theories. This is apparent in the n-diffuson contribution to the correlator (2.28), which blows up when n (t/τ th ) d/2 . It would be interesting to understand if processes involving these large numbers of diffusons can be tamed, or Borel resummed. In perturbative QFT, processes involving many particles can also lead to a breakdown of the perturbative expansion, which can however be saved by expanding around a different saddle [65,66] (see [67] for recent developments); diffusive systems are a natural venue to study multiparticle processes, and perhaps apply some of these techniques.

Semiclassical theory of heavy operators in CFTs
We found in the previous section that the late time thermal two-point functions of light neutral operators of any spin are governed by hydrodynamics in generic (thermalizing) CFTs.
Working in the microcanonical ensemble, this implies that off-diagonal heavy-heavy-light OPE coefficients C HH L are universal, at least on average. A priori, the averaging must be done over a microcanonical window of states. However, heavy operators in thermalizing CFTs are expected to look typical, so that much less averaging may be needed in practice.
This expectation is objectified by the ETH Ansatz [19,20,21,10] for the matrix elements of a light local operator O in energy-momentum eigenstates †P µ |H = |H p µ : 1) † More precisely, the energy of the heavy state on the cylinder R × S d is p0 and pi labels the spherical harmonic on the spatial sphere. We will mostly focus on regimes where the sphere can be approximated as (3.11) and comment below), so that pi = ki will denote regular spatial momentum.
where Ω(p) is the density of states at momentum p, and the R O HH behave like independent random variables with unit variance. Averaging Eq. (3.1) over a microcanonical window of heavy operators H, H simply states the equivalence between microcanonical and canonical ensembles; the non-trivial content of Eq. (3.1) is instead that microcanonical averaging is unnecessary: diagonal matrix elements directly produce thermal expectation values, and offdiagonal matrix elements probe out of equilibrium response, for example through symmetric two-point function OO . The appearance of OO in the variance above is required for the Ansatz to reproduce the two-point function [21,23] (note that the Wightman and symmetric two-point functions are approximately equal in the hydrodynamic regime (2.1)). In a CFT, the state-operator correspondence relates these matrix elements to OPE coefficients. For scalar operators (see e.g. [4]) The diagonal part of ETH (3.1) implies that diagonal heavy-heavy-light OPEs are controlled by equilibrium thermodynamics, as found in Ref. [10] (see also [4,68]). In section 3.1 their results are reviewed and extended to operators with spin. In section 3.2 we turn to the off-diagonal part of (3.1), and show how hydrodynamics controls the corresponding OPE coefficients.

Thermodynamics in OPE data
Consider a heavy operator H, with dimension ∆ ≡ ∆ H 1 larger than any other intrinsic number of the CFT (such as measures of the number of degrees of freedom). It will be useful to define the energy density of the state that H creates on the cylinder R × S d of radius R . We can then reach the macroscopic limit by taking R → ∞ while keeping fixed. The diagonal OPE coefficient is fixed by thermodynamics [10] where in the second step we used (3.1), and (3.3) in the last to eliminate the radius. We used (2.10) and (2.11) to express thermal expectation values as where b T = s o = sβ d is the dimensionless entropy density, and is related to the energy Eq. (3.4) can be straightforwardly extended to operators with spin. From Eq. (2.14) we find that the thermal expectation value of light operator of even spin takes the form where we again used scale invariance to write If the heavy operators are still scalars, the OPE coefficient C HHO is parametrized by a single tensor structure [69], which agrees with (3.6) (see Ref. [3] for similar checks in the large charge limit). Now if the heavy operators carry a spin J that is not macroscopic -i.e. J ∆ → 0 in the macroscopic limit -the states they create on the cylinder are homogeneous so that (3.6) still applies. However, many tensor structures can now appear [69], each with their own OPE coefficients. could depend on J and m (we are using SO(3) notation for simplicity, which can be easily generalized to SO(d + 1) for d > 2). Note that we are focusing on diagonal matrix elements here, so that both states have to have the same weight. Comparison with (3.6) shows that the leading answer is in fact independent of J and m as long as J ∆, so that (3.4) still . (3.8) Diagonal OPE coefficients involving heavy operators with macroscopic spin J ∼ ∆ are discussed in section 3.3.

Hydrodynamics in OPE data
We will use (3.1) and the correlators obtained in Sec. 2 to determine OPE coefficients (3.2) in the 'macroscopic' limit R → ∞, with a 'mesoscopic' difference in the dimensions of the in spacetime dimensions d + 1 ≥ 3, keeping the energy density and frequency ω finite.
When the mesoscopic difference in their dimensions is not large The lower bound comes from the fact that the hydrodynamic results will receive corrections from the finite size of the sphere of radius R at the Thouless energy ω ∼ D/R 2 -these could be obtained by generalizing the hydrodynamic correlators of Sec. 2 to the sphere, but we will not attempt to do so here; we expect the singular features that we find in OPE coefficients to be softened in that regime. The upper bound however is a fundamental UV cutoff of hydrodynamics (2.1), assuming τ th ∼ β.

The OPE coefficients C H
The simplest case of three scalar operators J = J = = 0 is somewhat subtle and will be discussed further below. Instead we start with a light spinning operator with ≥ 2, and take J, J 'microscopic', i.e. they are kept fix ∼ 1 in the macroscopic limit R → ∞.

Microscopic spin J, J
Keeping the spin J, J ∼ 1 fixed in the R → ∞ limit implies that the Green's function in (3.1) must be evaluated at spatial wave-vector k = J−J R = 0. In this case we found that the correlator is controlled by hydrodynamic loops: for components with¯ ≥ 2 spatial indices, † A CFT is expected to be weakly coupled when its twist gap γ ≡ min (∆ − J) − d + 1 ≥ 0 is small γ 1.
In this case the thermalization time is parametrically enhanced τ th ∼ β/γ, as can be observed e.g. in the O(N ) model by comparing its thermalization time [70] to its twist gap [71]. Generic CFTs are expected to satisfy γ 1 and τ th ∼ β.
the Green's function is (using (2.15) and (2.11)) (3.13) with α¯ = d 2 +¯ − 2 for¯ even, and α¯ = d +¯ − 3 for¯ odd (see (2.17)). Converting this into an expression for the OPE C H J H J O using (3.1) and (3.2), we see that the hydrodynamic answer predicts a tensor structure (and fixes the correspondig OPE coefficient) for each = 0, 1, . . . , . However if the heavy operators are both scalars, conformal invariance constrains the three-point function up to a single OPE coefficient (the illegal step was to apply ETH (3.1) before accounting for all symmetries). To accommodate the tensor structures obtained from hydrodynamics, it is sufficient to let one of the heavy operators have spin J ≥ -this leads to precisely¯ + 1 tensor structures in agreement with the CFT where σ = sgn(m) = ± denotes the spatial directions ± = x 1 ± ix 2 (x 1 , x 2 are the directions used to define the weight m, i.e. J 12 |H, Jm = |H, Jm m).
Combining Eqs. (3.1), (3.2) and (3.13) therefore gives Here we have replaced the random number with unit variance |R HH | 2 → 1. Strictly speaking this expression for |C HH L | 2 and those below should be thought as average statements, averaged over a few heavy operators H or H . We have taken the density of states at energy . (3.16) Eliminating the radius R in (3.15) gives (dropping numerical factors) We will not attempt to control subleading extensions to the Cardy formula (3.16), and therefore will not comment on the subexponential dependence on ∆. However, we attract the reader's attention to the non-analytic dependence on ∆ − ∆ , coming from hydrodynamic † In the notation of Ref. [69], the tensor structures in the sum are fluctuations. Corrections to the leading result are shown in the square brackets, and come from less relevant terms in hydrodynamics and finite volume corrections: , (3.18) both of which are parametrically small in the regime (3.11).
If both J, J ≥ 1, there may be more tensor structures allowed by conformal invariance than needed -the claim of thermality is that as in (3.8) the leading OPE coefficients will not depend on these extra indices.

Mesoscopic spin J, J
Let us now extend to heavy operators with 'mesoscopic' spin. More precisely, we want the spins to be non-macroscopic (so that the state on the sphere remains homogeneous), and the difference in spins to be mesoscopic, or (3.19) in the limit R → ∞. Let us first consider a light operator O with spin = 0. In Sec. 2 we found that its thermal expectation value leads to (see (2.13a)) Conformal invariance allows for many tensor structures for the three-point function and where the subexponential dependence on ∆ (which is degenerate with logarithmic corrections to S(∆)) was packaged in α ∝ . These OPE coefficients feature a 'resonance' at the sound mode ∆ − ∆ = ± 1 We are now ready to turn to the general case of the heavy-heavy-light OPE coefficient of three spinning operators. The hydrodynamic prediction was given in Eq. (2.15). To match with OPE coefficients we will need the precise index structure, which can be conveniently packaged by using the index-free notation with z 2 = z 2 = 0. In this notation, the full index structure of (2.15) is given in appendix A (see Eq. (A.9)). For a CFT (A.9) becomes with r = 0, 1, . . . ,¯ . The leading hydrodynamic result (3.24) only contains these structures for r =¯ ,¯ − 1 and 0; coefficients of the structures for other values of r will be controlled by subleading hydrodynamic tails (in d = 2 these additional tails are only log suppressed compared to the leading ones, see Sec. 2.2). One therefore obtains ¯ =0 (¯ +1) = 1 2 ( +1)( +2) structures. When J, J ≥ , Ref. [69] showed that the CFT three-point H|O|H function contains more structures: there are 1 2 ( + 1) ( + 2) min(J, J ) + 1 − 3 structures, which in their notation take the form where a, a , b run over all integers such that the powers above are positive. The hydrodynamic OPE coefficients (3.24) only depend on the contractions between the light operator and the heavy ones, hence on a, a but not on b (see Fig. 6). Since a, a = 0, 1, . . . , satisfy and with the numerical factors a d = 2(d−1) . Subleading corrections to these results are similar to those in Eq. (3.17).

Macroscopic spin
Let us now briefly comment on heavy operators with macroscopic spin Macroscopic spin has been treated in an EFT approach for large charge in CFTs with a U (1) symmetry [9,22]. It was found there that in the regime (3.28), the superfluid state forms a vortex lattice, such that the coarse-grained superfluid velocity is equal to that of a rotating body with angular momentum J. For a normal fluid, one expects a similar stationary solution to the Navier-Stokes equations † . Let us work in d = 2 spatial dimensions for simplicity, and search for a velocity profile u µ ≡ (u 0 , u θ , u φ ) = (1 + v 2 φ ) 1/2 , 0, v φ with an azimuthal velocity that only depends on the polar angle v φ = v φ (θ). Now a typical state with angular momentum on the sphere will equilibrate (preserving its angular momentum) to an equilibrium velocity profile v φ (θ) that does not dissipate; in particular it must be annihilated by the shear viscosity term in (2.4) -this leads to a differential equation which can be solved for v φ (θ). Energy eigenstates created by heavy operators are expected to look thermal and should have this velocity profile. The ideal stress tensor can then be obtained by imposing conservation

Global symmetries
It is straightforward to extend the results above to QFTs and CFTs with an internal symmetry group G; this section deals with the simplest example G = U (1). The additional Ward identity ∂ µ J µ = 0 protects a new slow excitation -charge density J 0 -whose fluctuations will give additional contributions to late time correlators.
A background chemical potential µ can be introduced for the internal symmetry. In  Fig. 2. This will allow us to connect to recent work on the large charge limit of CFTs [6,7,8,9,3,22].

Hydrodynamics of a charged fluid
The conservation laws ∂ µ T µν = 0, ∂ µ J µ = 0 must be supplemented with constitutive relations for the currents. In the Landau frame and up to first order in derivatives, the constitutive relation for the stress-tensor is still given by Eq. (2.4) and that of the U (1) † Dissipative superfluid hydrodynamics also describes 2+1d theories at finite temperatures 0 < T < T BKT , where strictly there is no spontaneous symmetry breaking; the protection of the long-lived superfluid phase can however be understood without reference to symmetry breaking [72].
current is [24] Three new parameters were introduced: ρ, κ and χ T . † These, along with those appearing in (2.4), are functions of both µ and β. Consistency with thermodynamics fixes ρ in terms of the equation of state ρ = ∂P/∂µ, and imposes κ ≥ 0 and χ T = 0. ‡ Hydrodynamic correlators can again be obtained by expanding around equilibrium (2.5) with µ(x) = µ + δµ(x). If we first take the background chemical potential to vanish µ = 0, then the background charge density ρ vanishes by CPT and we see directly from (4.1) that there is no mixing at the linear level between the new hydrodynamic degree of freedom δµ and the ones considered previously δu µ , δβ, at least to this order in derivatives. The stress-tensor correlator (2.6) is therefore unchanged, and the current correlator is given by where χ ≡ ∂ρ/∂µ is the charge susceptibility, and D c ≡ κβ/χ the charge diffusion constant.
It is easy to see that the new hydrodynamic field δµ does not allow the construction of more relevant operators -the results from Section 2 are therefore largely unchanged -except for odd spin operators with¯ = 0 or 1 spatial indices. The reason is that for these cases we found in Sec. 2 that the dominant hydrodynamic contributions to the correlators (2.13) involve the term λ 0 in (2.10), which was forbidden by CPT for odd-spin operators (see appendix A.2.1). However, thanks to the conserved U (1) charge this term is now allowed for odd spin as well where λ 0 (µ, β) is an odd function of µ by CPT. Expanding λ 0 in δµ, one finds that components with¯ = 0 spatial indices can overlap linearly with the density, so that (4.2) implies (4.5) † Another commonly used notation for the conductivity is σQ ≡ κβ. ‡ Note that χ T is only forbidden because Jµ is conserved. Generic non-conserved spin-1 operators will have terms like χ T in their constitutive relation. In holographic models, along the lines of Ref. [44], these could come from coupling a massive gauge field in the bulk to the Weyl tensor, e.g. through A µ ∂µC 2 Weyl . and components with¯ = 1 spatial indices are controlled by a hydrodynamic loop at k = 0 (4.6) A special case is the correlator of the current operator itself O µ = j µ . Then λ 0 = ρ so ∂λ 0 /∂µ = χ and (4.6) reproduces known results [28,29]. This correlator with O µ = j µ is the one that led to the original discovery of long-time tails [33].

Turning on a background µ = 0
A background chemical potential will allow the longitudinal hydrodynamic modes j 0 , T 00 and ∂ i T 0i to mix (the transverse sector is unaffected and still given by the second term in (2.6)). The longitudinal sector will still contain a diffusive mode and a sound mode, but these will be carried by linear combinations of j 0 and T 00 , see e.g. [24]. The correlators (2.15) of neutral operators therefore do not change qualitatively: the functional dependence on ω, k is unchanged, but the thermodynamic and transport factors are more complicated.
One exception is again for operators of odd spin . For example components with¯ = 1 spatial indices can now overlap linearly with hydrodynamic modes, even at k = 0. Indeed, λ(µ, β) in (4.4) is now expanded around µ = 0 so that the constitutive relation O i0···0 = λ 0 u i + · · · has the same form as (2.12b), and the two- given by (2.13b). More generally, in charged hydrodynamics at finite density, the results of Sec. 2 hold for both even and odd spin , because of the absence of any CPT constraint.

Dissipative superfluids
The hydrodynamic theory of relativistic, dissipative superfluids was thoroughly studied in Refs. [73,74]. Compared to normal charged fluids, superfluids contain an additional slow hydrodynamic degree of freedom carried by the Goldstone field φ that non-linearly realizes the internal U (1) symmetry. Here we will focus on conformal superfluids, and will not give an expectation value to the superfluid velocity, i.e. ∂ i φ = 0. This velocity can be thought of as the charge density associated with an emergent higher-form symmetry [72] -since the symmetry is emergent, heavy CFT operators creating superfluid states are not labeled by their representations under it. Working to linear order in the superfluid velocity, there is only one new thermodynamic parameter compared to (4.1) -the superfluid stiffness ρ s -and one new dissipative parameter ζ 3 [73] (see also [74]) : T µν = u µ u ν + P ∆ µν + 2ρ s µ s n (µ u ν) + ρ s µ s n µ n ν − ησ µν + · · · , (4.7a) where n µ ≡ ∆ µν ∂ ν φ/µ s and µ s ≡ −u µ ∂ µ φ. The projection ∆ µν ≡ η µν + u µ u ν , and σ µν is the shear viscosity tensor appearing in (2.4). The thermodynamic parameters , P, ρ s , ρ n and dissipative parameters η, κ, ζ 3 are inputs in the hydrodynamic treatment -however when the dissipative superfluid is obtained by heating up a T = 0 conformal superfluid (as in Fig. 2), these can all be expressed in terms of a single EFT parameter at low temperature, One new feature is that operators with finite charge q ∈ Z under the U (1) symmetry can now be matched in the IR using the Goldstone phase where as in Fig. 1  Expanding φ = δφ − µt shows that real time correlators contain an extra factor of e −iqµt , so that in frequency space one finds where the right-hand side simply refers to the general result (2.15) for neutral operators, evaluated at frequency ω − qµ. One may worry that the hydrodynamic features in this correlator appear at frequencies above the hydrodynamic cutoff ω ∼ qµ 1/τ th -however this is simply because the operator carries a phase e −iqµt which translates hydrodynamic features usually at ω ∼ 0 to ω ∼ qµ in the correlator of operators of charge q. We will see in the CFT application in Sec. 4.3 that ω − qµ measures the difference in dimensions ∆ − ∆ of the heavy operators (as did ω for neutral operators, see (3.9)).
Of course, superfluids at finite temperature also have well known static properties which control equal time correlators. For example terms like the first term in (4.8) will lead to (4.10) For d = 2, where the equal-time phase correlator φ(x)φ = 1 2πρs log |x| (where ρ s is the stiffness in (4.7)) one has instead (4.11) These also provide EFT constraints on the CFT data involving heavy charged operators that create a superfluid state.

Dissipative superfluids from the EFT
It is natural to expect that if a CFT exhibits a superfluid phase, this phase will be connected to a T = 0 superfluid, as in Fig. 2. At T = 0, a superfluid EFT describes the physics up to a cutoff, which in the case of a CFT must be proportional to the chemical potential µ [75,76,6,8]. Dissipative hydrodynamics can be seen to emerge from the EFT at finite temperatures 0 < T µ; in other words, all the thermodynamic parameters and dissipative parameters in the previous section can be computed in terms of the EFT parameters. Let us illustrate this to leading order in gradients, where the EFT is simply [76] with |∂φ| ≡ −∂ µ φ∂ µ φ. The dimensionless constant c 1 is non-universal and depends on the underlying CFT. The U (1) current is Expanding around the saddle φ = µt + π, one finds a zero temperature superfluid density (4.14) Thermodynamics δ = µδρ then fixes the zero temperature energy density as can be checked by computing the stress tensor directly from (4.12). The pressure for a CFT is given by P = /d. From the CFT perspective, c 1 can be defined by Eq. (4.14) or (4.15) and can be viewed as CFT data on a similar footing as the thermal expectation value of the stress tensor, b T in (3.5).
The action (4.12) can be expanded around the saddle φ = µt + π with π c ≡ c 1 µ d−1 π. Only the schematic form of the interactions S int will be needed -∂π c symbolizes either time or space derivatives and numerical factors have been dropped. The strong coupling scale of the EFT is given by the energy density Λ sc ∼ 1/(d+1) ∼ c µ.
Let us start with the thermodynamics, which to leading order can be studied from the free part Euclidean version of the action (4.16a). The simplest finite temperature quantity to compute is the entropy density, which can be obtained from the free energy with c s = 1/ √ d, from which we can obtain the dimensionless entropy density Terms that are higher order in gradients or field in the action (4.12) and (4.16) lead to corrections to the expressions above that are suppressed by powers of T /µ. The normal density is slightly more subtle: it comes from taking the thermal expectation value of nonlinear terms in the current (4.13), and showing the disalignment between the current and the expectation value of ∂ µ φ [72]. One finds with s given by (4.18). The non-relativistic limit c s 1 of this expression is well known [77].
A similar expression has appeared in a holographic context recently [78] -we see here that it is a universal prediction of the EFT. For a CFT, c 2 s = 1/d. Furthermore, the emergence of hydrodynamics at finite temperature leads to an additional sound mode (second sound)since the EFT is to leading order a free scalar and hence scale invariant at low energies, the speed of second sound is itself related to that of first sound as c 2 s,2 = c 2 s /d at low temperatures T µ (see e.g. [77,79]).
Finally, the dissipative parameters η, κ, ζ 3 that appeared in the previous section can also be computed from the EFT (4.16) by treating the weakly coupled phonons with kinetic (Boltzmann) theory. This was first done for non-conformal superfluids (which have two additional viscosities ζ 1 , ζ 2 ) in the non-relativistic limit in Ref. [77], see Ref. [80] for a more recent and detailed discussion. The calculation is quite lengthy so we only sketch it here, focusing on the shear viscosity η for illustration. The phonon differential cross section can be computed at tree level from the cubic and quartic terms in (4.16) (see e.g. [81]), the diagrams in Fig. 7 where is the energy density (4.15) and p symbolizes dependence on the individual phonon momenta p i , i = 1, 2, 3, 4. The dependence on the individual momenta can be important, in particular the total cross section σ diverges because of small angle scattering [77,81].
This divergence is regulated by more irrelevant terms in the action (4.12), so that the total cross-section is less suppressed by the cutoff than Eq. (4.20) suggests [77]. However it is large angle scattering that controls the shear viscosity [13], so that the naive expression (4.20) is sufficient for our parametric estimate. One can now estimate the thermalization time from the thermally averaged cross-section The thermalization time is large τ th β because the phonons are weakly coupled. The shear viscosity can then be estimated as The viscosity diverges rapidly as T → 0 because of the long thermalization time (4.21) of the superfluid.
It is interesting to contrast these results to holographic superfluids [82,83]. Because these theories have a large O(N 2 ) number of degrees of freedom, the superfluid sector only gives small O(1) corrections to thermodynamic quantities such as the entropy density s.
However, transport is more sensitive to the presence of the weakly coupled superfluid sector.
The holographic value of the low temperature shear viscosity .

(4.24)
A similar conclusion holds for more general hyperscaling-violating Lifschitz geometries [84,78], with a different exponent in (4.24). It would be interesting to understand if this non-commutativity of the T → 0 and N → ∞ limits signals a more important breakdown of low temperature finite density holographic solutions (such as extremal black holes) due to quantum effects [85].

Implications for heavy CFT operators with macroscopic charge
The ETH Ansatz (3.1) is slightly modified for systems with additional symmetries. For the case of an internal U (1) symmetry with generatorQ, the extension can be simply obtained by using the HamiltonianĤ →Ĥ − µQ without the need of using the grand canonical ensemble explicitly. One then obtains, for a few-body operator O q of U (1) charge q, In a superfluid phase, we found in Eq. (4.9) that the correlators of light charged operators O q are also controlled by hydrodynamics. Therefore, when the state created by the heavy operator H Q,J is a finite temperature superfluid we can use (4.25) to obtain hydrodynamic predictions for OPE coefficients of light charged operators. We find that the results in Sec. 3 for neutral operators (q = 0) are essentially unchanged, but now also hold for charged operators (with the obvious constraint of charge conservation). For example (3.17) becomes (4.26) the only difference with (3.17) being that this also holds for q = 0 and the relevant transport parameterη o is not simply the shear viscosity but a combination of the superfluid dissipative parameters η, κ and ζ 3 from Sec. 4 A natural possibility is that this state is a superfluid [6]. The superfluid EFT then predicts both the spectrum of low-lying operators, and OPE coefficients between these operators and light CFT operators [8,3]. As one moves away from the edge ∆ ≥ ∆ min (Q), the spectrum becomes dense and the many phonon state eventually start to look thermal. Notice that here the thermalization time is large (4.21), because the original EFT is weakly coupled.
The dimension of operators near the edge can be written as where δ 1 is related to the temperature by (this relation follows from Eq. (4.40) derived in the following section) with s sflu o given by (4.18). This implies that the hydrodynamic window (3.11) is parametrically smaller close to the edge of the spectrum (4.30)

Phase transitions in the spectrum
The equation of state of a CFT at finite µ and β is no longer fixed by scale invariance, but can depend on the dimensionless reduced chemical potential α ≡ βµ . (4.31) In the previous sections, we explored the hydrodynamic descriptions pertaining to two natural phases of CFTs at finite density -the superfluid phase that is expected for α 1 and normal phase for α 1 -and determined how hydrodynamics controls some of the CFT data. These phases should be separated by a phase transition. In this section, we explore how the nontrivial thermodynamic properties of the transition control the data of the underlying CFT, and leave for future work a hydrodynamic treatment of the system near the phase transition (this would require incorporating long-lived critical fluctuations, see e.g. Refs. [86,87]). In this sense, this section extends the work of Ref. [10], where thermodynamics was seen to control some of the CFT data, to situations where the thermodynamic equation of state and corresponding phase structure are non-trivial.
Expectation values of the currents now take the form reduce the equation of state to a single function of one variable, which we could take for example to be s o (α). However when studying the operator spectrum in a CFT, it is most convenient to work in the microcanonical ensemble and to think instead of α (or µ) and β as functions of the densities, say and ρ. In particular, it will be convenient to study a slice of Fig. 2 at fixed ∆ 1, i.e. fixed energy density , and vary charge. Since is fixed, we can use it to define a dimensionless charge density and temperature . The potentials α(n) andβ(n) are dimensionless functions of the dimensionless charge density n. The thermodynamic relations (4.33) imply that these functions satisfy so that only one function is independent, sayβ(n), and can be thought of as the equation of state characterizing the thermodynamic properties of the CFT. Thermodynamic stability further implies so that both α andβ are positive, monotonically increasing functions of n.
The asymptotic properties of the equation of state can be related to familiar parameters of the CFT. For example, as n → 0 one has The first term simply comes from (3.5), and the subleading term follows from (4.33) and (4.35) and features the dimensionless charge susceptibility (χ can also be expressed as a thermal 2-point function of the current at zero chemical potential). The monotonicity ofβ (4.36) for n 1 is equivalent to χ o ≥ 0.
The equation of state is also fixed in the opposite limit if we assume, following Ref. [6], that the state at is a zero-temperature superfluid † . Using again the thermodynamic identities (4.35), one finds that the equation of state near the zero-temperature superfluid takes the form (as n → n max ) , (4.40) where s sflu o is given by (4.18). Note that the two asymptotic behaviors (4.37) and (4.  across a U (1)-restoring thermal phase transition † . That the transition is in the mean-field universality class in this case [82], with η = 0 and ν = 1/2, is likely an artefact of large N ; mean-field critical exponents are not expected for generic CFTs.
Consider for example a (3+1)d CFT with a global U (1) symmetry, and assume following Ref. [6] that the lightest operator of charge Q creates a superfluid state when n = n max .
When n is decreased past n c , the symmetry is restored and we expect the transition to be in the 3d Wilson-Fisher universality class, with ν 0.672 and η 0.038 ‡ . These exponents control correlators near or at the critical n c , which like the hydrodynamic long-time tails will lead to predictions for some of the CFT data. For example the anomalous dimension η will control the equal-time correlator of light, charged operators (4.41) The correlation length critical exponent can be obtained from the vanishing of the thermal † See Fig. 3 , (4.44) and the equation of state β(n) is very smooth, with an essential singularity at n = n c .
The phase diagram can be considerable enriched by considering operators with spin J = j∆, with 0 ≤ j ≤ 1 (still in the ∆ → ∞ limit) † . The corresponding states at zero temperature, i.e. keeping the charge density as large as possible n = n max ( , j), were studied in d = 2 and d = 3 in Refs. [9,22], where it was found that the angular momentum of the state is carried by different objects (on top of the superfluid background) depending on j: vortex-antivortex pair (4.45b) As j = J ∆ → 1, the superfluid EFT breaks down and the spectrum is instead governed by the light-cone bootstrap. Departing from the manifold of maximal charge at fixed dimension and spin, these 'phases' will be embedded in a larger phase diagram with finite temperature phases. At large enough temperatures, the U (1) symmetry will be restored, and the vortex lattice will melt. In Fig. 9, we show a tentative operator spectrum 'phase diagram' for heavy operators ∆ 1 of a CFTs with a global U (1) symmetry. also turn into a first order transition at larger j, as is observed in holographic superfluids [83].
Similar phase diagrams have been observed in liquid helium [92], Bose-Einstein condensates [93], thin film superconductors [94], and quantum Hall systems; in the last spin per charge is mapped to the filling fraction ν = J/Q ∼ ∆ 1 d+1 j/n (see e.g. [93,95,96]). Comparison with these systems suggest a number of possible exotic features in the phase diagram in Fig. 9. For example in (2+1)d, the spinning operators studied in [9] can lead to states with opposite vorticity on each poles, and vanishing vorticity along the equator. Gapless edge states are then expected to live along the equator. Since these are supported in (1+1)d, their hydrodynamic interactions are relevant and dissipation anomalous [25,26,27,53].
The CFT spectrum may also probe the melting of the vortex lattice in Fig. 9. In (2+1)d this transition is infinite order like BKT (4.44), but with different exponents [97]. Finally dynamical response and transport near the equilibrium critical point is also singular [86].
We leave a more thorough exploration of this phase diagram for future work.

Conclusion
We showed that hydrodynamics controls a large portion of the CFT data, namely OPE coefficients of any two heavy operators close enough in dimensions (see (1.5) [108,109,110,111,112], see in particular [113,114,115] for discussions on the off-diagonal part of ETH in this context. Far from equilibrium techniques and turbulence may also be useful in higher dimensions to determine OPE coefficients C HH L away from the hydrodynamic linear response regime regime (1.5), e.g. to study There are a number of possible interesting extensions, which we leave for future work.
We list a few below: • It should be possible to extend our results to CFTs with anomalies or non-trivial current algebras by studying hydrodynamics with anomalous Ward identities [116,117,118,119].
• We have mostly focused on local operators. Certain nonlocal operators, for example in gauge theory, have signatures in the corresponding hydrodynamic theories as higherform charges [120,121].
• Operators that are odd under parity (or inversion) can be considered as well, with hydrodynamic tails that depend non-trivially on dimensionality. One can also study heavy operators in CFTs without inversion symmetry, using parity-violating hydrodynamics e.g. in 2+1d [122].
• Boost symmetry plays only a minor role in Sec. 2 -hydrodynamic tails control late time correlators in non Lorentz-invariant QFTs as well. The CFT implications in Sec. 3 rely on a state-operator map. We expect similar results to exist in non-relativistic CFTs (with Schrödinger symmetry), since these also enjoy an operator-state correspondence [123]. The large charge bootstrap has already been extended in this direction [124,125].
The present work revealed hydrodynamic constraints on CFTs. It is our hope that the favor may one day be returned, with techniques such as crossing and unitarity leading to constraints on dynamics in thermalizing CFTs, e.g. in the form of bounds on transport and thermalization [17,70,126,127,23].
It would also be interesting to explore if the novel features in late time thermal correlators discussed here have implications for cosmology, where thermal physics enters both in the thermal desription of de Sitter space and through the actual temperature of the universe.
We end with an amusing observation: since reflection positive Euclidean CFTs can be continued to a unitary Lorentzian CFTs [128,129], the equilibrium properties of certain statistical mechanical systems at their critical point know about hydrodynamics in one lower dimension! †

A Detailed hydrodynamic correlators
Correlators involving many indices can be treated by using an index free notation (see e.g. [69]). Consider a spin-operator O µ 1 ···µ ; its elements involving¯ spatial components can be packaged as where z lives in d-dimensional space (not d + 1-dimensional spacetime), and satisfies z 2 = 0.
This projects on the spatially traceless part (spatial traces are related to components with more time indices since η µ 1 µ 2 O µ 1 µ 2 ···µ = 0). We are interested in thermal 2-point functions in the hydrodynamic regime ωτ th , k th 1.

A.1 Hydrodynamic loop computation
When k = 0, we found in Sec. 2 that correlators (A.2) are dominated by a hydrodynamic loop. For¯ even this comes from the following term in the constitutive relation (2.14) (q · z)¯ −2 (q · z )¯ −2 (q · z)(q · z ) q 2 cos(c s |q||t|)e − 1 2 Γsq 2 |t| + z · z − (q · z)(q · z ) q 2 e −Dq 2 |t| 2 , (A.4) with q ≡ d d q (2π) q . Terms in the integrand that oscillate with q lead to exponentially decaying terms ∼ e −|t|/τ th (see e.g. [24]). Dropping these gives (A.5) forbidden by CPT † . Higher derivative terms in constitutive relations are also constrained by CPT (see e.g. [37]), however these constraints allow for all the λ i in (2.10) as long as O is not itself a conserved current.
The result (2.15) (or more precisely (A.9)) therefore holds for odd and¯ even as long as λ 0 is not involved, i.e. it holds for¯ ≥ 4 even. Where λ 0 is involved, it is replaced by a subleading hydrodynamic tail. We detail here the cases¯ = 0, 1, 2 (¯ ≥ 3 odd will be treated later). When¯ = 2, the first two lines in (A.9) are unchanged since they do not involve λ 0 . The second line came from a hydrodynamic loop through λ 0 , the most relevant term when k = 0 in the constitutive relation is now and scales as k(k d/2 ) 2 = k d+1 (which is indeed less relevant than the forbidden λ 0 term ∼ k d+¯ −2 = k d ). This leads to a long-time tail contribution to the correlator i T 00 . The resulting correlator will still take the form (2.13), but with additional ω or k suppression.
Since the loop contributions will only dominate terms in the first line when k → 0, we have dropped total derivative terms such as λ¯ −2 . The two most relevant loop contributions come from λ ¯ −1 (1-loop) and λ¯ −3 (2-loop). These terms scale as λ ¯ −1 : k¯ −1+d , λ¯ −3 : k¯ −3+ 3d 2 , (A. 13) so that λ¯ −3 dominates for spatial dimensions d ≤ 4 (and λ ¯ −1 dominates in higher dimensions). The leading correlator then behaves as (A.14) One final special case is when¯ = = 3. Then as shown in A.2.1, λ¯ −3 = λ 0 is forbidden by CPT, so that the top line is replaced by a subleading hydrodynamic tail.

A.3 Subleading tails
In Sec.2, the leading hydrodynamic contribution to correlators O (¯ , ) O (¯ , ) (ω, k) were found by matching the operators O (¯ , ) to composite hydrodynamic operators as in Fig. 1, and then Gaussian factorizing the hydrodynamic fields. Gaussian factorization however only holds at the lowest energies (or smallest ω, k) and irrelevant interactions give subleading corrections to the correlators. These corrections can be captured systematically by using a dissipative effective field theory for fluctuating hydrodynamics, see e.g. [31,25,36]. This was performed to next to leading order in [42] for simple diffusion.
In this appendix, we will illustrate the structure of these subleading corrections with a specific example. We will do so without using the full effective action, and therefore miss certain subleading contributions to the correlator. However, important qualitative features of the answer (such as the analytic structure) will be captured. where · · · denotes higher derivative terms that we will ignore. The two point function O (¯ , ) O (¯ , ) will receive a single tree-level contribution from the linear term δT 00 δT 00 , which is given in (2.13a). It will receive 1-loop corrections from δT 00 δT 00 , δT 00 δT 2 00 and δT 2 00 δT 2 00 . The first two come from interactions in the action and will not be captured here -we will focus on the last. At leading order it can be factorized δT 2 00 δT 2 00 (ω, k) 2 d d xdt e ik·x−iωt T 00 T 00 (x, t) 2 .
(A. 16) The hydrodynamic correlator appearing in the integrand can be obtained from (2.7) T 00 T 00 (x, t) s 2β 2 c 2 e −(x+c|t|) 2 /2Γs|t| + e −(x−c|t|) 2 /2Γs|t| (2πΓ s |t|) d/2 , (A. 17) so that performing the integral yields δT 2 00 δT 2 00 (ω, k) (A.18) The quantity in square brackets is formally divergent -the divergence can be treated in dimensional regularization by expanding around integer d and throwing away the divergent piece (this UV divergence can be absorbed in the bare transport parameters [42]). This gives which should be contrasted to the pole in the tree-level part of the correlator (2.13a), at ω = ±c s k − i 2 Γ s k 2 . This is the sound analog of the two-diffuson branch cut at ω = − i 2 Dk 2 found in [42]. The analytic structure of hydrodynamic correlators is shown in Fig. 4. For generic ω ∼ k, (A.18) is suppressed compared to the tree-level contribution (2.13a) -however it dominates as k → 0 where one finds O (0, ) O (0, ) (ω, k = 0) ∼ (∂ 2 β λ 0 ) 2 ω d 2 −1 . Higher-loop corrections will lead to additional branch points at the threshold for production of n diffusons ω = − i n Dk 2 or n sound modes ω = ±c s k − i 2n Γ s k 2 , but the discontinuities across the cuts are increasingly suppressed at small ω and k. However, the Fourier transform ω → t picks up these non-analyticities, and the leading behavior of G(t, k) is controlled by multi-diffuson decay at late time, as shown in Sec. 2.3.