$U(1)$ current from the AdS/CFT: diffusion, conductivity and causality

For a holographically defined finite temperature theory, we derive an off-shell constitutive relation for a global $U(1)$ current driven by a weak external non-dynamical electromagnetic field. The constitutive relation involves an all order gradient expansion resummed into three momenta-dependent transport coefficient functions: diffusion, electric conductivity, and"magnetic"conductivity. These transport functions are first computed analytically in the hydrodynamic limit, up to third order in the derivative expansion, and then numerically for generic values of momenta. We also compute a diffusion memory function, which, as a result of all order gradient resummation, is found to be causal.


Introduction and Summary
Hydrodynamics [1][2][3] is an effective long-distance description of most classical and quantum many-body systems at nonzero temperature. Within the hydrodynamic approximation, the entire dynamics of a microscopic theory is reduced to that of conserved macroscopic currents, such as expectation values of energy-momentum tensor or charge current operators computed in a locally near equilibrium thermal state. An essential element of any hydrodynamics is a constitutive relation which relates the macroscopic currents to fluid-dynamic variables (fluid velocity, conserved charge densities, etc), and to external forces. Derivative expansion in the fluid-dynamic variables accounts for deviations from thermal equilibrium. At each order, the derivative expansion is fixed by thermodynamics and symmetries, up to a finite number of transport coefficients such as viscosity and diffusion coefficients. The latter are not calculable from hydrodynamics itself, but have to be determined from underlying microscopic theory or experimentally.
The most simple example of constitutive relation is the diffusion approximation for electric current J = −D 0 ∇ρ (1. 1) where ρ is the corresponding conserved charge density. Diffusion equation which follows from the current conservation has a well-known major conceptual problem-it violates causality. Relativistic Navier-Stokes hydrodynamics is plagued by similar pathology. In the case of the latter, this problem reveals itself in the form of numerical instabilities, rendering the framework into practically useless.
To restore causality one has to introduce higher order gradient terms in the constitutive relation. Moreover, inclusion of a finite number of terms is insufficient, and causality is restored only after all (infinite) order gradients are resummed. Generically, this is equivalent to a non-local constitutive relation of the type: where D is a memory function, which is most generally non-local both in time and space. Causality implies that D has no support for t < t . In practice, the memory function is typically modelled: the simplest model imposes an exponential relaxation time approximation D(t − t ) ∼ e −(t−t )/τ θ(t − t ). Müller-Israel-Stewart [4][5][6][7] second order relativistic hydrodynamics, commonly used in simulations of heavy ion collisions, is an example of such a modelling 1 . AdS/CFT correspondence [9][10][11] opens a possibility to study transport properties exactly and partially analytically, at least for a class of quantum gauge theories for which gravity duals can be constructed. This holographic duality maps hydrodynamic fluctuations of a boundary fluid into long-wavelength gravitational perturbations of a stationary black brane in asymptotic AdS space [12][13][14]. Viscosity and many other transport coefficients could be computed from the gravity side of the correspondence. The ratio of viscosity η 0 to the entropy density s was computed in [12,13,15] η 0 s = 1 4π (1. 3) and was found to be universal for all gauge theories with Einstein gravity duals [16][17][18].
Remarkably, the fluid/gravity correspondence is not limited to linear response theory for small perturbations of the velocity field u µ . Navier-Stokes equations are completely encoded in the Einstein equations. Particularly, the formalism of [19] provides a systematic framework to construct nonlinear fluid dynamics, order by order in the fluid velocity derivative expansion, with the transport coefficients determined from the gravity side. The study of [19] was subsequently generalised to conformal fluids in higher dimensions [20], weakly curved background manifolds [21], and to forced fluids [22]. We refer the reader to [23][24][25] for comprehensive reviews of fluid/gravity correspondence. In [26,27], two of us built upon previous work [28] and constructed a flat space allorder linearly resummed relativistic conformal hydrodynamics using the fluid/gravity correspondence. By linear in the fluid dynamic variables, we mean small amplitude terms like (∇∇ · · · ∇u) µν ; but with all nonlinear structures such as (∇u) 2 µν neglected. This is mathematically a well controlled approximation. The relativistic hydrodynamics with all order derivatives resummed was found to have a rich structure, absent in a strict low frequency/momentum approximation. The fluid stress-energy tensor was expressed using the shear term with η 0 replaced by a viscosity functional η of space-time derivative operators, and a new viscous term, which emerged starting from the third order in the gradient expansion. In Fourier space, both transport coefficient functionals became functions of frequency and spatial momentum, η(∂ t , ∂ 2 ) → η(ω, q 2 ). We performed a similar study for the Einstein-Gauss-Bonnet gravity in [29]. Ref. [30] generalised computations of [26,27] by including small non-dynamical metric perturbations in the boundary stress tensor. It was found that four new transport coefficient functions are needed in order to linearly couple the stress tensor to external gravitational perturbations. In [28] these transport coefficient functions were called gravitational susceptibilities of the fluid (GSF). They play a role similar to electrical conductivity, which will be discussed below. After the resummation, the memory function associated with the viscosity term was found to be consistent with the causality constraints.
Our strategy in [26,27,29,30] was to further develop the method of [19,22] by including all order linear structures in a self-consistent manner. Interestingly, we found that the dynamical components of the bulk Einstein equations are sufficient to compute an "off-shell" stress-energy tensor of the boundary fluid, with the transport coefficient functions fully determined. Additionally imposing four remaining constraints from the set of all bulk equations is equivalent to putting "on-shell" thus obtained stress-energy tensor.
In the present work, we extend our study into dynamics of a conserved vector current associated with a global U (1) symmetry, focusing on the properties of all order charge diffusion and conductivity. The most general linear in the fields (small amplitude) covariant form constitutive relation extending (1.1,1.2) is 2 where u α = (1, 0, 0, 0) and P µν = η µν +u µ u ν is a projection operator onto spatial directions. The field strength tensor F µ is associated with a non-dynamical external electromagnetic potential A (0) µ (see equation (4.4)). More transparently, the constitutive relation (1.4) reads The transport coefficients D, σ e and σ m are scalar functionals of spacetime derivatives In the Fourier space, via replacement (∂ t , ∂ 2 ) −→ (−iω, −q 2 ), these functionals of the derivative operators turn into complex functions of frequency ω and momentum q: 2 A very similar constitutive relation was deduced from an effective action approach in the most recent paper [31], where the chemical potential (instead of the charge density ρ) was used to express J µ . However, no values for the transport coefficients were determined there.
In the hydrodynamic limit, the transport coefficients are expandable in Taylor series, each term corresponding to a fixed order gradient (1.9) D generalises the usual concept of diffusion constant D 0 and will be referred to as diffusion function. σ e is an extension of electric DC conductivity σ 0 e . The second conductivity σ m reflects a response to external magnetic field. This second derivative structure (the curl of magnetic field) appeared already in [32,33] but no information about σ m was provided there. τ D ,e,m are corresponding relaxation times. The constitutive relation (1.5) is an offshell construction. Nevertheless, for a specific underlying microscopic theory, we are able to determine these transport coefficient functions. Putting the current on-shell reproduces the usual linear response theory, such as the Ohm's law 3 . While the magnetic term drops from the continuity equation σ m does contribute to current-current two-point correlators (see (5.6)). Our construction is quite analogous to that of all-order stress-energy tensor derived in [30]: D is analogous to the viscosity η, while σ e , σ m play the role similar to the GSFs in the constitutive relation for the stress tensor. Our goal below is to compute D, σ e , and σ m for a specific 4d theory defined holographically: the bulk theory is a probe Maxwell field in the background of Schwarzschild-AdS 5 geometry. This theory is a consistent reduction from probe D-brane models such as D3/D7 construction of [34]. In quenched approximation, i.e., when the number of probe D-branes is much smaller than that of colour branes responsible for the bulk geometry, it is reasonable to neglect the back-reaction of the bulk Maxwell field on the metric. Thus we study charge dynamics decoupled from any energy-momentum flow. Coupling the dynamics of the U (1) current to that of the stress tensor is certainly interesting and is likely to reveal new transport structures. We plan to explore this in the future.
R-current of the boundary CFT is a particular example of the U (1) current we are interested in. Its hydrodynamics was previously studied via holography in [13,[35][36][37][38][39][40][41][42][43][44]. For various supergravity backgrounds, ref. [36] derived a simple expression for the diffusion constant D 0 in terms of the metric components. In [17,45,46], black hole membrane paradigm was used to explore universality of the diffusion constant and DC conductivity. For N = 4 SU (N c ) SYM plasma in the large N c and large 't Hooft coupling limit, the Rcharge diffusion constant D 0 = 1/(2πT ), DC conductivity σ 0 e = πT and τ e = log 2/(2πT ). For convenience, below we will adopt dimensionless ω and q 2 normalised to π times the temperature, πT = 1. The physical frequency and spatial momentum are πT ω and πT q.
It is important to stress that the holographically defined boundary theory does not have a dynamical U (1) gauge field. Hence the background external fields ( E and B) are by construction non-dynamical. However, one could imagine gauging this global U (1) symmetry [37,40] with a small coupling. Then the current J µ can be treated as an induced current in a conducting medium described by a self-consistent macroscopic electrodynamics (see section 6).
Our derivation consists of several steps. We first solve the bulk Maxwell equations in a static AdS 5 black brane geometry. An important element of the procedure is that we solve the Maxwell equations for a set of unspecified boundary conditions A (0) µ imposed at the conformal boundary. The boundary field A (0) µ acts as a source of the current J µ . All order transport coefficient functions are then read off from the near boundary behaviour of the bulk Maxwell field. In the hydrodynamic limit ω 1, q 1, our analytical results derived in subsection 5.1 are −π 2 ω 2 + q 2 (6 log 2 − 3π) + · · · , σ e = 1 + log 2 2 iω + 1 24 π 2 ω 2 −q 2 (3π + 6 log 2) + · · · , σ m = 0 + 1 16 iω 2π − π 2 + 4 log 2 + · · · . (1.11) These results fully agree with all previous studies (and extend beyond): the constant term in D (D 0 ) was first reported in [13]; the DC conductivity σ 0 e is consistent with that of [37] obtained using the low-frequency results of [13]; the expression for σ e at q 2 = 0 is in agreement with [42,43] (see also subsection 5.2 for a more comprehensive comparison extending to generic values of ω). The underlined terms are new results unavailable in the literature. Quite interestingly, the hydrodynamic expansion of σ m starts from a term linear in ω, rather than from a constant.
For generic values of momenta, the transport coefficient functions are computed numerically. The results are presented and discussed in subsection 5.2. Particularly, the diffusion function vanishes at very large momenta, which is a necessary condition for causality restoration. The memory function introduced in (1.2) is defined as Fourier transform of the diffusion function D D t, D is also computed in subsection 5.2 and it is found to be causal D(t) ∼ Θ(t). For positive times, D(t) displays damped oscillations, which reflects the presence of a (infinite) set of complex poles-the bulk's quasi-normal modes. In the next section 2, the meaning of all-order gradient resummation and its implication on the dynamics are briefly clarified, essentially quoting the discussion presented earlier in [30]. The holographic model is introduced in section 3. Section 4 presents our main calculations. Maxwell equations in the bulk are discussed first. We search for solutions which are functionals of the boundary data ρ(x α ) and A µ (x α ). Decomposing these solutions in a basis constructed from the boundary fields, we are able to cast the problem into a set of ordinary differential equations for the decomposition coefficients. Then, we read off J µ from the near-boundary expansion of the decomposition coefficients. The latter are computed first analytically and then numerically in section 5. In subsection 5.3 we derive the continu-ity equation from the constraint component of the bulk Maxwell equation. Current-current two-point correlators are also expressed in terms of the transport coefficients D, σ e,m .
In section 6, we gauge the external field, turning it into a dynamical Maxwell field on the boundary. Combined with the current constitutive relation (1.5) this becomes a selfconsistent electrodynamics of a conducting medium. For completeness of our presentation we compute the medium's dielectric functions (transverse and longitudinal) and remark on negative refraction phenomena. In section 7 we make a brief summary. The calculations presented in section 4 are derived in Landau frame. In Appendix A we demonstrate that in fact the results are frame independent. Perturbative solutions for the bulk Maxwell fields are given in Appendix B.

All order gradient resummation
In this section we provide a clarification about what we actually mean by all-order gradient resummation, which along spacial derivatives includes an infinite number of time derivatives. This particularly means that the effective dynamical equations require an infinite set of initial conditions or, equivalently, the dynamics at hand is a theory of infinite number of degrees of freedom. Indeed, these are the quasi-normal modes of the bulk Maxwell theory. An alternative way to see the origin of the problem is to realise that exact dynamics in the bulk cannot be mapped onto a single degree of freedom on the boundary.
Yet, it turns out that it is possible to formulate this all-order hydrodynamics as a normal initial value problem with all the time derivatives absorbed into the memory function. For the sake of argument we will ignore any space dependence and will set B = 0. The causal current (2.1) andσ e is a memory associated with the conductivity function. The external electric field E is normally turned on at negative times, so to create an initial charge density profile at t = 0, and then turned off at t = 0 (E(t > 0) = 0), letting the system to freely relax to its equilibrium at infinite future. For such experimental setup, for positive times the current Generically, J H is not vanishing and it accounts for the entire history of the system at negative times. This is in contrast to a typical memory function-based approach, where one introduces constitutive relation (2.2) assuming J H = 0 and then also models D [2,3]. Our construction is so far formally exact. However, in order to solve the dynamical equation (1.10), it is not sufficient to provide the initial condition for the density only, but we also need the "history" current J H at all times, equivalent to providing infinitely many additional initial conditions.
We are now to discuss under what conditions we can nevertheless set J H to zero, casting our theory into a well-defined initial value problem. The response functions D and σ e are some given functions defined by underlying microscopic theory. Thus the equation J H (t > 0) = 0 is in fact an equation for the electric field E at t < 0. It is, however, not obvious that there exists a solution for generic D andσ, because we want the current J H = 0 at all positive times. Even though we cannot guarantee vanishing of the current identically, it is safe to assume that its effect could be rendered negligibly small. Particularly, it is obvious that J H vanishes at late times and its only potential influence could be at very early times when the current J J H .

The holographic model
We would like to study a finite charge density/chemical potential system having a global U (1) symmetry exposed to an external electromagnetic field. The holographic model is a Maxwell field in the background of Schwarzschild-AdS 5 geometry.
Conventions: the upper case Latin indices {M, N, · · · } denote the bulk while the lower case Greek indices {µ, ν, · · · } the boundary coordinates; the lower case Latin indices {i, j, · · · } indicate the spatial directions on the boundary; g M N is the bulk metric; γ µν is the induced metric on a hyper-surface Σ of constant radial coordinate r, and η µν is a Minkowski metric at the boundary. The gauge coupling of the bulk Maxwell field is set to unity.
The regularised bulk action is is a counter-term action to be specified below. From the action (3.1), the Maxwell equations are Radial gauge A r = 0 is used throughout this work. In the ingoing Eddington-Finkelstein coordinates, the Schwarzschild-AdS 5 geometry is where f (r) = 1 − 1/r 4 . The horizon is set to be at r = 1, which normalises the Hawking temperature to πT = 1.
Holographic renormalisation of the U (1) gauge field in the asymptotically AdS 5 space was considered previously using the Fefferman-Graham coordinates [47][48][49]. We are going to revise these studies for the case of the Eddington-Finkelstein coordinates. First, the Maxwell equations (3.2) are assumed to be solved. Near the conformal boundary r = ∞ the solution is expandable in a series where A (1) µ and B (2) µ are fixed by the Maxwell equations (3.2) The values of A (0) µ and A (2) µ are left unspecified. These are the equation's integration constants. Second, the near-boundary expansion for A M (3.4) is substituted into the Maxwell action of (3.1), resulting in a logarithmical divergence near the conformal boundary r = ∞, where "· · · " denotes a finite piece when Λ → ∞. S reg. is a regularised action at the cutoff Λ. Its covariant form makes up the counter-term in (3.1) where the indices are contracted with the metric γ µν on Σ: In writing down the counter-term (3.7) only the minimal subtraction scheme is implemented, that is (3.7) has no finite piece. Keeping some in the counter-term would induce certain contact terms (see e.g. [50,51] for discussion of ambiguities in the renormalisation scheme).
The boundary current is defined as a variation of the on-shell action S with respect to the source: The variation of S is computed using equations of motion (3.2) and assuming A µ falls off sufficiently quickly at infinity where ∇ µ is compatible with γ µν , and n M is the unit out-pointing vector normal to Σ, The boundary current is Expressed in terms of the near-boundary asymptotic expansion, (3.4) becomes (3.13) -8 -

Boundary current from the Maxwell dynamics in the bulk
The five Maxwell equations (3.2) split into four dynamical components EQ µ = 0 and one constraint EQ r = 0. Following the strategy developed in [26,27,29,30], we will initially solve the dynamical equations only. Substituting the solution into (3.13) will result in an "off-shell" current J µ with values of transport coefficients being fully determined. The constraint EQ r = 0 is then found to be equivalent to the continuity equation (1.10). We first notice that the equations EQ µ = 0 admit the most general static homogeneous solutions where A (0) µ and ρ are, for the moment, assumed to be some constants. Substituting (4.1) into (3.13), the boundary current corresponding to this solution is just Hence we identify ρ with the expectation value of the charge density operator. The boundary theory is a static uniformly charged plasma with no external fields present (F (0) µν = 0). Next, following the spirit of [19], we promote the constant parameters A (0) µ and ρ into unspecified functions of the boundary coordinates, A µ of (4.1) would cease to solve the Maxwell equations for arbitrary functions A µ (x α ) and ρ(x α ). The solution has to be amended: For given A µ (x α ) and ρ(x α ), the Maxwell equations (3.2) will be solved for the correction a µ . The equations are subject to specific boundary conditions to be discussed next. Near r = ∞, we require that a µ (r → ∞, x α ) → 0, (4.5) so that A µ (x α ) could be identified with the external electromagnetic field sourcing the current J µ . We would like to keep identifying ρ(x α ) as the charge density profile, equivalent to Landau frame convention The Landau frame convention can be argued to be equivalent to residual gauge fixing. In Appendix A we relax this condition and demonstrate that the results on all the transport coefficients are still uniquely fixed 4 . Finally, since we are working in the ingoing Eddington-Finkelstein coordinate, we further require that a µ (r, x α ) is regular over the whole range of r. Particularly, the regularity requirement at the horizon r = 1 is equivalent to the in-falling wave condition of [52] in the Poincare patch. The above boundary conditions are sufficient to fix a µ uniquely.
The Maxwell equations (3.2) reduce to the following equations for a µ , where the inhomogeneous terms are induced by the fields A µ (x α ) and ρ(x α ). Thus, the solutions a µ are linear functionals of these fields. Following a similar procedure in [26,27,29,30], a µ are decomposed in a basis made of A The decomposition coefficients S i and V i are functions of r and also functionals of the boundary derivative operators ∂ t and ∂ 2 i . In Fourier space, S i and V i become functions of (r, ω, q 2 ). The equations (4.7) for a µ map into a system of second order ordinary differential equations. We summarised these equations by grouping them into decoupled sub-sectors. (4.10) (4.11) The boundary conditions for S i and V i are derivable from those for a µ . The condition (4.5) translates into The regularity condition for a µ (r) is equivalent to that for the coefficients S i and V i . Finally, the Landau frame choice becomes a non-trivial constraint on the pre-asymptotic behaviour of the coefficients S i (see below).
Since the boundary current J µ (3.13) is defined in terms of the near-boundary behaviour of the field A µ , our next step is to re-express it in terms of the pre-asymptotics of the coefficients S i and V i . The coefficients S i and V i near the boundary r = ∞ have the form (4.13) Eqs. (4.9, 4.10, 4.11) were used in order to obtain (4.13). The boundary condition (4.12) was also imposed. From (3.13), the current J µ reads (4.14) The Landau frame convention (4.6) requires (see Appendix A where this condition is relaxed) The remaining coefficients (v i ) cannot be fixed from pre-asymptotic considerations only and have to be found through integration of the dynamical equations (4.9, 4.10, 4.11), from the horizon to the boundary. We, however, were able to find two relations among the coefficients S i and V i , which help to reduce the number of unknown coefficients. Consider the combinations which satisfy two coupled homogeneous ordinary differential equations (4.17) X and Y have to vanish under our boundary conditions. Thus, two relations emerge At r → ∞, the first relation is already imposed by the Landau frame choice (4.15). The second one is which will be used to eliminate v 1 in favour of the others. Finally, the current J i expressed in a gauge invariant form reads The transport coefficients D, σ e and σ m are So far, we have only expressed the transport coefficients in terms of the coefficients v i . The latter are yet to be determined by solving the ODEs (4.9, 4.10, 4.11). This will be our goal in the next section. We will first solve these equations analytically in the hydrodynamic limit ω 1, q 1 and then numerically for generic values of ω and q.

Boundary current: Results
In this section, we solve the equations (4.9,4.10,4.11) and determine the transport coefficients, first perturbatively and then numerically.

Perturbative analytic results: hydrodynamic expansion
In the hydrodynamic limit ω 1, q 1, the equations (4.9,4.10,4.11) could be solved perturbatively. To this goal, we introduce an expansion parameter λ via ω → λω, q → λ q. Then, S i and V i are formally expandable in series of λ, λ keeps track of the gradient expansion order and, eventually, will be set to unity. The expansion (5.1) is substituted into (4.9,4.10,4.11), and the equations are solved iteratively, order by order in λ. At each order in λ, direct integration results in the expansion coefficients S When substituted into (4.21), these are the results quoted in (1.11).

Numerical results: all-order derivative resummation
To all order in the gradient resummation, the equations (4.9,4.10,4.11) have to be solved for generic values of ω and q 2 , which we have been able to do numerically only. We are facing a boundary value problem for a system of second order ODEs. Much like what was done in [26,27,29,30], we apply the shooting technique. We briefly sketch the numerical procedure. First, take any trial initial values for functions V i and S i at the horizon r = 1 and integrate the equations up to the boundary r = ∞. The solution generated in this way has to satisfy the requirements (4.12) and (4.15) at the conformal boundary. If it does not, the trial functions must be adjusted and the procedure is repeated until the boundary condition is fulfilled with a reasonably good numerical accuracy. This fine-tuning process of finding the correct initial data is reduced to root-finding problem, which we solve by the Newton's method.

Diffusion
We start presenting our numerical results with 3D plots for the diffusion function D (Figure 1). In figures 2, we plot 2D sections for ω = 0 and q = 0, which are more transparently revealing the asymptotic behaviour of D. Imaginary part of D vanishes at ω = 0. As seen from the plots, D is a very smooth function of q 2 . The weak dependence on q 2 reflects quasi-locality of the diffusion process. The ω dependence is much more interesting. The imaginary part of D first approaches a maximum at ω 1.5, which could be interpreted as an emerging scale in the problem. Both the real and imaginary parts display damped oscillations as a function of ω. This behaviour indicates the presence of a set of complex poles, which are the quasi-normal modes of the bulk theory. Importantly, D eventually vanishes at large momenta. This behaviour is consistent with restoration of causality.
It is interesting to notice that the diffusion function D behaves much like the viscosity η. The latter was computed in [26,27]. For the sake of comparison, we quote the results here ( figure 3). Roughly speaking, η ∼ 2D. It is not exact, however, as is obvious from the  small momenta expansion for the viscosity η which is not proportional to one of D. In principle, there seems to be no reason why these two transport coefficient functions should be proportional. Yet, the overall similarity in the behaviour is quite striking suggesting a sort of universality in dissipative transport.
Vanishing of D at large ω signals recovery of causality but does not imply it automatically. To deeper explore the problem, we follow [30] and turn to the memory function formalism. The memory function D is defined in (1.12). It is displayed on figure 4. As anticipated, causal D indeed does not have support at negative times (up to numerical noise).

Conductivity
The results on the conductivities are shown in figures 5, 6 (3D), and 7 (2D). Similarly to the case of diffusion, the dependence on the spatial momentum q is much less pronounced compared to the frequency ω. The conductivity σ e is a growing function of ω. It will become clear from the next subsection that at q = 0, σ e coincides with the current-current two-point correlator computed analytically in [42,43],

Continuity equation and two-point correlators
In section 4, the dynamical components of (3.2) were solved for the correction a µ . The resulting solutions for a µ were expressed in terms of unspecified boundary fields A (0) ν and ρ. The fifth Maxwell equation is the constraint EQ r = 0 ν . In Fourier space, the current-current correlators G µν are The transverse component of the correlators G µν was analytically computed in [37] for q = ω. The analytically computed spectral function χ µ µ ≡ −2Im[η µν G µν ] is confronted against our numerical calculation in figure 8, clearly demonstrating high accuracy of the latter. The longitudinal response function has infinitely many poles, solutions of the dispersion  In the bulk, these poles correspond to the quasi-normal modes. They are absent in the transverse part of the correlator. In the hydrodynamic limit, substituting the expansion (1.11) and solving the dispersion relation (5.7) perturbatively, we reproduce the lowest diffusive pole which is in perfect agreement with [13]. Perturbative expressions for the correlators (5.6) are also in agreement with [13].
6 Remark on electromagnetic properties of the medium: dielectric tensor and refractive index So far we have treated the external fields as non-dynamical. As was mentioned earlier, we can gauge these fields turning them into fully dynamical. Coupled with the induced current, for which the constitutive relation was constructed, the theory at hand becomes a self-consistent electrodynamics of a conducting medium. Electromagnetic properties of this medium are related to two-point correlators and in principle could be deduced from the literature. For self-completeness of our study we summarise them below.
For medium with spatial dispersion, the more appropriate way to describe its response to external fields is the ( E, D, B) approach, where H is set equal to B. The macroscopic Maxwell equations are supplemented by a constitutive relation for the matter whereˆ is a generalised dielectric tensor, taking into account both the dielectric and magnetic response. In Fourier space,ˆ is a function of both ω and q. A comparative discussion between the ( E, D, B) approach and the formalism based on four fields-( E, D, B, H) can be found in the textbooks [53,54].
For an isotropic medium,ˆ is decomposed into longitudinal and transverse components 3) The constitutive relation (1.5) is the relation for the induced current Using the on-shell condition (continuity equation) and the Faraday law 5 (the third Maxwell equation), the constitutive relation (6.4) becomes the Ohm's law for the induced current The longitudinal and transverse components of the conductivity are proportional to the current-current correlation function G µν and are trivially related to the transport coefficient functions of (6.4): The dielectric functions are related to the conductivity via the relation: Numerical results for L,T are shown in figure 9. The electromagnetic waves propagating in the medium should satisfy dispersion relations Transvese mode : q 2 ω 2 = T , Longitudinal mode : L = 0. (6.8) The transverse mode defines the refractive index n 2 = T , which has two branches for n.
Depending on the sign of Im[n], only one branch is physical. The first three propagating waves for both transverse and longitudinal modes were computed numerically in [55], revealing a rich structure of the electromagnetic response of such holographic medium. Remarkably, the authors of [55] found that one of the transverse modes displays negative refractive phenomenon. Previous studies [56,57] noticed this phenomenon in the hydrodynamic regime of charged fluids which holographically are described by R-charge black holes. The criterion proposed in ref. [58] for negative refraction was used in [56,57]. The latter criterion is, however, valid for perturbatively small q only. The negative refraction mode of [55] is gapped and cannot be treated perturbatively. This is the reason the authors of [55] relied on the condition Re[n] < 0 (Im[n] > 0) to establish the phenomenon.
We would like to remark that requiring negativity of Re[n] might not be sufficient to make a conclusive statement about the existence of propagating modes with negative refraction. One must carefully study the energy propagation in the medium and establish direction of the Poynting vector. The time-averaged Poynting vector is [53] where U is the electromagnetic wave energy density. While refs. [55][56][57] have made a significant progress towards establishing the negative refraction in holographic fluids, we feel that a further exploration of this phenomenon is worth pursuing.

Summary
The most general off-shell constitutive relation for a globally conserved U (1) current (1.5) can be parametrised using three independent momenta-dependent transport coefficient functions: the diffusion D(ω, q 2 ) and two conductivities σ e,m (ω, q 2 ). The latter couple the current to external electromagnetic fields. At small momenta, these functions are power expandable reflecting the gradient expansion of the hydrodynamic regime. In contrast, causality of the theory is manifested in the large momentum asymptotics. Precise determination of all three coefficient functions using the fluid-gravity correspondence is our main result. In the hydrodynamic regime, we analytically reproduced all known results on the gradient expansion and extended them to the third order (1.11). Particularly, for the first time an analytic expression was given for σ m . The rest of the results were obtained numerically, beyond hydrodynamic regime. Of particular interest is the large ω behaviour, which for the first time was shown to be consistent with the causality constraints. We also studied causality using a numerically computed memory function. The memory function D(t) was demonstrated to have no support for negative times thus further revealing causality restoration due to all order derivatives resummation.
We have noticed an interesting similarity between the diffusion D(ω, q 2 ) and the viscosity η(ω, q 2 ): both functions demonstrate a strikingly similar behaviour as functions of ω and q (even though they are not exactly equal). We speculate that this might point to a certain universal behaviour among various dissipative transport coefficient functions.
The transport coefficients were deduced using the holographic prescription. In the dual gravity side, we constructed a general solution for a probe Maxwell field in Schwarzschild-AdS 5 black hole spacetime. The Maxwell equations in the bulk were split into dynamical and constraint components. Interestingly, the off-shell constitutive relation (1.5) and the corresponding transport coefficient functions were found to be fixed uniquely through solutions of the dynamical components only. The constraint component was shown to be equivalent to the continuity equation ∂ µ J µ = 0.
As yet another check of our formalism, we used the continuity equation (on-shell condition) to reproduce some known results on two-point correlators (5.6). It is, however, worth emphasising that in linear response there are only two independent two-point correlators (transverse and longitudinal), whereas our formalism reveals three fully independent transport coefficients. Our off-shell structure emerges as an important element of hydrodynamic effective action [31] and is clearly more rich than what is encoded in the on-shell correlators.

A Relaxing the Landau frame condition
In this Appendix, we demonstrate that all the transport coefficients can be uniquely fixed without imposing the Landau frame convention (4.6). The latter is equivalent to fixing the residual gauge, whereas the physical quantities (transport coefficients) are gauge invariant.
The gauge invariance of J µ (eq. (4.14)) under boundary gauge transformation gives the relations Along with the condition (4.5) and regularity requirement for S i , V i at the horizon, the relations (A.2) are equivalent to the constraints (4.18). The constitutive relation for J i of (4.14) in a gauge invariant form is Here ρ was exchanged in favour of J t and the Bianchi identity ∇ × E = −∂ t B was used to substitute E i for B i . The transport coefficients of interest are The coefficients σ e and σ m are linear combinations ofσ T andσ L . The equation for V 1 is decoupled from all the other equations. Moreover, the asymptotic boundary condition (4.5) and regularity at horizon are sufficient to fully determine V 1 . Thus its pre-asymptotic value v 1 and hence the transport coefficientσ T are uniquely fixed and just the ones quoted in the text.
The situation with the coefficients D andσ L is more involved. Below we will prove that D andσ L are also uniquely fixed even without the Landau frame condition (4.15) imposed.

(B.18)
The boundary data v i are straightforwardly read off from the above perturbative solutions.