The timbre of Hawking gravitons: an effective description of energy transport from holography

Planar black holes in AdS, which are holographically dual to compressible relativistic fluids, have a long-lived phonon mode that captures the physics of attenuated sound propagation and transports energy in the plasma. We describe the open effective field theory of this fluctuating phonon degree of freedom. The dynamics of the phonon is encoded in a single scalar field whose gravitational coupling has non-trivial spatial momentum dependence. This description fits neatly into the paradigm of classifying gravitational modes by their Markovianity index, depending on whether they are long-lived. The sound scalar is a non-Markovian field with index (3-d) for a d-dimensional fluid. We reproduce (and extend) the dispersion relation of the holographic sound mode to quartic order in derivatives, constructing in the process the effective field theory governing its attenuated dynamics and associated stochastic fluctuations. We also remark on the presence of additional spatially homogeneous zero modes in the gravitational problem, which remain disconnected from the phonon Goldstone mode.


Introduction
Perturbations of asymptotically AdS black holes relate to the thermal stochastic dynamics of the dual boundary conformal field theory plasma [1][2][3][4]. The classical decay of perturbations via quasinormal modes is accompanied by the thermal and quantum fluctuations, manifested in the form of Hawking quanta. The physics of such dissipative systems should be encoded in a real-time open effective field theory (EFT). Such a theory can be constructed as the dynamics induced onto a probe system coupled to the thermal environment after integrating out the plasma degrees of freedom. Recently, such open EFTs were constructed from holography: for instance, [5] describes the non-linear dynamics of quasinormal modes associated with probes that have short lifetimes of order the thermal scale, and [6,7] (cf., also [8][9][10]) describe the physics of diffusion which originates from the fact that the conserved currents of plasma have long-lived modes. Other related work includes [11][12][13][14][15]. Progress in this direction was spurred by [8] who gave a simple geometric proposal for computing real-time, Schwinger-Keldysh ordered observables in holography (cf., [16][17][18][19] for earlier work in this direction).
In the previous works [6,7] the focus was on constructing an open effective field theory of diffusive modes focusing on momentum diffusion. This occurs through shear modes, which carry transverse momentum and shear the fluid elements (hence their name) without compressing the plasma as they propagate. Hence, the shear hydrodynamic mode is present in compressible as well as incompressible fluids alike. The shearing causes transverse viscous drag, resulting in the diffusion of the transverse momentum. These shear modes are long-lived and diffuse slowly through the medium.
On the other hand, compressible fluids have an additional degree of freedom: the sound mode. Sound travels by carrying a longitudinal momentum, applying pressure on the fluid elements, which in turn results in a wave of compressions and rarefactions. As the fluid gets compressed, there is a local change in pressure and energy density, unlike in the case of shear modes (where the local pressure and energy density remain unperturbed). We remind the reader that the relativistic fluids are always compressible since incompressibility requires an instantaneous transmission of pressure which is forbidden within special relativity. Thus, relativistic fluids always have sound modes. The physics of shear and sound are qualitatively different. While shear modes are diffusive and obey parabolic PDEs, sound modes are oscillatory and obey hyperbolic PDEs. The difference owes to the fact that a fluid at rest already has an energy density and pressure. A perturbation over this background results in the sound mode. 1 In contrast, a fluid at rest has no momentum density, so the transverse momentum diffuses, resulting in the shear mode.
We have described until now the sound and the shear mode as being governed by second order PDEs. In a realistic system, higher derivative corrections appear and the effects of the thermal fluctuations need to be considered, resulting in a higher derivative stochastic PDE. In the case of the sound mode, higher derivative corrections describe the sound relaxation whereas fluctuations describe the noise background. More precisely, the sound mode (and the shear mode) should be thought in terms of an open effective theory incorporating both fluctuation and dissipation. Such a theory is obtained by integrating out the fast modes in the plasma. This is more easily said than done: one first needs to systematically separate out the fast modes from the slow modes. In weakly coupled theories, we need to further deal with the fact that sound relaxes over non-perturbatively long time scales. We will sidestep these issues by considering a strongly coupled CFT plasma and study it using holography. Sound and shear modes in holography appear as linearized metric perturbations over planar Schwarzschild-AdS d+1 black holes [3,4]. Non-linearly, the plasma dynamics is captured by the fluid/gravity correspondence [20,21] (which as constructed to date however only captures the dissipative but not the fluctuation effects).
The holographic set-up using gravitational dynamics to model the dynamics of a strongly coupled plasma has some important features highlighted in [6]. Firstly, it makes clear that the dynamics of short-lived and long-lived modes, dubbed Markovian and non-Markovian, respectively, are qualitatively different. At the linearized level, each such mode is described by a scalar field propagating in the AdS black hole background. These fields were called designer scalars as they are non-minimally coupled to the gravitational background, with the coupling modeled as a (auxiliary) background dilaton. This dilatonic coupling modulates the interaction the field has in the radial direction. Heuristically, the short-lived modes are repelled from the boundary of the spacetime, while the long-lived modes are floppy and have a large wavefunction support near the boundary, which prevents them from decaying away rapidly. Technically, the Markovian nature is captured by the asymptotic behaviour of the dilaton and can be encoded in a single number, the Markovianity index, M.
The Markovian nature of a bulk field is characterized by the boundary conditions imposed on it in order to compute correlation functions via the GKPW dictionary.  This definition mildly updates the definition given in [6]. The main distinction is that fields in the window M ∈ (−1, 1) can either be Markovian or non-Markovian, depending on the boundary conditions imposed. This window is similar to the usual discussion of relevant operators close to the unitarity bound in AdS/CFT [22], when near-boundary fall-offs are slow enough that one can switch to general multi-trace boundary conditions [23]. Secondly, the linearized metric perturbations which contain both short and long-lived modes can be decoupled by working with diffeomorphism/gauge invariant degrees of freedom. This provides a clean separation of fast and slow modes, allowing one to integrate out the former. It also makes clear what the natural gauge choices are for analyzing perturbations in the real-time gravitational Schwinger-Keldysh (grSK) geometry of [8]. For instance, use of standard gauge fixing in the AdS black hole background leads to solutions which have spurious singularities at the horizon, necessitating artificial boundary conditions in the interior of the spacetime, in tension with the rules of the AdS/CFT correspondence (this is an issue for instance in the solutions discussed in [8] and [10]). In contrast, the gauge invariant variables chosen in [6,7] leads to a smooth solution on the grSK geometry with appropriate boundary conditions on the two asymptotic boundaries (corresponding to the bra and ket pieces of the Schwinger-Keldysh evolution).
Let us understand these statements by examining the dynamics of the plasma stress tensor. Given a direction of propagation we can identify polarizations labeled by the little group in the transverse spatial geometry. The stress tensor operator has traceless spin-2 polarizations which are short-lived, and hence are Markovian with index M = d − 1. The transverse vector spin-1 polarizations are non-Markovian with index M = 1−d and correspond to the shear modes of the plasma. This leaves the single longitudinal mode which is the focus of the present paper -it relates to the sound mode resulting from energy transport. These statements translate directly in the dual gravity picture, for the selfsame decomposition can be applied to the linearized gravitons [24,25].
In [6] and [7] we derived the effective dynamics from the dual gravitational perspective for the spin-2 transverse traceless tensor and spin-1 transverse vector modes. The non-Markovian shear mode was captured by a designer scalar with index 1−d. This scalar, which is weakly coupled near the AdS boundary due to the dilatonic modulation, should be quantized with Neumann boundary conditions for purposes of computing the generating function of dual stress tensor correlators. This is not an ad-hoc choice, but rather one enforced by the underlying Einstein-Hilbert dynamics which, when distilled through the field redefinitions necessary to arrive at the decoupled designer scalar dynamics, ensures that the variational principle arises with suitable boundary terms for this choice.
Since the shear mode is long-lived, constructing a local effective field theory requires that we treat it in a Wilsonian fashion. Rather than compute the generating function of Schwinger-Keldysh correlators one therefore chooses to compute a Wilsonian influence functional (WIF) parameterized by the long-lived boundary modulus field -the momentum flux operator for shear modes. In practice, one can obtain this WIF by Legendre transforming the generating function of correlators. Happily, the Legendre transformation is simple and can be obtained by quantizing the designer non-Markovian field with (renormalized) Dirichlet boundary conditions.
Based on the description of sound propagation in relativistic plasma, one expects that the longitudinal modes of the gravitons are similarly captured by an effective designer non-Markovian scalar, which is dual to the boundary energy flux operator. While this observation is morally correct, the technical details are significantly more involved. For one, the field dual to the energy flux operator by itself (which we denote as Θ below) does not have a nice dual geometric description. Rather, a linear combination of this field and the conformal mode of gravitational perturbations has simple autonomous non-Markovian dynamics. Even more curiously, this designer field, denoted Z, has its coupling to the background geometry modulated non-trivially as a function of momentum. While the technical reason for these statements can be traced through the derivation, we do not understand the physical reason for why this should happen. On the contrary, the field Θ, while exhibiting no such pathologies, does not have simple autonomous dynamics. All of this is in stark contrast with the physics of diffusive modes; even for charged plasmas where the shear mode mixes with the Markovian charge current mode, the dynamical system allows for clear decoupling and relatively simple dynamics for the resulting non-Markovian degrees of freedom [7].
We will primarily focus on the spatially inhomogeneous modes of the plasma. Physically, we imagine cutting-off spatial momenta as k ≥ k IR and examining the dynamics of propagating sound modes above this cut-off. This is sensible for the plasma on R d−1,1 to mitigate IR effects. 2 This perspective will be important for us, since the dynamical system we analyze exhibits a somewhat discontinuous behaviour as a function of momentum. Spatially homogeneous modes (zero spatial momentum) have a qualitatively different behaviour. For one, their dynamics appears to be Markovian, and further there are zero modes that do not merge into Z. One such is the mode which corresponds to a homogeneous static heating of the plasma, which as explained in footnote 1 is unphysical (it changes the background solution). Since these modes do not directly affect the dynamics of sound, we will mostly not discuss them in the main text. Nevertheless, for the sake of completeness, we include an analysis of the homogeneous solution space. We demonstrate that it can be understood as the space of large diffeomorphisms, and moreover demarcate the part which is captured by our designer field in Appendix E (see [26] for an earlier analysis).
Part of the complication in the sound mode sector has to do with the fact that there are many degrees of freedom in the dual gravitational description. One has seven metric functions which can be combined into diffeomorphism invariant combinations. Four of these can be eliminated a-priori, three by gauge fixing and one by using an algebraic constraint arising from the dynamical equations of motion (latter for k = 0). The remaining three fields satisfy three linearly independent equations of motion, of which two can be solved by introducing the field Θ and enabling elimination of two functions. The final step is to show that this field combines with the conformal mode to produce the designer field Z for spatially 2 The Goldstone mode for sound in Minkowski spacetime has soft modes which may be tamed by considering the plasma on a large compact spatial sphere, i.e., using global AdS to provide a regulator. While this is an interesting problem to study, working with a momentum cut-off will suffice to extract the physics of fluctuating phonons.
inhomogeneous modes. Somewhat amazingly, the off-shell Einstein-Hilbert action simplifies considerably when we parameterize the linearized gravitons using Z.
The complexities are all pushed to pure boundary terms (a consequence of the large degree of redundancy in the classical system). These also simplify significantly when we consider asymptotically locally AdS boundary conditions, allowing one to show that for purposes of computing correlation functions, Z should be quantized with Neumann boundary conditions. Earlier analyses of scalar sector quasinormal modes in AdS black hole geometries have not carefully analyzed the variational principle, leading to some inaccurate statements in the literature.
The reduction to the single scalar field Z was first ascertained by [24] at the level of equations of motion. Using their results and examining the asymptotics in global Schwarzschild-AdS 4 , [27] argued that one should impose a Robin boundary condition for the global analog of our field Z. During their study of planar black hole quasinormal modes, [28] argued that the field Z should have Robin boundary conditions in d = 3, 4 but should have Dirichlet boundary conditions in d ≥ 5. 3 Since we follow the field redefinitions and examine the variational principle, we have a clear prescription to obtain an unambiguous answer: Z should be quantized with Neumann boundary conditions to obtain boundary stress tensor correlation functions. 4 One way to see the issue is to note that the field Z, and not its conjugate momentum, gets renormalized by higher order counterterms. The renormalization of Z starts with the counterterm at quartic order in boundary gradients (the leading cosmological constant and boundary Einstein-Hilbert counterterms do not renormalize Z). Usual AdS/CFT dictionary with Dirichlet boundary conditions renormalizes the field momentum and not the field, which is held frozen as the source.
Once we have isolated the dynamics in terms of a single designer scalar degree of freedom, the rest of the analysis follows along the lines described in [6]. We first obtain the ingoing boundary to bulk (inverse) Green's function by solving the dynamical equation of motion, order by order in boundary gradient expansion. We discover in this process another surprise: the solution for Z can be written in terms of the solution to the minimally coupled scalar wave equation (equivalently the tensor mode solution) with some simple replacement rules up to cubic order in gradients. At quartic order there are new functions which reflect the change in the nature of the dynamics.
Armed with the boundary to bulk (inverse) ingoing Green's function we can construct the full boundary to bulk solution on the grSK geometry, parameterized by the expectation value of the corresponding boundary plasma operator. This information then suffices to obtain the 3 A different combination of the linearized metric components with decoupled dynamics was constructed in [25] in radial gauge. This combination should be quantized with Dirichlet boundary conditions as noted by the authors. We have not made direct contact with their variables, but note that the relation to the field Z can be recovered from Eq. (4.7) of [28]. 4 This is true even in low dimensions d = 3, 4 where the index for Z lies in the window M Z ∈ (−1, 1) mentioned earlier, allowing for both sets of boundary conditions. So while general boundary conditions are technically allowed, they don't compute stress tensor correlation functions.
WIF at the Gaussian order by evaluating the on-shell action (which gives the saddle point semiclassical answer). This expression can be equivalently written in terms of the energy flux operator by a suitable field redefinition. This is the primary result of the paper: we reproduce in the process the expected correlation function of the energy-momentum tensor isolating the locus of the sound pole and obtaining thus the sound dispersion relation of the holographic plasma. Our results are consistent with the earlier computations of [29] and [30], who obtained the non-linear sound dispersion for N = 4 SYM (d = 4) and ABJM plasma (d = 3), respectively, and generalize them to arbitrary dimensions.
The outline of the paper is as follows. We begin with a quick summary of the grSK geometry which forms the arena of the computation in §2. In §3 we argue that the dynamics of Einstein's equations can be distilled into one designer field. We primarily present the final result of the analysis and give a bird's-eye view of the arguments leading up to capturing the dynamics into a single field. We describe in §4 the solutions of the designer field and its on-shell action on the grSK geometry. We then use this data in §5 to write down the open effective action for the energy flux operator, completing thereby the task initiated in [6] for the neutral plasma.
To avoid cluttering the text with technical details, we have tried to summarize the essential points as much as possible. Readers interested in understanding the steps leading to our statements in detail are invited to consult the appendices. Appendix A gives a detailed argument of how to assemble the gauge invariant data and use the dynamics implied by Einstein's equations to deduce the field Z for spatially inhomogeneous modes. While the final result is not original, with a close relative of Z having already been motivated in [24], we have tried to make more transparent the origins of the designer field. On the other hand, Appendix B contains a careful examination of the Einstein-Hilbert dynamics at the level of the off-shell action, which has hitherto not appeared in the literature. For completeness, we give two presentations, one in terms of metric functions after suitable gauge fixing to make connection with the dynamical equations, and another in terms of the designer field Z. The latter is important for our analysis, since we wish to establish the boundary conditions for Z -we prove that it satisfies Neumann boundary conditions asymptotically for computing the generating function of correlations. In Appendix C we give the expressions for the boundary observables which we use in § 5. Finally, Appendix D compiles details of the solution we obtain in gradient expansion. In particular, we give the near-boundary asymptotics of the functions which enter the computation of the boundary observables.
For completeness, we also include Appendix E where we characterize the dynamics of the spatially homogeneous modes and relate them to the large diffeomorphisms of the background geometry. We demonstrate that there are additional zero modes in the problem, which do not smoothly connect to the solution space parameterized by Z.

Background
We are interested in analyzing the real-time dynamics of energy current and the associated propagation of sound in a neutral plasma. Holographically, a neutral conformal plasma is dual to the Schwarzschild-AdS d+1 geometry, which is a solution of the Einstein-Hilbert action (with a negative cosmological constant). In ingoing Eddington-Finkelstein coordinates the background metric reads For analyzing real-time dynamics in the d dimensional field theory we want to use the grSK geometry. In the conventions introduced in [5,6] the metric reads Here ζ is the mock tortoise coordinate, a real-time contour in the bulk complexified radial coordinate Fig. 1 and is defined by subject to the following boundary conditions at the cut-off surface r = r c For further details we refer the reader to the aforementioned papers. 5 Re(r) Figure 1: The complex r plane with the locations of the two boundaries and the horizon marked. The grSK contour is a codimension-1 surface in this plane (drawn at fixed v). As indicated the direction of the contour is counter-clockwise and it encircles the branch point at the horizon.
When we analyze the dynamics of gravitational fluctuations, it will be convenient to work with suitably dressed radial derivative operators that make Schwinger-Keldysh time reversal properties manifest. Time-reversal is a Z 2 involution acting as v → iβζ − v and thus a natural 5 As in those references, uppercase Latin indices are used for the bulk AdS spacetime, with lowercase Greek indices reserved for the timelike boundary. Spatial directions along the boundary are further indexed by lowercase mid-alphabet Latin indices. Furthermore, gAB is the bulk metric, γµν the induced boundary metric, and n A the unit normal to the boundary. basis tangent and cotangent spaces are {D + , ∂ v , ∂ i } and { dr r 2 f , dv − dr r 2 f , dx i }, respectively. We defined following [6] the operators D ± (2.5) in the time and frequency domain, respectively. In presenting our results we find it helpful to scale out dimensions using the horizon scale set by b −1 and work with dimensionless parameters

Dynamics and the designer sound field
As discussed in §1, to understand energy transport and the sound modes in the plasma, it suffices to focus on metric perturbations involving scalar plane waves and their derivatives. These scalar polarized gravitons will be the only set of modes we analyze here; the discussion of tensor and vector polarizations was previously described in [6,7] for neutral and charged plasmas, respectively. For this decomposition we pick a direction for the spatial momentum k and define harmonics on R d−1,1 to be the SO(d − 2) harmonics of the corresponding little group. Concretely, we consider metric perturbations of the form (3.1) Here S = e ik·x−iωv is the scalar plane wave on R d−1,1 and S i and S T ij are derived harmonics, defined as It is thus distinguished from the derived harmonic S ij defined in [6], which is neither transverse nor traceless. We have also defined a short-hand for the measure on the Fourier domain to keep expressions compact.
There are seven metric components in the perturbation above. A-priori we expect that these seven functions obey seven coupled radial ODEs arising from the linearized Einstein equations. Owing to diffeomorphism invariance not all of these dynamical equations are independent. One must first identify the pure gauge modes from the physical perturbations and focus on their dynamics. This problem has been analyzed in detail in [24], whose discussion we can adapt for our purposes. Following their analysis we will refer to the two-dimensional {v, r} spacetime as the orbit space. The physical sound mode ends up being described as a non-Markovian scalar field in this orbit space.
We now explain how to efficiently distill the dynamics into a single gauge invariant degree of freedom inspired by the previous analysis of gauge dynamics in [6,7]. 6 The key point is that there are seven metric functions while the classical phase space, parameterized by gauge invariant data, has only one physical degree of freedom in the scalar sector of gravitational perturbations. To arrive at this result we make use of the following observations which are elaborated upon in Appendix A: • There are seven diffeomorphism invariant combinations of the metric perturbations Ψ AB which can be organized as an orbit space traceless tensor, an orbit space vector, and two orbit space scalars. One first deduces that the orbit space vector (Ψ vx and Ψ rx ) and one scalar (Ψ T ) can be gauge fixed to vanish, leaving four functions, which are essentially Ψ S , Ψ vv , Ψ vr , Ψ rr .
• Time-reversal involution is an orbit space diffeomorphism: Ψ S , Ψ vv , Ψ vr + 1 2 r 2 f Ψ rr are time-reversal even, while Ψ vv + r 2 f Ψ vr is time-reversal odd. This information is useful to constrain the structure of the equations of motion and the action.
• Among the equations of motion we find an algebraic constraint, which allows elimination of Ψ rr for non-zero spatial momentum. Two other equations are the momentum constraint equation (which is the boundary energy-momentum tensor conservation) and a first order radial equation. These two equations can be used to solve for Ψ vv and Ψ vr in terms of a function Θ, which is related to the boundary stress tensor component (T CFT ) i v .
One therefore finds that the metric can be parameterized by two functions: an overall Weyl rescaling of the background geometry and a function which encodes the physical data of the boundary stress tensor. We will label these fields as Φ W and Θ, respectively, and judiciously define them with suitable factors of r to simplify the dynamical equations; The metric parameterized by these fields, subject to the gauge fixing where all the metric components involving derived scalar spherical harmonics are set to zero, is said to be in the Debye gauge. Our scalar perturbations are then captured by 7 We have written the metric in the basis of the cotangent space that is adapted to the timereversal involution of the background. The last step involves analyzing the remaining dynamical equations of motion and discerning that they can be solved if one further introduces a field Z to parameterize Φ W and Θ as The function Λ k is curious. It is a non-trivial function of spatial momentum (indicated by the subscript). It will turn out to be a designer dilaton for the field Z and is given by Because of the momentum dependence, the field Z should be seen as residing in the orbit space and not in the entire background geometry. 8 Note that there is a linear relation between Θ and Φ W from (3.6) Our main claim is the following: linearized Einstein equations for scalar perturbations of the Schwarzschild-AdS d+1 geometry are satisfied provided the field Z obeys a second order linear differential equation This equation is the 'master field equation' obtained in [24], which we have obtained in a somewhat different parameterization. 9 We have written the effective dynamics in a manifest time-reversal invariant form. This implies that for analyzing the solution on the grSK geometry (2.2) we only need to solve (3.9) in a single copy background (2.1) with ingoing boundary conditions. One can then use the time-reversal properties to construct the outgoing solution and thence the full grSK solution with boundary conditions specified on the two asymptotic boundaries at r → ∞ ± i0 7 We have written the metric directly in position space as the fields Φ W and Θ are simply Fourier transformed with the scalar harmonic S. 8 The origin of Λ k is analogous to the Ohmic function h(r) which appears in the analysis of vector perturbations of the Reissner-Nordström-AdS d+1 black hole [7]. The modulation there was due to the background charge whereas here it directly relates to the momentum carried by the perturbation. 9 The authors of [24] prefer to write the equations in Schwarzschild coordinates and express it as a Schrödinger equation in orbit space for a field S KI . To bring our equation to their form, one first transforms from our ingoing coordinates to Schwarzschild coordinates and implements a field redefinition: as described in [5]. Arriving here was a key step in the analysis: had one worked with other gauge fixing methods traditionally employed in AdS/CFT such as radial gauge, while one would have been able to solve the ingoing problem, the corresponding outgoing solution would be singular. The issue is similar to the problems encountered with gauge fields and momentum diffusion analyzed in [6,7]. Not only does one end up with a single second order differential equation to solve in the scalar sector of gravity, but there is also a remarkable simplification of the Einstein-Hilbert action. Plugging in the parameterization, we find after a series of algebraic simplifications the dynamics of Z to be governed by a simple dynamical system 10 (3.10) The dilaton χ s which modulates the gravitational interaction is while the normalization is fixed by a parameter ν s , which in turn is This is a remarkable simplification given the complexities inherent in the scalar perturbation; our field redefinitions in (3.6) imply that the metric in Debye gauge is a function of {Z, D + Z, D 2 + Z}, so the truncation of two derivative dynamics is indeed a welcome surprise. As one might expect, much of the complication is hidden in the boundary term S bdy [Z] in (3.10). It captures all the contributions from the Gibbons-Hawking term, additional boundary terms encountered while writing the action in terms of Z, and counterterms. In Appendix B, we explain how to obtain the action and the variational principle for the field Z starting from the Einstein-Hilbert dynamics. This action (3.10) is of the general form of a non-Markovian designer scalar introduced in [6], albeit with some additional novelties. Firstly, the dilaton χ s modulates the gravitational interaction non-trivially as a function of spatial momentum k. The Markovianity index, 11 which depends on the asymptotic behaviour of the dilaton, is determined to be lim r→∞ e χs = 1 (3.13) We have made explicit our spatial momentum cut-off for clarity. Thus, we see that for nonzero momentum the field Z is a non-Markovian designer scalar with index M = 3 − d. Note 10 The effective central charge is defined as c eff = d−1 AdS 16πG N . 11 The Markovianity index was defined in [6] with a minimally coupled massless scalar having index d − 1.
that while the field Z has a non-trivial potential which we will need to take into account, the potential does not, as in the analysis of [7], modify the Markovianity properties. The latter is purely governed by the radial kinetic operator.
While our distillation of the dynamics into the designer field Z is only valid for spatially inhomogeneous modes, for purposes of finding the solution, we can examine the behaviour of the dilaton and the wave equation (3.9) at zero spatial momentum. On this locus Markovianity index changes to M = d−1, which suggests that the zero mode sector comprises of short-lived modes (sound of course doesn't propagate without momentum). The full dynamics of k = 0 modes is a bit more involved, though a part of the solution space is indeed captured by Z with an effective Markovian dynamics as (3.9) suggests.
From a pragmatic standpoint, this switch between Markovian and non-Markovian behaviour in §4.2 will prove very useful when we solve (3.9) in a gradient expansion. Thanks to this observation, we will be able to write down the solution directly in terms of known solutions [6] of the Markovian wave equation with M = d − 1 at low orders in the gradient expansion. As mentioned earlier, an analysis of zero modes can be found in Appendix E, where we describe how Z(r, ω, 0) connects onto the set of large diffeomorphisms.

Sounding out the grSK geometry
For the rest of the discussion we will focus on the dynamics of the gravitational system encoded in Z and recover the physics of sound propagation with attenuation. We begin with the solution the equation (3.9) with ingoing boundary conditions. Subsequently, using the time-reversal involution of the grSK geometry we will construct the full linearized solution parameterized by the expectation value of a boundary stress tensor component at the L and R boundaries.

Solving the designer equation
To better understand the nature of the designer scalar Z(r, ω, k), we first analyze the asymptotics of the wave equation (3.9). As described above, while the derivation of the Z equation is valid for k = 0, we can consider k ≥ 0, since part of the homogeneous mode solutions merge into Z.
Let us first focus on zero frequency solutions, but examine both the zero momentum and non-zero momentum behaviour separately, owing to (3.13). At zero spatial momentum, we have the asymptotic behaviour determined by where the constant mode c a is analytic, but the monodromy mode c m typically has a logarithmic branch cut emanating from the horizon. On the other hand at non-zero spatial momentum one finds the expected non-Markovian behaviour: Per se, this is not unexpected given the general analysis of Markovian and non-Markovian degrees of freedom as described in [6]. This is an explicit realization of the observation above that naively there is a change in character as a function of momentum. Let us pause here to comment on two special cases: d = 3 and d = 4. The situation in d = 4 is marginal; the monodromy mode behaves logarithmically. When we present solutions we will be careful to normalize them as r d−4 d−4 so the limit d → 4 can be easily taken by replacing this function by log r. On the other hand, in d = 3, it appears from the fall-offs that the field Z is Markovian since the index vanishes. This is misleading since, as noted in §1, the field Z is actually non-Markovian for all d ≥ 3. This is not manifest from just analyzing fall-offs, for we also need to take into consideration the variational principle and the boundary conditions for Z, which we discuss in Appendix B.3. We demonstrate there that Z has Neumann boundary conditions imposed on it for purposes of computing the generating function of correlators and hence it is non-Markovian for all d ≥ 3.
Returning to the dynamical problem we can solve (3.9) by disentangling the physical non-Markovian behaviour from the auxiliary Markovian one. We first rewrite the differential equation for Z in terms of the designer scalar wave operator [6,7]: When expanded in powers of momentum we find a remarkably simple form for (3.9): where the parameter ν s was defined earlier in (3.12).
Since the operator D d−1 annihilates a Markovian scalar of index d − 1, which we write as the field ϕ d−1 in the notation of [6], one can subtract out this piece from Z and write a general solution as We are almost done: Z satisfies an inhomogeneous linear differential equation of the form D d−1 Z = source. The source is simply determined in terms of a Markovian field ϕ d−1 in the background Schwarzschild-AdS d+1 geometry, which has already been solved for in [6,7]. That analysis has already inverted the operator D d−1 , so one simply needs to add in a particular solution to determine the full behaviour.
Not only do we have access to the solution space, but the nature of the change in asymptotics is now transparent. At zero momentum the second term in (4.4) vanishes and simply has D d−1 annihilating Z(r, ω, 0). At non-zero momentum however, the source provided by the Markovian solution is scaled up because of the 2 r(1−f ) D + term. This accounts for the change in behaviour. We will parameterize the solutions in a particularly convenient manner so that the non-Markovian behaviour in fact only sets in at quartic order. 12 This will turn out to be another artifact of the field Z; translating back to the metric function Φ W , or even Θ, we will see a change at quadratic order in momenta.
We now present the ingoing solution for the field Z and use it to determine the full solution on the grSK geometry. Our analysis will be accurate to quartic order in the boundary gradient expansion.

Ingoing solution in gradient expansion
Working with the dimensionless variables (2.6), we obtain the ingoing inverse Green's function for the field Z, normalized to satisfy lim ξ→∞ G in Z (ξ, ω, 0) = 1: (4.6) We have written the solution using the notation employed in [6] since the functions appearing below are essentially those already encountered there. 13 We have introduced two new functions of frequencies and momenta in parameterizing (4.6). The first, K s (ω, k), which we will later confirm to be the sound dispersion function, is defined as (4.7) Up to quadratic order K s captures the propagation of sound while Γ s encodes its attenuation. The second parameter, p s , may be viewed as a 'deformed momentum' parameter arising from the spatial modulation of the dilaton and is In fact, this behaviour is quite similar conceptually to that observed for transverse vector polarizations of photons and gravitons in [7]. There it was found that one had a Markovian mode which mixed with a non-Markovian mode and picked up additional divergent terms. In that discussion there were two independent degrees of freedom, which were decoupled to isolate the long-lived and short-lived modes. Here there is only a single mode whose character changes owing to the underlying gauge invariance. 13 In [6] we had only obtained the solution to cubic order in gradients. Solutions accurate to the quartic order were obtained in [7] both for the neutral and charged black holes. In the latter reference the gradient expansion used a slightly different grouping of terms in the solution; we describe this form for completeness in Appendix D.3.
We have judiciously combined terms from the solution for ϕ d−1 and Z to write the result in this compact form (introducing the parameter p s in the process). 14 There is one important subtlety which the reader should be aware of in the way we have presented the solution for Z. The solution to (3.9) up to the quartic order in gradients does not determine Γ s owing to the explicit factor of q 2 multiplying V k (ξ) in (4.6). This would be unfortunate since Γ s will turn out to be the sound attenuation function. However, if we compute the functions Φ W and Θ from Z (or any other metric function) then this overall factor of q 2 disappears and we obtain Γ s accurate to quadratic order by examining the coefficient of the non-normalizable mode in the solution. A physical way to say this is to note that the boundary source is the conjugate momentum with a factor of k 2 stripped off, see (5.5). This behaviour is manifest in the ingoing Green's function of the field Φ W which is given in (D.10) (the information in Θ is similar, but not independent owing to the linear relation (3.8)). To obtain Γ s directly from Z as written above, we would have to compute the solution accurately to sextic order, another peculiarity of the fact that it is related to radial derivatives of the metric functions.
In fact, we would like to conjecture that one can actually get Γ s to quartic order in gradients. From the leading divergent mode in Φ W (or Θ), which scales as r d−2 , we get an expression for K s to quartic order. However, we can also examine the constant mode, scaling as r 0 , near the boundary. It has a coefficient which is a non-trivial function of ω and k. By judiciously parameterizing the solution for Φ W , one finds that this function starts off as Γ s and gets corrected at cubic and quartic corrections (the same quantity appears as the constant mode of G in Z ). If one uses a parameterization of the solution in terms of Γ s , simply noting that it starts at linear order, then the on-shell action is finite up to sextic order. Based on this observation and the nature of the explicit solution, we predict an expression for Γ s up to the fourth order in boundary gradients, which we record for completeness in (D.9) (where for clarity it is denoted as Γ s (ω, k)).
All of the functions that appear above, except for V k (ξ), are defined using the ingoing solution for massless Klein-Gordon scalar ϕ d−1 in the Schwarzschild-AdS d+1 background. More precisely, solves D d−1 ϕ d−1 = 0 with ingoing boundary conditions, lim r→∞ ϕ d−1 = 1, and thus is the ingoing bulk to boundary Green function for the massless scalar described in [6]. The asymptotic boundary condition along with the regularity at the horizon (ingoing boundary condition) 14 One way to observe that p 2 which has the form of a Markovian wave operator with analytically continued momenta.

Function
Source Function Source uniquely fixes the eight functions {F, The asymptotic boundary condition implies that they all vanish as r → ∞. The conversion from the massless scalar to the field Z involves replacing factors q 2 with p 2 s in various places up to the quartic order: this in particular means that G in Z remains Markovian up to that order (for reasons explained in footnote 14).
The key function that controls non-Markovianity is V k (ξ), which is defined by the integral For d > 4 this function grows as r → ∞, unlike the other functions. We have extracted the leading divergence as r d−4 , which is precisely what one expects for a non-Markovian field of index M = 3 − d as described around (4.2). Since we are effectively inverting the Markovian operator D d−1 to find the functions above, they can all be given formal integral representations, assuming that the sources are regular on the horizon and do not grow too fast at infinity. We can write in general following [6] with F(ξ) = F(ξ) − F(1) defined to measure function values relative to that on the horizon. We tabulate the data for the functions appearing in the gradient expansion originating from the minimally coupled scalar in Table 1. Further details of the solution, including asymptotic expansions and expressions for the gravitational data from our solution for Z, are compiled in Appendix D.

The grSK solution for the designer field
The ingoing Green's function for Z suffices for us to determine the full solution on the grSK geometry thanks to the time-reversal invariance of (3.9). We want to impose suitable boundary conditions at the two boundaries of the grSK geometry at r → ∞ ± i0. It was argued in [6] that a non-Markovian fields should be quantized with Neumann boundary conditions if we wish to compute their correlation functions. Equivalently, the asymptotic field value does not correspond to the boundary source, but rather gives the dual boundary operator (akin to the alternate quantization of low lying operators). This was cleanly formulated for probe non-Markovian fields and established to be the case for diffusive modes in the aforementioned reference and [7]. The question we should first ask is what the Einstein-Hilbert dynamics with its usual Dirichlet boundary conditions for computing correlation functions of dual energy-momentum tensor correlation functions implies for the field Z. It turns out that while the dynamics of Z is pretty simple in the bulk, as evidenced from (3.10), it has a pretty involved set of boundary terms owing to the redefinitions in (A.27). The boundary terms are a general quadratic form in Z, D + Z, and D 2 + Z, implying that generically we need to fix a combination of these three quantities for stationarity of the action. However, one can do better: armed with the asymptotic fall-offs of the field, one learns that the leading set of boundary terms in a near-boundary expansion are simpler, and one indeed finds that Z ought to be quantized with Neumann boundary conditions to obtain the dual stress tensor correlators. The detailed argument analyzing the variational principle for the designer field Z is given in Appendix B.3.
While the field Z is quantized with Neumann boundary conditions in order to compute the generating function of boundary correlators, the resulting correlators have a sound pole, reflecting the non-Markovian nature of the field. For this reason, it was proposed in [6] that one should, for purposes of computing a real-time Wilsonian influence functional, parameterize the solution for Z in terms of the normalizable mode, which for a non-Markovian field corresponds to the conjugate momentum.
With this understanding we will now parameterize the general solution on the grSK geometry in terms of the sound modulus, which is the expectation value of the dual boundary operatorȎ Z (thinking of Z probe field in the fixed Schwarzschild-AdS d+1 background). We will denote the modulus associated with the auxiliary field Z asZ and write In terms this boundary modulus fieldZ we can write the full grSK solution in the average-difference basis as being the time-reversed propagator. It is useful to note that Z has mass dimension d − 2 and henceZ a,d are likewise boundary fields with this dimension. We will use this information in the next section to present the effective action for the phonon modes in the relativistic plasma.

Effective dynamics of sound and energy transport
With the grSK solution for the designer sound mode, we can evaluate the on-shell action parameterized by the fieldsZ a,d , which then gives us the open effective field theory of sound propagation in the holographic plasma at the Gaussian order in amplitudes. We first outline the details of this effective action and then turn to the boundary stress tensor, which we will express in terms ofZ and background polarization terms involving the sourcesζ.

The sound Wilsonian influence functional
Since the background Schwarzschild-AdS d+1 geometry has a non-vanishing free energy, we expect to see two contributions to the on-shell action. One is a piece from the background, which based on fluid/gravity intuition should correspond to the ideal fluid free energy. In addition, there will be the true dynamical data corresponding to the Wilsonian influence functional (WIF) of the field Z. These two contributions are cleanly separated in the shear and tensor sectors of a neutral fluid [6] because, in those sectors, we do not have a propagating mode.
A relativistic ideal fluid has a propagating phonon mode, so the split in this sound sector is not a-priori manifest. We will therefore not identify the ideal fluid contribution at this stage, but simply separate the on-shell action into the WIF and contact terms. Subsequently, working out the stress tensor will enable us to understand which pieces should be regarded as part of the ideal fluid contribution. With this preamble, let us write the full action as a sum of two pieces, viz., We now summarize the final result for the Schwinger-Keldysh effective action obtained by computing the on-shell action for Z with the gradient expansion solution given in §4.2. The reader can find details of the evaluation in Appendix C.3. Let us start with the contribution to the WIF, which is ascertained to be From this expression we solve for the boundary source of the field Z in terms of the moduli fieldZ and obtainζ One can give an alternate expression for the source directly from the conjugate momentum of the field Z, after stripping off a factor of k 2 , viz., It can be checked that our identification of the source in (5.5) agrees with the expectation from the WIF (5.4) (which we expect on general grounds from [6]) and is verified in Appendix B. Isolating a factor of k 2 in the WIF results in stripping off a similar factor from the conjugate momentum (5.5). 15 The rescaling by a factor of k 2 is only allowed since we are focusing on spatially inhomogeneous modes and are implicitly working with k > k IR . The source for the designer field can be given several equivalent expressions in terms of the metric functions. Using the asymptotics of the solutions obtained in Appendix D.2 (see in particular (D.10) and (D.12)) one can show These source terms are basically capturing the deformation of the boundary metric order by order in the gradient expansion. Indeed, upon examining the induced metric γ µν on the boundary we see that spatial and temporal components of the perturbed boundary metric can be viewed as sources for Z. Since we have only one physical degree of freedom in the longitudinal sector, we do not have independent metric perturbations, but rather see that the red-shift factor captured by the temporal term is related up to a dimension dependent constant to Weyl rescaling of the background. With the identification of the sources we can now present the contact term part of the action, which is a functional of these sources. Ignoring the background free energy term for simplicity, we find the following result at linear and quadratic order in amplitudes: In §5.4 we will argue that this contact term can be understood as arising from the on-shell action of an ideal fluid propagating on (5.7). This includes the somewhat counter-intuitive numerical factor in the quadratic term, which vanishes in d = 6.
We will explain the individual contribution to the action, in particular the split between ideal and non-ideal parts, after we discuss the stress tensor. With that understanding we will be able to cleanly identify the ideal fluid contribution. Along the way, we will also argue that the non-dissipative part of the hydrodynamic action, the Class L terms in the terminology of [31], can be extracted from the WIF. For the scalar sector we will for example see the curvature coupling of the fluid at quadratic order in gradients. These statements will be elaborated in §5.4.
Let us take stock of some physical implications from the grSK solution for Z, and the results for the Wilsonian influence functional. We see that the inverse Green's function computed from the WIF is proportional to K s (ω, k). In other words, the retarded Green's function ofȎ Z has a pole at the vanishing locus of K s . This function also appears as the coefficient of the divergent (non-normalizable) mode in Z. As described in [6], non-Markovian fields have a completely normalizable solution on a codimension-1 locus in the boundary Fourier domain. For Z this is the vanishing locus of K s ; it will end up defining the dispersion function for sound. This also implies that Γ s defined in (4.7) is the rate of attenuation of sound speed due to viscosity. We will elaborate on this further below when we discuss the physical sound degree of freedom and compute boundary energy-momentum tensor correlators.

The boundary stress tensor
With the Wilsonian influence functional parameterized by the boundary value of the field Z at hand, we can now turn to the physical boundary observables, which are the scalar polarizations of the boundary energy-momentum tensor density. The stress tensor has both a background contact piece and a contribution given in terms of the fluctuating fieldȎ Z . In the scalar sector, the presence of a non-trivial boundary metric (5.7) means that the result we quote for the contact terms depends on the index positions and whether or not we work with stress tensor densities. We will work with tensor densities, quoting the mixed components for the stress tensor operator, but compute the correlator for the operator with both indices raised. 16 For the sake of notational simplicity, we define Here T µν CFT is the counterterm corrected Brown-York stress tensor density given in (C.4). We will take the viewpoint that the standard rules of the extrapolate dictionary in holography (which can be extended to grSK geometry [5]) as applied to non-Markovian operators dictates a canonical split between operator (or vev) and source contributions. As explained in Appendix C.2, this can be seen by examining the Brown-York stress tensor (corrected by the counterterms) directly in terms of the metric fields Φ W and Θ. This tensor must satisfy conservation and be traceless (up to the conformal anomaly, which we do not access in 16 In the AdS/CFT context, it is natural, as from any effective action, to extract the boundary stress tensor density, since it only requires variation with respect to the boundary metric and no removal of metric determinants (recall that the stress tensor operator is . This statement follows naturally from the momentum constraint equation in the bulk geometry and underlies the original identification in [32].
With this choice the Schwinger-Keldysh stress tensor operator has the following representation: We recognize the background contribution which says that the unperturbed planar Schwarzschild-AdS d+1 black hole is a conformal plasma. The terms linear inζ andȎ Z are the terms we should understand.
As it is written, the stress tensor is not manifestly traceless. Nor is the conservation Ward identify obvious on the induced boundary geometry (5.7). The two do hold, and are, in fact, equivalent to the relation between sources and vevs (5.4), which picks out the sound dispersion locus: 17 We should view (5.11) as giving us the on-shell condition for the sound mode which occurs when the solution is purely normalizable, i.e., when the source contribution is set to vanish. This confirms that the operator K s gives us the dispersion relation for sound in the plasma and identifies Γ s (ω, k) as the attenuation function.
With this identification we can now give an alternate presentation of the stress tensor. Let us use the dispersion relation (5.11) to shift the source and vev contributions in the spatial part of the stress tensor, i.e., use the replacement rule Re-expressing K sZ in terms ofζ we find that the stress tensor in (5.10) can be equivalently 17 The covariant derivative in the conservation equation is the one appropriate for the stress tensor density, cf., footnote 16.
presented as The representation of the CFT stress tensor (5.13) renders the trace Ward identity manifest. Additionally, it also isolates the sound attenuation contribution captured by Γ s solely into the spatial part of the stress tensor. In fact, the contribution is governed by the longitudinal trace-free tensor structure S T ij . In this presentation, the tracelessness Ward identity is manifest, but conservation now implies the on-shell condition (5.11).
One could go a step further and replace the Γ s term in the spatial part of the stress tensor once again using (5.12), i.e., express it as In this manner of presentation, both the conservation and tracelessness Ward identities of the stress tensor are identically satisfied. The operator contribution encoded inZ has the right momentum and frequency dependent factors for the conservation to be rendered trivial, while the terms involvingζ, one can check, respect ∇ γ µ T µν = 0 by themselves. While there appear to be three distinct parameterizations, the form given in (5.13) is the one that separates the ideal fluid contribution from the dissipative part. We will demonstrate this in §5.4 after writing down the stress tensor correlators.

Correlation functions
From the Wilsonian influence functional (5.2) we can read off the correlation functions of thȇ O Z field operator on the boundary. The retarded Green's function is given by .
The structure of the WIF respects the KMS condition, as explained in the earlier works, implying that the Keldysh correlator satisfies the fluctuation dissipation condition, viz., We can use this information to write down the stress tensor correlators given the explicit expressions for the components in (5.10). From this expression it is clear that the result is given by the two-point functions ofȎ Z , suitably dressed to account for the derivative operators present (the functions of ω, k in frequency/momentum domain).
Before doing so however, we should note that we have the boundary metric determined by a single source function (5.7), which implies relation between sources of various components (and an absence of source for the spatial-temporal component). We will view the boundary metric as defining the source for the energy density. Thus, we can obtain the energy density two-point function naturally. Once we have this piece of data, we will use flat spacetime Ward identities, as explained in [4], to fix the remaining correlation functions. In particular, we demand that the contact terms in the retarded correlation functions are fixed so that flat spacetime momentum conservation holds -this implies that correlation functions with at least one temporal index vanish at zero momentum. To fix the purely spatial components we utilize full energy-momentum tensor conservation.
With this understanding, we now quote the result for the physical retarded correlator for the energy-momentum tensor density operator T µν CFT . We parameterize these as To avoid writing involved expressions we pick a spatial direction for sound propagation, setting We choose a representation where G µν,ρσ (ω, k) are polynomials in ω, k and take the form: Any other representation of G µν,ρσ would differ from the above by factors of K s . The Keldysh propagator follows naturally from the fluctuation dissipation theorem. Our choice for nonanalytic part of the retarded Green's function recovers the result for the analytically continued Wightman function given in [4] for d = 4. The physical aspect of the result which is interesting is the fact that the stress tensor correlators have a sound pole at the dispersion locus characterized by the vanishing of K s . As noted around (5.11) this function gives us the on-shell condition for the phonon mode. At the leading order in gradients it enforces the expected equation of state condition, which fixes the speed of sound in a conformal plasma to be 1 √ d−1 . At higher orders the Γ s pieces serve to attenuate the propagation and predicts its half-life.
Solving (5.11), and using the expression for Γ s given in (4.7), we find 18 We have indicated only one branch of the solution in the above expression.
As particular cases, note that with d = 4 we recover the results for the N = 4 SYM plasmas, 19 20) which was first obtained to cubic order in [29] (extending the original result of [4]). Setting d = 3 we obtain the corresponding result for the ABJM plasma obtained in [28,30]: As we describe in Appendix D.2, from our solution we can actually extract higher order corrections to the dispersion, and quote the result accurate to quintic order in (D.13).

A fluid dynamical perspective
We have now reproduced the physical results expected for sound propagation as encoded in the stress tensor correlators. Let us therefore try to analyze features of the solution from a hydrodynamic perspective and, in particular, attempt to understand the contributions to the action (5.8). The gravitational calculation gives us the answer tout ensemble, but we can attempt to decompose it using hydrodynamic intuition. The physical question of interest is how one delineates the ideal and non-ideal parts of the effective action. Addressing this question will make clear that one should think of the action as comprising of a Schwinger-Keldysh factorized part corresponding to sound propagation, in addition to the physical influence phase, the part that governs sound attenuation, as we indicate in (5.31) below.
The non-ideal part, by definition, includes all the gradient contributions in the stress tensor parameterized in hydrodynamic variables, whether or not they lead to dissipation. We will use the stress tensor in the form parameterized in (5.13), which judiciously isolates the non-ideal contributions into Γ s . As has been argued earlier [31], not all higher order transport 18 We where Har (x) is the Harmonic number function. 19 We work with dimensionless frequencies and momenta defined in (2.6). This definition differs from normalizations used in earlier references. Our normalization is twice that used in [29], while [30] uses a normalization set by temperature and not (inverse) horizon size, which in d = 3 differs by a factor of 3 4π .
is dissipative. While dissipative transport leads to entropy production, in general, there exists non-dissipative transport which is adiabatic and leads to no entropy production. While at leading order Γ s (ω, k) captures sound attenuation, which originates from the dissipative shear viscosity term, it also includes contributions from higher order non-dissipative gradient terms.
To keep the discussion transparent, we will first identify the ideal fluid contribution, which by definition only captures zeroth order terms in the gradient expansion of the stress tensor. Having understood this part, we will then attempt to address higher order (spatial) gradient terms, which capture non-dissipative transport.
Consider an ideal fluid with energy density and pressure related by the conformal equation of state = (d − 1) p, viz., the tensor density (nb: T is dimensionless for simplicity) 20 We claim that the part of our solution which is insensitive to sound attenuation, i.e., with Γ s → 0, describes the dynamics of such an ideal fluid on the boundary geometry (on either L or R boundary). Specifically, we assert We can confirm this by solving the conservation equations arising from (5.22) on the induced boundary geometry (5.7). We can parameterize the temperature and velocity by a fieldZ, which is a-priori unrelated to the stress tensor expectation value Imposing the conservation equation one finds thatZ andζ must satisfy the relation (5.11) with Γ s = 0. 21 With this constraint recovered, we may identifyZ as the physical phonon mode, i.e., as the boundary value of the non-Markovian field. We have effectively isolated 20 We are not keeping track of the normalization factor translating between the horizon size parameter b and the physical temperature T , cf., (2.3). 21 This parameterization can be motivated by considering a phonon mode for a relativistic plasma in flat spacetime. All one needs is the statement that the dynamics is captured by conservation of the stress energy tensor. It is not hard to check that for a linearized perturbation about an equilibrium plasma in flat spacetime, the dynamical sound mode, which importantly does exist even in the absence of dissipation, on the inhomogeneous dynamical boundary spacetime (5.7). One can, furthermore, use the on-shell relation to write the temperature more suggestively as In this presentation, we see that the sourceζ contribution to the temperature is just the red-shift effect for a fluid propagating on (5.7). The contribution fromZ parameterizes the response of the fluid. Beyond the ideal fluid, the first correction comes from the dissipative shear viscosity term, which physically leads to the damping of the sound in the medium. This is the leading iω contribution to Γ s . To isolate any non-dissipative contributions we must therefore switch off time-dependence and focus on equilibrium data. In stationary equilibrium, Γ s (0, k) = d−3 ij part of the spatial stress tensor can be identified as the coupling of the fluid to background curvature. For a conformal plasma it takes the form Here C γ µανβ is the Weyl tensor of the boundary geometry (5.7), ω vor is the fluid vorticity, σ sh is the shear tensor of the fluid, and the angle brackets indicate transverse projection. 22 The value of the transport coefficient κ = 2 c eff T d−2 is known for Schwarzschild-AdS d+1 black holes [33] (it was initially derived in d = 4 in [29]). Using the temperature and velocity profiles identified above, one can directly check that our result captures this contribution to the stress tensor.
Having understood the contributions at the level of the stress tensor we can now explain how to interpret the contact terms in the action. Furthermore, we can also isolate terms corresponding to the Class L adiabatic action for holographic plasmas conjectured in [31], by reverse engineering a boundary action from transport data. In terms of the Weyl covariant hydrodynamic variables this action reads Here R W is the Weyl covariant curvature scalar on the boundary. Let us begin with the ideal fluid part which is the leading contribution in (5.27). To understand this we first note that the contribution to S[Z] can be understood directly from the variational definition of the stress energy tensor. Specifically, the on-shell action with Γ s → 0, which prior to our Legendre transformation is the usual generating function of correlators, is given by contracting the ideal stress tensor (5.22) with the change in the background metric from flat spacetime, viz.,

(5.28)
We have dropped the background constant free energy part and focused on the pieces arising from the solution to the linearized equations of motion. In particular, the contact term contribution in its entirety originates from the propagation of an ideal fluid on (5.7). It should now be clear that the curious factor of (d − 6) is just a numerical accident; it arises due to the relation between the metric components in γ µν . There is nothing special about relativistic conformal fluids in six spacetime dimensions, nor are Schwarzschild-AdS 7 black holes (and the dual (0, 2) SCFT plasma) in any way singled out.
Having understood the connection between the gravitational on-shell action and the stress tensor, we can connect to the adiabatic effective action (5.27). Prior to Legendre transformation the ideal part is simply the free energy evaluated on the sound mode solution. On the grsK solution we can represent it using the rescaled thermal vector In order to compute the quadratic part of the action it will suffice to know rescaled thermal vector accurate to linear order in amplitude, consistent with our identification using the stress tensor. The reason is that the amplitude expansion of the ideal fluid action results in two terms. One is the contribution which leads to the ideal fluid stress tensor in (5.28), and the other the 'heat current' term, which originates from the change of the rescaled thermal vector (the δb µ = b µ − 1 b ∂ v variation), cf., [31]. Using the on-shell relation for the ideal fluid one can check that these two contributions nicely sum up (up to the aforementioned Legendre transform) into the terms arising from (5.30). This is structurally similar to the Wilsonian influence functional in the vector sector, which captures the shear modes driving momentum diffusion in [6]. The main difference in that case was that since there was no propagating mode; the ideal piece was purely expressible in terms of a contact term, and moreover could be isolated directly from the boundary terms of the Einstein-Hilbert action.
For the sound mode we can re-express the on-shell action as an ideal piece and a term that captures sound attenuation. To wit, Here the 'LT' term in the subscript is present to remind us that one should Legendre transform the ideal fluid part to account for the fact that we have a sound pole. From the Γ s term we can also reproduce the Weyl curvature contribution in (5.27) above, as promised. One finds We have used the on-shell relation (5.4) betweenZ andζ to make clear that this term arises from the spatial curvature of the boundary geometry.

Discussion
In this paper, we have extended the analysis of open quantum systems with memory to include Goldstone modes with a decay width. The earlier analysis of [6,7] focused on diffusive modes in thermal plasmas, which does not incorporate such propagating modes. Specifically, we analyzed the dynamics of energy transport and the physics of associated sound modes in a relativistic thermal plasma. While there is a physical difference in the nature of the longlived modes, our analysis confirms the general paradigm articulated in these earlier works. Namely, the long-lived modes can be captured into a Wilsonian Schwinger-Keldysh effective action, which we derived for conformal relativistic plasma with a holographic dual. The key result is the Wilsonian influence functional parameterized directly in terms of the boundary expectation value of the energy flux operator (T CFT ) i v . While the field theory result shares many characteristics with the corresponding effective action for diffusive modes, there are interesting technical peculiarities in the gravitational description. For the diffusive modes one was able to distill the bulk dynamics into non-minimally coupled designer scalar fields (one per polarization), where the non-minimal coupling was captured by an auxiliary dilaton, whose primary characterization was its asymptotic fall-off rate (the Markovianity index). This auxiliary dilaton however modulated only the interactions in the radial direction, i.e., as a function of energy scale in the field theory, but was spatially homogeneous. This no longer holds for the bulk dual of the sound mode; the auxiliary dilaton has a non-trivial modulation along the spatial directions of the boundary. It nevertheless remains true that the dual field has a non-Markovianity index 3−d for spatially inhomogeneous modes.
Furthermore, the gravitational analysis gives a beautiful picture for the dynamics of energy transport. The physical phonon degree of freedom is part of the ideal fluid, and thus should be captured by the hydrodynamic sigma model (Class L) actions of [31]. Owing to the presence of a gapless mode, one should not construct the sigma model action directly, but rather write the Wilsonian analog, which effectively captures the two derivative kinetic term of the Goldstone mode. Since this part is conservative, the resulting Schwinger-Keldysh effective action is factorized into L and R pieces, as indicated in (5.31). The dynamical information, viz., the dispersion relation (5.4), is obtained from this action by the constrained variational principle outlined in the aforementioned reference. The Class L action also captures higher order non-dissipative contributions, like the background curvature coupling (5.32). 23 Once we have separated out the propagating mode, what is left is the physics of sound attenuation. Since this is driven by the leading order dissipative terms, the shear viscosity of a conformal plasma, the structure is isomorphic to that found for momentum diffusion in [6]. In other words, the physical influence phase of the sound mode is the non-factorized part of the Wilsonian influence functional, with a physical kernel Γ s (ω, k). The dissipative part of this kernel is not captured by the Class L sigma model actions, as it should be; it is these frictional effects which drive the plasma to behave as an open quantum system. So in a sense, Γ s (ω, k) is the physical influence functional for phonons, though their complete dynamics also requires the kinetic operator arising from the conservative part of the action.
Owing to the spatial modulation of the designer field dynamics in the gravitational description, one finds there to be an interesting discontinuity in the dynamics between vanishing and non-vanishing spatial momentum. As explained at the outset, in order to isolate the physics of sound, it suffices to imagine there being an infra-red cut-off in momenta and to study modes which are long-wavelength above this cut-off scale. To understand the physics of the soft zero modes however, needs a bit more work. These modes can be understood as large diffeomorphisms of the background, but we have not attempted to quantize this system. It would be interesting to do so. Alternately, one could work with a physical cut-off, say by placing the plasma on a compact spatial volume, e.g., on S d−1 .
While our analysis was focused on sound modes in a neutral plasma, it can be readily extended to include additional conserved charges. For example, in a charged plasma we have fluctuations of both energy and charge density; while the former leads to sound, the latter leads to charge diffusion in the scalar sector. In this case, one has intricate dynamics where two long-lived modes are coupled to each other. It is nevertheless possible, as in the case of vector perturbations discussed in [7], to decouple these two modes and construct the desired effective action. Our preliminary investigations suggest that the general paradigm explained here continues to hold; we hope to report on this in the near future [34].
Finally, it is worth contrasting the analysis of real-time fluctuating hydrodynamics with the earlier work on the fluid/gravity correspondence [20,21]. The focus in that work and extensions thereof was to construct the gravity dual of a fluid flow of the boundary CFT. In particular, given a holographic system whose stress tensor one-point function can be parameterized in terms of hydrodynamic variables, viz., temperatures and velocities, obeying the conservation equations, the fluid/gravity paradigm constructs a spacetime geometry characterized by this hydrodynamic data. By virtue of focusing on thermal one-point functions, that analysis had a technical advantage of being able to work with SO(d − 1) tensor decomposition, but more importantly was fully non-linear in amplitudes of departures from thermal equilibrium. 23 There is a specific prediction for fourth order (in gradients) transport data contained in the q 4 terms of (D.9). We have not attempted to classify the terms in the Class L action that are responsible for it.
The open effective field theory paradigm however addresses a slightly different question: "What is the gravitational dual of a fluctuating plasma"? More precisely, realizing that the plasma consists of both short-lived and long-lived modes, we seek to parameterize its dynamics in terms of the sources for the former and the operators (or fields) corresponding to the latter. This was the philosophy outlined [6] for the study of the Wilsonian influence functional of the plasma. One has to not only keep track of the dissipative pieces which relate to infalling quanta in the dual gravity, but also the Hawking quanta that correspond to stochastic fluctuations. But this is precisely what has been achieved in terms of the designer fields, which now parameterize the fluctuating bulk metric. While they are not manifestly SO(d−1) covariant, and our analysis thus far has been restricted to linear order in amplitudes, the close resemblance of the ingoing part of our solutions to those obtained in the fluid/gravity literature makes it highly suggestive that it should be possible to bootstrap onto a non-linear solution. It would be fascinating if this goal can be realized.

A Dynamics of scalar gravitons
Our starting point for analyzing the action is simply the Einstein-Hilbert action with its Gibbons-Hawking variational boundary term and appropriate counterterms. We are going to be working to quartic order in gradients. A-priori we expect that we would need counterterms accurate to that order. However, as we shall see there are some additional subtleties in this system which will allow us to obtain certain finite results from the quadratic counterterms alone. Irrespective of this we will quote here the full counterterm action accurate to fourth order in boundary derivatives. 24 The gravitational dynamics we consider is governed by 25 We will first examine the equations of motion which we write as and then proceed to analyze the variational principle.

A.1 Gauge invariant data and time-reversal
To understand the dynamics of the scalar gravitational perturbations and deduce that the dynamics can be captured by a single field Z we will proceed in a series of steps. Our first task will be to identify the diffeomorphism invariant combinations of the metric perturbations for the ansatz (3.1). A natural way to capture this information is to look at the curvature tensors which we write in terms of orbit space tensors. It will be convenient to define a connection on this part of the geometry: 26 Some useful identities which we have used to simplify the expressions are: We start with the metric parameterized in terms of Ψ AB as presented in (3.1). For this geometry the gauge invariant combinations organized into the orbit space tensors are [24]: • An orbit space vector X a whose dual one-form has components • An orbit space symmetric traceless rank 2 tensor, Y ab , with components (A.6) 25 We eschew the overall normalization by 1 16πG N to keep the expressions simple. Boundary quantities will be obtained by multiplying by c eff at the end. 26 We will use lowercase early alphabet Latin characters to indicate orbit space tensors in addition to the conventions specified in footnote 5.

TR Parity Gauge invariants Metric components Debye gauge data
Even Odd • And finally, we have orbit space scalars To understand the time-reversal properties of these combinations we use the observation that on the orbit space time-reversal is just a diffeomorphism. Hence we conclude that Y a a and Y S have even time-reversal parity. For the reminder we use the fact that the orbit space vectors can be decomposed into the basis adapted to time-reversal introduced above (2.5), and use the fact that dv− dr r 2 f is odd under time-reversal and dr r 2 f is even. Similar decomposition for the tensors leads to For the purposes of analyzing the equations of motion it is helpful to define some rescaled combinations of fields which have definite time-reversal parity. We introduce: We summarize the essential information from this analysis in Table 2.

A.2 Dynamics in the Debye gauge
In Appendix A.1 we introduced the gauge invariant combinations of metric perturbations. One can however fix some of these metric functions by using the diffeomorphism freedom. A-priori we can gauge fix three functions, leaving behind four of the seven perturbation functions appearing in (3.1). We will implement this by working with a set of 4 functions {Ψ S , Ψ vv , Ψ vr , Ψ rr } by first rescaling out a factor of k 2 from the gauge invariant scalar and tensor data, Y S and Y ab , i.e., setting X a = Ψ T = 0. Equivalently, we have the gauge conditions We can interpret Y S and Y ab in terms of metric components in the scalar sector in a Debye gauge. 27 This is a coordinate chart such that metric has no derivatives of scalar plane waves under plane wave decomposition, viz., the perturbation can be recast into the form This was the gauge choice adopted in [24]. We will now present the linearized Einstein equations in terms of these scalars by decomposing (A.2) into plane waves. Employing the definitions in (A.10) and further introducing the combination: we end up the metric which at linear order takes the form: (A.14) With this choice of gauge the vector gauge invariants vanish X a = 0, while the remaining scalar and tensor invariants simplify in the parameterization (A.14) to It will be helpful to assemble the equations of motion (A.2) into time-reversal invariant orbit space tensor combinations as above. We first have the scalar equation, which involves only Φ B and takes a simple form: 28 This equation is actually a simple algebraic constraint on the invariants: The orbit space tensor equations assembled again into time-reversal invariant combinations take the form: We introduced here the quantity T, which is defined to be This leaves us with the vector equations which being coefficients of S i have an explicit momentum factor, keeping track of which will be important for understanding spatial zero modes. We find: (A.19) 28 With Ψ T = 0, this equation gets modified to The tilded equations strip out the momentum factor which is convenient to do. The remaining equations which are orbit space scalars picking out the trace and the S ij part of the spatial harmonics can be naturally expressed in terms of them as Finally, a natural way to combine the equations involves taking a particular combination of E B and E 5 : (A.21)

A.3 Parameterizing the solution space: k = 0
Since there are only four physical functions, we should only have to use four of the equations of motion. One can check that not all the equations given above are independent (explicitly visible for example in E 6 and E 7 ), which suggests that a judicious choice of four equations should suffice to distill the dynamics into a manageable form. For k = 0 an efficient choice turns out to be the set {E T , E 1 , E 2 , E 3 }, satisfying which will ensure that the reminder are also upheld. We will now analyze the equations introducing Θ and Z to simplify the dynamics in the process. Let us begin with E T = 0 which says that Φ B (r, ω, k) = 0, as long as we focus on non-zero k, spatially inhomogeneous modes. We will use this to set Φ B = 0 this subsection and return to the case where we have a spatially homogeneous function in Appendix E.
In the rest of the section we will give a brief discussion of how one simplifies the dynamics to that of a single scalar field. We have three independent linearized Einstein's equations in the set (A.17), (A.19) for the fields {Φ E , Φ O , Φ W }. It will be convenient to pick the following combinations as our independent Einstein's equations, setting Φ B = 0 in the process to simplify our expressions: We will see shortly that E 4 in (A.19) will be accounted for (actually it can be eliminated algebraically using an algebraic identity). The combination E 3 above will simultaneously take care of E 5 and E B by definition.
The Weyl factor and momentum flux fields: To solve these equations we adopt a strategy similar to the one employed in the analysis of gauge field equations in [6,7]. One notes that E 2 is the energy conservation equation; in fact Φ O is the only time-reversal odd field which is related to the momentum flux. This suggests we should algebraically solve this equation by letting Φ O ∝ ω. We express T in terms of the same variable and then fix Φ E using the first equation. To wit, This choice ensures that the first two equations in (A.22) are satisfied. We are then left with third equation E 3 , which can be viewed as a relation between Φ W and Θ. This gives a constraint on the parameterization, isolating the true dynamical equation.
The designer field for sound: At this point, based on the experience with vector polarizations and diffusive mode, one would expect that Θ is the physical variable that should parameterize the designer field dual to the sound mode in the plasma. While this is physically correct (as we will justify) there is however a technical obstacle. The parameterization (A.23) does not immediately give an autonomous equation for Θ but rather leads to a coupled system between Θ and Φ W from E 3 . One can however isolate a new field Z by realizing that Φ W and Θ are not independent but related to each other through the relation The first equality follows from (A.18) and the second from (A.23). This can be solved by introducing an auxiliary field Z and solving for Φ W and Θ in terms of it. 29 This results in the expression (3.6) quoted above. This explains the origin of the designer field Z and the momentum dependent factor Λ k , which originates during decoupling the Θ and Φ W dynamics. Once we arrive here, it is straightforward to check that the final constraint equation on the system is the equation of motion for Z given earlier in (3.9). It is easy to see that the resulting equation is second order once one appreciates that Θ, Z, and Φ W satisfy a linear relation from (3.6) To do so we use the observation that equations of the form (∂ + A)X = (∂ + B)Y can be solved by setting This allows us to write (A.26) A slight simplification of (A.26) using the explicit expression for Λ k and (A.4) leads to the equation of motion (3.9) quoted in the main text. We emphasize that the dynamics is governed by a time-reversal invariant equation, which as explained in [5], allows one to construct smooth solutions on the grSK geometry.
The parameterization of the metric functions in terms of Z is easily obtained to be This suffices to determine the linearized geometry (3.5) once we know the solution for Z. The change of variables involves Z, D + Z, and D 2 + Z and appears to be necessary to ensure that the classical phase space is only two-dimensional, parameterized by an effective source and a corresponding dual CFT plasma operator. We will focus on parameterizing the phase space for the present by boundary values of Z and subsequently argue that the physical solution space is best parameterized by the stress tensor component (T CFT ) i v or equivalently by Θ.

B Variational principle in the scalar sector
We have indicated in (A.1) that L will refer to the Lagrangian including the measure factor. Unless explicitly indicated, we will write the terms in the action in a series of steps below, quoting at each stage this Lagrangian in momentum space. Integrations over momenta and over the bulk radial coordinate can thus be avoided in the expressions, which themselves tend to be pretty long. We also use † to indicate the frequency and momentum reversed field, viz., Φ † (ω, k) = Φ(−ω, −k), and thus use + cc to account for symmetrization. This analysis is restricted to k = 0 as we seek to establish the variational principle for Z at the end of the day.

B.1 Action for time-reversal invariant fields
Since the background Schwarzschild-AdS d+1 solution has a non-vanishing on-shell action, when we expand the Einstein-Hilbert action with the perturbation ansatz, we will have terms starting at the zeroth order in the amplitudes of the perturbation. We will separate out the zeroth and first order contributions out ab-initio -they do not contribute to the dynamics of the linearized modes. Rather, these terms correspond to the background free energy and represent the ideal fluid contribution of the boundary action. 30 We therefore will write: We work with the fields {Φ E , Φ O , Φ W } having chosen to eliminate Ψ rr using (A.13) (and use (A.16) to set Φ B = 0). To begin, let us look at the contribution from the background and the linear terms in the fluctuations takes the form Up to this order, we can express this result as an ideal fluid action on the induced boundary geometry. This is similar to the earlier discussion in the vector and tensor cases [6]. One can equivalently write Here b µ is the (rescaled) thermal vector that picks out the inertial frame. We will later see that on our solution there are corrections to the thermal vector, which we will need to account for, to get the correct ideal fluid action at quadratic order as noted in footnote 30.
Turning to the quadratic part, we will proceed in a series of steps, outlining independently the contributions of the bulk Einstein-Hilbert and Gibbons-Hawking terms. This will suffice to demonstrate that the dynamics is governed by the familiar Dirichlet boundary conditions for the aforementioned fields (and hence for Ψ AB ). Finally, we will outline the contribution from the counterterms that render the on-shell action finite.
The bulk Einstein-Hilbert term can be shown to decompose into a bulk piece, a boundary term, and a total temporal derivative term, viz., 30 In the analysis of the tensor and vector modes in [6,7] we in fact even extracted a part of the quadratic terms which assembled nicely to give ideal fluid contribution at the outset. In the present case, given the relative complexity of the dynamics, we find it useful to keep the quadratic pieces together and only isolate the part which involves terms at most linear in the fluctuation fields.
We start with the bulk term which can be simplified to the form: This part of the action is obtained by direct evaluation and integrating by parts to isolate the boundary terms. It is interesting to observe that only Φ W has a quadratic kinetic term and Φ E appears in the kinetic part only coupled to Φ W . Since Φ O is the only time-reversal odd field its appearance in the action is highly constrained (and it only shows up with explicit time-derivatives). Note also that the field Φ W has a wrong sign kinetic term (from the last line) reflecting the familiar issue with the conformal mode in gravity. This already suggests that despite appearances Φ W is not the physical field. The temporal boundary term L dot EOW does not enter the analysis and can be dropped ab-initio. The radial boundary term cancels against a similar contribution from the Gibbons-Hawking term. The precise form of these terms will therefore not be necessary for us.
This is a good consistency check for our computation ensuring that the Einstein-Hilbert action together with the Gibbons-Hawking term has a good variational principle. The remaining part of the Gibbons-Hawking term turns out to be: 31 In [6,7] this contribution was referred to as the ideal piece, since in those cases, it corresponds to the bare ideal fluid action. We will refrain from employing that notation here; the ideal fluid contribution to the on-shell action is a bit more involved in the scalar sector, owing to the presence of a propagating mode.
The variation of (B.5) gives us three independent equations. One of these is the momentum conservation equation E 1 from (A.22) which comes from varying Φ O . The other two equations are linear combinations of the ones we have given above. It is interesting to note that the variation does not directly produce the E 3 equation which was crucial to derive the autonomous second order equation for Z. One aspect that is clear from the variational analysis is that the fields {Φ E , Φ O , Φ W } obey Dirichlet boundary conditions. This is manifest from the structure of the Gibbons-Hawking term which is a quadratic form in these three fields.
Finally, the counterterm action is given as (B.8) We will quote results accurate to quartic order in the gradient expansion for which it suffices to include the boundary counterterm that is quadratic in derivatives (i.e., it only includes the boundary Einstein-Hilbert term in (A.1)). We have included here the contribution from the quartic counterterm for completeness.

B.2 The designer scalar action
We would now like to distill the action in terms of the designer scalar field Z. To do so we can directly compute the terms the bulk Lagrangian and boundary terms from the Einstein-Hilbert dynamics defined in (A.1) with the metric ansatz (3.5) and expand to quadratic order. We could equivalently begin with the action given in terms of the {Φ E , Φ O , Φ W } fields and use the substitutions given in (A.27). This is a bit more useful, since we have already removed in the process redundant boundary terms. We can therefore focus on just the three terms computed earlier: L bulk EOW , L GH EOW , and L ct EOW . Direct substitution of (A.27) into the bulk term L bulk EOW leads to formidable expression, denoted as L EH [Z]. There is however a nice structure beneath this mess. Lets first see why L EH [Z] reduces into a two-derivative action in the bulk, with the complications relegated to the boundary terms. Since Φ E appears with a single radial derivative in (B.5) and the change of variables to Z involves a double radial derivative (A.27), we end up with a action with higher derivative terms. The highest derivative term is d 3 If we naively vary this action with respect to the field Z we expect to get a quintic order equation of motion which should be implied by (A.26), i.e., the resulting equation must be some combination of E Z and derivatives thereof. Carrying out the exercise we find however The fact that the higher order action leads to a second order equation of motion is a sign of the hidden simplicity. Once one knows this it is a matter of corralling the higher derivative terms and showing they are total derivatives. With some effort one can show that The bulk action in momentum space is given in (3.10) which we reproduce here The complicated boundary terms can be understood as follows. Firstly, the leading The Gibbons-Hawking term L GH EOW does not have a corresponding term with this high derivative order, but the counterterm does (from Φ E Φ † E ). The cleanest presentation of the boundary terms turns out to be to combine the contributions from L ∂ [Z] and L GH EOW and express the result as a general quadratic form in the variables D + Θ ∼ D 2 + Z, Φ W ∼ D + Z, and Z itself. We will refer to this total collection of boundary terms as the variational boundary terms of Z and write: We find with coefficient functions (B.14) One can check directly that Having dealt with the bulk and boundary terms let us turn to the counterterms. These can be evaluated by direct substitution, though we note that the presence of the Φ E Φ † E does mean that the counterterms are functionals of D 2 + Z, which is somewhat unusual. The counterterm action can be corralled into:

B.3 The variational principle for Z
We have all the pieces at our disposal to deduce the variational principle for the field Z. Our first task will be to work out the variational principle that leads to the generating function of boundary correlators, viz., the usual boundary conditions in the AdS/CFT parlance. We will then work out the appropriate Legendre transform that computes the WIF of the boundary theory from the grSK contour.
For the variational analysis we will treat the factor k 4 Λ 2 k as an overall pre-factor that we will account for at the end of the day. Equivalently, we work with an auxiliary system for Z where the action has this factor scaled out.
With this understanding let us first record the momentum conjugate the field Z. From (B.11) we find: Further, using the Green's function (4.6) we can check that the asymptotically the conjugate momentum is constant and parameterized by the inverse Green's function K s . With our conventions for G in Z we have We see that the conjugate momentum Π Z is finite; the ellipses is (B.19) denote the subleading terms of O ξ 2−d . This behaviour of the field Z and its conjugate Π Z is indeed what one expects from a non-Markovian field based on the analysis of [6]. At this point we can even guess that the boundary conditions for Z are Neumann for the purposes of computing correlation functions. We will however want to compute the Wilsonian Influence Functional (WIF) parameterized by the boundary moduli fields for Z, which will turn out to be computed by quantizing the field with renormalized Dirichlet boundary conditions.
We will now argue for this directly by analyzing the variational principle for the action S [Z]. Recall that we have organized the classical action for the designer field as (B.20) Let us begin with the bulk term whose variation is simply: The variation of the boundary term L var [Z] produces various terms which can be expressed a combinations of δZ, δD + Z ∝ δΠ Z and δ 2 D + Z. To deduce this we note that Θ and Φ W contain D + Z while D + Θ has a piece that behaves as D 2 + Z. Putting this together we expect that the stationarity of the action demands where A i depends on (Z, D + Z, D 2 + Z) and background metric data through the coefficient functions defined in (B.14). We factored out the k 4 Λ 2 k piece as advertised which will be helpful when we use our knowledge of the solution in the gradient expansion (at this stage it was not strictly necessary). This complicated second order boundary condition is what is necessary for ensuring the stationarity of the action at a generic radial hypersurface. We are however interested in an asymptotically locally AdS geometry, so we should understand what the behaviour of the boundary condition is in the r → ∞ limit (see similar discussion in [7]). For this we need not only the asymptotics of the coefficient functions c i in (B.14) (which is clear from their definition), but also the large r behaviour of the functions Θ, Φ W , D + Θ. These are of course easy to extract given the solution (4.6) for Z.
Generally, for the purposes of the variational principle it suffices to focus on the bulk and the boundary terms as we have done above. One ignores the counterterms -they are important to ensure that we have a finite norm on the phase space. Crucially, they should be expressed in terms of the intrinsic data on the boundary that are being held fixed by the boundary conditions. Since we inherit the counterterm action from those for the Einstein-Hilbert dynamics (A.1) they are naturally expressed in terms of variables that are held fixed. This is easy to see from our discussion of the dynamics in terms of the triple {Φ E , Φ O , Φ W } in Appendix B.1. We have a complicated boundary action as evidenced above, and if we wish to discern the asymptotic behaviour of (B.22) it will be helpful to work with a regularized phase space. Consequently, we will include the counterterms in our variation and discern what they tell us about the variational principle asymptotically.
Before we evaluate the variation of the boundary terms and counterterms it will helpful to corral them into a nicer form. Using (B.13) and (B.16) we find: We now simplify this in the following steps: • The coefficient functions are explicitly known, so we can use the background data to estimate the leading large r behaviour. To do so we must carefully extract the pieces that scale as r −2(d−1) since the functions are highly divergent. Carrying out this exercise we notice that some terms drop out (eg., the Φ † W Φ W term). • Next using the asymptotics of the solution given in (D.10), (D.12), and (B. 19) we learn that at large r This implies that we can for purposes of the large r behaviour carry out the replacements above. We will not directly eliminate all occurrences of D + Θ since the relation to Π Z involves a factor of k 2 .
Implementing this we find the boundary action (B.23) reduces to The two terms that are indicated above are the leading contribution near the boundary while the terms we have dropped have subleading pieces that serve as counterterms. For the purposes of ascertaining the variation principle the terms we have retained suffice. With this simplified boundary term in hand, we can now understand the full variational principle. Taking a variation of (B.25) and including the bulk contribution (B.21) we find: Finally, using (B.24) we conclude that the first term is subdominant and the total variation of the action is proportional to Z δΠ † Z . 32 Thus, the conventional boundary conditions used to 32 Strictly speaking this is not necessary; from (B.24) we note that δD+Θ ∼ δΠ Z so the total variation is indeed proportional to δΠ † Z implying the Neumann boundary condition deduced above.
compute the generating function of stress tensor correlators is a Neumann boundary condition for Z. We emphasize that this is a highly non-trivial statement relying on the nature of the asymptotic fall-offs for the fields involved. At finite radial position we would have a complicated mixed boundary condition fixing some relation of the schematic form given in (B.22). In fact, this is what would have been suggested if we examined the field redefinitions (A.27) and used the fact that {Φ E , Φ O , Φ W } have Dirichlet boundary conditions imposed on them. It is also not something we have imposed by hand, but rather it is completely inherited from the original gravitational dynamics (A.1). 33

C Boundary observables
We now have understood the dynamics of the gravitational system encoded in Z. This information can be used to decipher the boundary observables directly. We describe the computation of the on-shell action, the boundary stress tensor, and then turn to a brief discussion of the relative merits of the field Z versus Θ. We will continue to drop an overall factor of c eff which we have restored in the main text.

C.1 The boundary sources and operators
We have determined that the Z should be quantized with Neumann boundary conditions for purposes of computing correlation functions. This suffices for us to see that Z does behave like a regular non-Markovian field introduced in [6]. The field Z should, by virtue of this boundary condition, limit to the dual boundary operatorȎ Z . Taking the limit on the grsK geometry leads to the statements asserted in §4.3, in particular, (4.12) and (4.13).
The field Z is only divergent starting at O r d−4 , which is lower than what one would expect for a field that is supposed to encode the physics of the boundary stress tensor. A consequence of this is that the renormalized field Z is not modified by the boundary counterterms up to the quartic order (i.e., it is uncorrected by the boundary cosmological constant and Einstein-Hilbert counterterm). It only gets renormalized by the quartic R 2 counterterm in (A.1) which is contained in the coefficient b 0 at O k 4 in (B.17). Taking this into account we learn that the correction comes from the k 4 contribution to b 0 in (B.17). The renormalized field Z ren can be determined to be (C.1) Having identified the boundary conditions for Z and its renormalized counterpart Z ren , let us turn to identifying the source. Since Z is quantized with Neumann boundary condition, the source should be defined in terms of the conjugate momentum Π Z . However, the relations in (B.24) suggests that the conjugate momentum can be traded for the fields D + Θ r d−1 and Φ W r d−2 . This is consistent with the fact that these fields determined the induced boundary geometry. Indeed, including the background piece the induced metric on the boundary evaluates (on either boundary) to It therefore makes sense to identify the temporal and spatial components of the boundary metric as the sources. These are however not independent from each other or from that conjugate momentum Π Z owing to (B.24). We define therefore the boundary source as in (5.5) and note that This definition of the boundary source for Z in terms of D + Θ has a proper gradient expansion, unlike Π Z which has an additional k 2 . As noted above, it is also the physically correct variable; the temporal component of the boundary metric Ψ vv that couples to the energy density is indeed Φ E = D + Θ up to a factor of r d−1 . Fixing all the normalization factors we find that the induced boundary metric can be expressed as in (5.7).

C.2 The boundary stress tensor
The boundary stress tensor density is given by varying the boundary Gibbons-Hawking term and the counterterms given in (A.1). This leads to the following expression accurate to quartic order in gradients We will now present the result for (T CFT ) ν µ which makes it easier to see the traceless condition by inspection. At the first order in amplitudes one evaluates the components of the stress tensor from the Brown-York analysis supplemented with counterterms (C.4). We quote the results for the individual components in turn.
to simplify the answer. The first summand in T P simplifies to Φ E 2 r d−1 , which is a source term, while the designer field Z assembles from the subsequent pieces involving Φ O and Φ W which constitute part of the operator contribution. One gets a term proportional to ω 2 Z from them. The contribution from Φ E − (d − 2) r Φ W is however also an operator contribution, the leading order terms in this difference cancels, as the reader can verify from (D.10) and (D. 12).
With this information we can evaluate the stress tensor on our solution parameterized by the designer field. The result of this exercise is what is reported in (5.10) where we have adhered to the identification of source and vev terms. Specifically, contributions of the form Φ E r d−1 and Φ W r d−2 are written in terms ofζ using (C.3).

C.3 The on-shell action
Once we have identified the boundary conditions for Z, we can evaluate the on-shell action on the grSK solution. At quadratic order this is just a boundary term and can be easily evaluated on the grSK geometry. One will obtain from this the generating function of boundary correlation functions with sources on the two boundaries of the grSK geometry, viz., We would however like to evaluate the Wilsonian influence function, for which we need to perform a Legendre transformation of the generating functional. This requires one to evaluate the on-shell action after including a suitable boundary term to carry out the Legendre transform: As noted in [6] this amounts to quantizing Z without the additional variational boundary term, i.e., we quantize Z with (renormalized) Dirichlet boundary conditions as have noted in (C.1). We will find two distinct contributions to the on-shell action: one arises from terms of the form Π Z Z † and originates from a combination of bulk action, various boundary terms, and the Legendre transform. The other contribution will turn out to be purely a functional of the source originating from the D + Θ † D + Θ in (B.25). Using (5.5) this piece can be written as a factorized source contribution on the two boundaries of the grSK geometry. We end up with the result quoted in (5.1) in the form of a contact term (the source contribution) and the genuine influence functional.
The explicit evaluation of the on-shell action is straightforward since we have already deduced the asymptotic behaviour of the field and the conjugate momentum. Let us start with the non-contact term and record the influence functional. The grsK solution for Z is given in (4.14). Accounting for the contribution from the Legendre transform we find this evaluates to where we introduced a large radius regulator at r = r c . This is the result quoted in (5.2). The contact term may likewise be evaluated directly. From the definition of the sources (5.5) and the form of the boundary terms in (B.25) we see that this is simply given in terms of the boundary source. The contributions furthermore factorize leading us to the following expression at quadratic order: This contact term contribution is quite peculiar. However, as we describe in §5.4, it is nothing but the contact term part of an ideal fluid on the boundary geometry (C.2). We also explain there how one can isolate various hydrodynamic transport data from our answer and connect to the discussion of Class L adiabatic fluid lagrangians of [31].

C.4 On field redefinitions and boundary operators
The solution for the designer sound field Z on the grSK contour can be repackaged directly in terms of field theory data. Recall that we expect a single mode of the boundary stress tensor that captures the effective dynamics in the low energy limit. The conservation equation of the stress tensor turns out to be a constraint on the field Θ, cf., the discussion around (A.22) and (A.23). Equivalently, examination of the induced boundary stress reveals that the leading contribution comes directly from Θ, The full expression for the stress tensor including the counterterms can be found in (C.5), where the √ f Φ W term can be seen to arise from counterterm contributions. These arguments suggest to us that the holographic dual of the sound mode in the plasma should be identified with Θ. However, Θ does not by itself have simple dynamics. 34 Consequently, we rely on Z as an intermediate auxiliary field to analyze the problem and translate the physical data back onto Θ therefrom.
The relative choice between Z and Θ is effectively a field redefinition in the boundary. To appreciate this, let us obtain the grSK solution for Θ by first constructing the inverse Green's 34 We were able to derive a third order radial differential equation for Θ directly. At each order in the gradient expansion this equation turned out to be a second order inhomogeneous equation for d dr Θ which suggests again that there is further simplification possible by passing onto Z. function G in Θ for Θ which is reported in (D.11). We find Θ SK (ζ, ω, k) = G in ΘZ a + n B + The key point to note is that the coefficient of radially homogeneous mode G in Θ from (D.11) is which suggests thatȆ One could, if one wished to do so, convert the expressions in the main text directly to expressions involvingȆ, but we have refrained from doing so to avoid complicating the already involved discussion.

D Further details of the gradient expansion solutions
The solution for the designer field Z which satisfies (3.9) with ingoing boundary conditions was given in §4.2. This form was chosen to make direct contact with the functions appearing in the fluid/gravity literature [33] and the earlier analysis in [6]. We compile some useful results about the functions appearing in the expansion in this appendix. The functions {F, H ω , H k , I ω , I k , J ω , J ωk , J k } that parameterize the solution for Z in (4.6) have compact integral expressions tabulated in Table 1. In Appendix D.1 we collect several useful facts about them and determine their asymptotic behaviour. Using this data we record in Appendix D.2 the asymptotic expansions for the fields Z, Θ, and the metric functions Φ E , Φ W which will prove useful for the phase space and boundary condition analysis. Subsequently, in Appendix D.3 we present the solution in the alternate form parameterized in [7] for ease of comparison.

D.1 Asympotics of the solution
As noted in §4.2 the functions parameterizing the solution for Z in gradient expansion are cleanly written in terms of a double-integral transform of a source function J, cf., (4.11). The sources for the various functions are collated in Table 3. Examining this data we immediately see that there are some useful relations: • First, the sources for the functions H ω and H k determine an useful identity for I k : H ω (ξ) + H k (ξ) = d − 2 2 I k (ξ) . (D.1) with T g g Asymptotics Table 3: The functions appearing in the gradient expansion of the Markovian ϕ d−1 up to the fourth order in gradients, given in the form of an integral transform defined in Eq. (D.17). We also present the leading asymptotic behaviour of the functions which is used for computing boundary observables. computed earlier, but had to introduce two additional functions: (D. 20) Of these, only ∆ 0,4 Z (ξ) has a divergent behaviour at large ξ. It is defined using the earlier functions in (4.3) by analytic continuation. In this parameterization, both ϕ 2,2 Z (ξ) and ϕ 0,4 Z (ξ) have divergent corrections (to quartic order), while that used in § 4.2 the divergence was isolated into V k (ξ).

E Spatially homogeneous modes
This appendix is somewhat outside the main line of development of the paper and is included for completeness. As we saw above in Appendix A.3 the dynamics of spatially inhomogeneous modes can be distilled into that of a single non-Markovian scalar Z satisfying a second order non-Markovian differential equation. We wish now to analyze what the equations of motion imply for spatially homogeneous modes.
We will carry out the analysis in two steps. First we examine the large diffeomorphisms of the Schwarzschild-AdS d+1 solution that respect the Debye gauge choice. Subsequently, we look at the solution space for spatially homogeneous modes looks like and parameterize it in terms of the most general allowed data compatible with asymptotically locally AdS asymptotics.
We will find that the two sets of analyses lead to the same set of zero modes. The surprise will be that there are more zero modes than those that can be lifted to physical moduli captured by Z.

E.1 Large diffeomorphisms of the background
The background Schwarzschild-AdS d+1 geometry is parameterized by b which is a measure of the black hole temperature or mass. One has additionally chosen a particular Weyl frame, by making a suitable choice of the radial coordinate r. Consider the now the following diffeomorphism and parameter shift on the orbit space, leaving spatial homogeneity intact: We have two functions on the orbit space and two constants C x and C b which the reader will recognize is precisely the freedom to rescale spatial length scales and the boundary temperature homogeneously.
To check this explicitly, we simply implement this change on the background solution, We give now a short complementary perspective on the spatially homogeneous modes from the dynamical equations of motion. As one might anticipate the physical solution space is already fully characterized by the large diffeomorphism modes, so we expect to see the same degree of freedom in the solution space, as we shall verify below. The manner in which this happens is that the dynamics of the system is modified at k = 0 resulting in additional moduli in the problem. Technically, at k = 0 some of the equations degenerate. For one, the scalar equation E T , the spatial vector equations E 4 and E 5 , and the spatial tensor equation E 7 are trivially satisfied, each being explicitly proportional k i . We are then left with a simpler set of equations, which we write as 36 where we have written the last equation succinctly using the parameterization defined in (A.19) for convenience. The first two equations in (E.7) imply that the combination T 0 k − r Φ 0 k B must be a constant, We are then left with two equations for effectively three variables. However, the remaining two equations are not independent, for The one remaining equation can be written as The solution space is parameterized by two functions Φ 0 k W and Φ 0 k O . The other two functions {Φ 0 k E , Φ 0 k B } can be solved in terms of them using (E.8) and (E.10). Inspired by (E.4), 36 We use the superscript 0 k to remind us that we are looking at spatially homogeneous modes of our fields.
we can w.l.o.g parameterize the functions Φ 0 k W and Φ 0 k where we demanded that the solution be asymptotically AdS.
Here C x is a coefficient that will related to integration constants momentarily. We recover the other two functions as given in (E.4) as 37 (E.12) One can readily see that we may identify F 1 = r d−3 χ r and F 2 = −χ v + χr r 2 f recovering the background diffeomorphisms. This implies that the time-independent pieces in Φ 0 k B and Comparing with (E.4) we conclude that we can parameterize the integration constants by defining them in terms of C x and C b : (E.14) We thus recover the full set of large diffeomorphisms from the analysis of the equations of motion.
Zero modes and designer field: One can ask if this solution is related to the zero momentum solution of the designer field Z. The ingoing Green's function given in (4.6) requires C T = C B = 0, which can be inferred by noting that the source terms vanish at ω = 0. In addition, we also note that Φ B vanishes identically in our parameterization by Z, which demands a relation between F 1 and F 2 , viz., This is consistent with the relations (A.27) between Φ E , Φ O , Φ W and Z at k = 0, which in turn require that Z is given in terms of the diffeomorphism functions F 1 and F 2 as 16) 37 The operator D+ annihilates e − 1 2 w ζ but this solution does not satisfy ingoing boundary conditions and hence we restrict to allowing an integration constant CB.
Solving for F 1 and F 2 using (E.15) and the first equation of (E.16), we can then write an autonomous equation for Z 0 k . The resulting equations turns out to be implied by the Markovian wave equation (with k = 0) which is of course the zero momentum limit of our designer field equation (3.9). This is also the zero momentum limit of a minimally coupled massless scalar field in Schwarzschild-AdS d+1 .
We thus see that a part of the large diffeomorphisms is indeed the homogeneous solution for the designer field. In particular, given a solution for Z(r, ω, k) satisfying (E.17), we can determine large diffeomorphism functions F 1 and F 2 . However, one does not recover the full set of large diffeomorphisms. Specifically, the part parameterized by the constants {C T , C B } (or equivalently {C x , C b }) is not recovered from the designer field dynamics.