Hydrodynamics of a relativistic charged fluid in the presence of a periodically modulated chemical potential

We study charged hydrodynamics in a periodic lattice background. Fluctuations are Bloch waves rather than single momentum Fourier modes. At boundaries of the unit cell where hydrodynamic fluctuations are formally degenerate with their Umklapped copy, level repulsion occurs. Novel mode mixings between charge, sound, and their Umklapped copies appear at finite chemical potential -- both at zero and finite momentum. We provide explicit examples for an ionic lattice, i.e. a periodic external chemical potential, and verify our results with numerical computations in fluid-gravity duality.


Introduction
In considering the quantum mechanical wave function of a single electron in a lattice of atoms Bloch had the insight that one should expand the wavefunction in a manner consistent with the discrete periodicity1 The novel part of Bloch was its application to quantum wavefunctions rather than waves in general.How waves propagate in periodic structures was already considered by Newton, and that waves in periodic structures exhibit peculiar interference phenomena that we now know as level repulsion/Umklapp/gap opening at Brillouin zone boundaries or Bragg reflection from point-like lattices was already recognized by Kelvin in the 1880s [1].In electrical engineering the propagation of electromagnetic waves in periodic structures was [2], and is an important topic, see e.g.[3]. 2 Also sound waves in lattices were considered from the earliest days up to today, see e.g.[4].Sound waves, however, are hydrodynamic fluctuations -a long-time long-wavelength perturbation around thermodynamic equilibrium of a conserved charge associated to a global symmetry -and in that sense differ from electromagnetic waves or single particle wavefunctions in that the fundamental equations of motion, i.e. the hydrodynamic conservation laws, are non-linear.The wave-like fluctuations propagate on a background that is itself a full (equilibrium) solution to the non-linear set of equations, and through the non-linearity the properties of the fluctuating waves depend on this background solution.Though gradients are energetically disfavored, through external forcing the equilibrium background can be imprinted with a spatially varying temperature T (x), pressure P(x), or chemical potential µ(x).Due to the non-linear coupling between fluctuations and the background in hydrodynamics, the wave propagation properties can be self-consistently determined from the (spatially varying) background.This was elucidated particularly clearly in recent years in the context of electron hydrodynamics in systems with random charge impurities [5,6].Such charge disorder is encoded in a spatially varying chemical potential with average [µ(x)] = µ 0 and variance [µ(x)µ( y)] − [µ(x)] [µ( y)] = σ 2 µ δ(x − y).Quantum mechanical single particle electron motion in the presence of random impurities is a classic condensed matter problem.As Anderson showed, the random wavefunction interference is essentially uniformly destructive; at low temperatures all motion is inhibited and the system becomes an insulator.In the hydrodynamic regime, however, i.e., in a situation where many electrons collectivize to a classical fluid rather than a quantum mechanical wave, the conductivity rather strikingly remains finite indicating the existence of an "incoherent metal" state [6].Observing this electron hydrodynamics in sufficiently pure 2D systems is currently actively pursued, see e.g.[7] or [8], references therein and the recent review [9].
Here we study not hydrodynamics with random spatial disorder but with strictly periodic modulations of the background, i.e. a lattice.Moreover, we also consider hydrodynamics of a charged rather than a neutral fluid with an eye towards condensed matter systems.Compared to the many existing studies on sound waves in periodic structures, the presence of electromagnetic charge as an additional conserved quantum number changes the fluctuating wave response fundamentally.This is again due to the non-linear nature of the hydrodynamic equations.At finite chemical potential sound mixes with charge diffusion.In a companion article we focus on the significant consequences of this cross-coupling of Bloch modes in a lattice for the measurable DC and AC conductivities in condensed matter systems where this hydrodynamics approach may apply [10].In this article we provide the deeper hydrodynamic analysis of the full fluctuation spectrum of charged hydrodynamics in a periodic background.

Hydrodynamics: Set-up and brief review of homogeneous fluctuations
The principal reason that linearized hydrodynamic fluctuations in a lattice background should also be expanded in Bloch modes has already been emphasized: the essence is wave propagation in a periodic structure.Waves are described by coupled first order differential equations of the form 3 (∂ t + M (x))φ(x) = 0 . (3) If M (x) is periodic M (x + 2π G ) = M (x), then φ(x) can be decomposed in Bloch waves 4 φ(x) = n G/2 −G/2 d k φ n (k)e i(k+nG)x .Taking M (x) = −M 0 ∂ 2 x +Acos(G x) as canonical example, one can solve Eq. (3) perturbatively in A. Diagonalizing M in terms of φ n (k)= p A p φ (p) n (k)/p!, the lowest eigenvector to first order in A is n+1 (k) + . . ., (5) in terms of the unperturbed eigenmodes.This mixing between the different Bloch waves is Umklapp.In this article we shall only focus on these perturbative solutions for small lattice amplitudes. 3The standard wave equation (∂ 2 t − M 12 M 21 )φ 1 = 0 follows from (2) 4 The Bloch theorem essentially states that the plane wave decomposition φ(x) = ∞ −∞ dqφ(q)e iqx can be segmented into unit cells q n ∈ [(n − 1/2)G, (n + 1/2)G] where n ∈ labels each cell -or Brillouin zone - dq n φ(q n )e iq n x .The wavevector in each Brillouin zone can be shifted The advantage of this decomposition is that discrete periodic shifts x → x + 2mπ/G relate modes in different Brillouin zones at the same Bloch momentum k ∈ {−G/2, G/2}.
We also already noted that what is special about hydrodynamics is that the fluctuation equations are themselves a linearization expansion of the fundamental non-linear equations.The principle behind the theory of hydrodynamics is local equilibrium and encoded in the local conservation laws of macroscopic charges, i.e., of a slowly spatially varying energy-momentum tensor T µν (x) and in the presence of a U(1) charge, a current J µ (x).In turn this implies that one can also describe fluid behavior in the presence of a slowly spatially varying external potential whether temperature T (x), pressure P(x), or chemical potential µ(x).
For simplicity -as well as for the experimental supposition that strongly correlated condensed matter systems can have an emergent Lorentz symmetry at low energies -we shall use d = 2 relativistic charged hydrodynamics in this article.In principle all we state also applies to arbitrary d non-relativistic charged hydrodynamics, even if the precise expressions may be subtly different.In relativistic charged hydrodynamics the dynamical equations are simply the conservation equation of the energy-momentum tensor and the charge-current Here we have allowed for an external electromagnetic field strength in terms of a local external vector potential.In this paper, we will be interested in taking A µ,ext = (µ ext (x), 0, 0) with µ ext (x) a periodic function.Though again, in principle our results also hold for a spatially varying (external) pressure (see e.g.[11]), or a spatially varying (external) temperature. 5he dynamical variables of the fluid are the temperature T , the unit timelike velocity vector and the chemical potential µ.Away from equilibrium, the conserved currents in our theory -which we assumed to be parity-invariant, see [12] for more general cases -are given by the constitutive relations at first order in gradients in Landau frame Here d = 2 is the number of spatial dimensions and the shear viscosity η, the bulk viscosity ζ, and the microscopic conductivity σ Q are hydrodynamic transport coefficients -in principle set by the microscopic details of a given theory, see e.g.[13,14], in practice phenomenologically determined.∆ µν = η µν + u µ u ν is a projector orthogonal to the fluid velocity.The Landau frame choice is such that at any order in gradients, we have J t = n and T t t = ε.The above constitutive relations also hold in a static equilibrium background.In a system with Galilean or relativistic Lorentz boost invariance -which we use in this paper -it is convenient to choose the reference frame for which the equilibrium fluid is at rest.In absence of contact to a spatially varying heat bath, the temperature must then also be constant and independent of position.In the presence of a spatially varying external chemical potential µ ext (x), the equilibrium solution to the hydrodynamic equations Eqs. ( 6) is then parametrized as In the grand canonical ensemble, the hydrostatic equilibrium yields moreover Throughout this paper we will use the bar notation X to denote such static background quantities.For a homogeneous background, they will be spatially constant and we will use a subscript 0 as X 0 to denote them.For spatially varying quantities, we will use superscripts Y (n) to describe higher order (Bloch wave) moments The hydrodynamic equations need to be supplemented by an equilibrium equation of state relating the energy density ε, the pressure P and the charge density n to solve in terms of the equilibrium values of T and µ.In this paper, we will be considering a general fluid whose equation of state P(T, µ) determines the thermodynamic equilibrium of the theory.In Sec. 4, we will specialize to conformal systems.
On top of this background, we now consider perturbations X (t, x) = X (x) + δX (t, x).The conservation equations ( 6) then take the form which can further be simplified into In the previous expression, we have defined a renormalized viscosity η ≡ ζ + We also introduced the "potential"-variations δλ ε ≡ δT T 0 and δλ n ≡ δµ − μ T 0 δT conjugate to the energy and charge densities.The velocity perturbations δv i are conjugate to the momenta δπ i .These are not independent due to the hydrodynamic local equilibrium condition.The momenta δπ i are related to the velocity perturbations δv i through the constitutive relations δπ i = δT t i = χππ δv i at this order. 6Similarly, the charge and energy densities δn, δε are related to the sources δλ n , δλ ε through the thermodynamic susceptibilities derived in Appendix A. Using that χππ = ε + P, we can use the fundamental thermodynamic relation ε + P = T 0 s + μn and the first law of thermodynamics to relate Finally, we introduced an external electric field δE x ≡ ∂ t δA x,ext and in (12) we introduced the total electric field δE tot x ≡ δE x +∇ x δλ n .In what follows, we will be interested in the hydrodynamics response of the modes {δε, δn, δπ x , δπ y } obeying the equations (12).In the following sections, we will use the static susceptibilities relating the potentials {δλ ε , δλ n , δv x , δv y } to the densities {δε, δn, δπ x , δπ y } to express the equations in terms of the latter only, i.e., we work in the microcanonical ensemble. 6Formally the susceptibility χππ (x 1 , t 2 ; ) is influenced by a (chemical) potential v i (x 2 , t 2 ) at a different space-time point.Here Z(v i ) is the partition function in the presence of a chemical potential (velocity) v i for the charge density (momentum) π i .In the hydrodynamic limit, however, one assumes that all equilibrium (t 1 + t 2 = 0) static (t 1 − t 2 → ∞) charges depend only locally on the potentials π i (v i (x)).In a homogeneous equilibrium background where χstatic ππ = χ static (x 1 − x 2 ) this is equivalent to approximating the static susceptibility with its constant part χstatic ππ (x 1 − x 2 ) = χ ππ,0 + (x 1 − x 2 )∂ x χππ (0) = χ ππ,0 + . ... We discuss this in more detail below and in Appendix A.

Hydrodynamic fluctuations in a homogeneous background
In this section, we first review the hydrodynamics of a long wavelength perturbation above a homogeneous conformal charged fluid.Further details can be found in [15,16].Since the background is homogeneous, this means every barred quantity will be a constant X = X 0 .The equations (12) decouple into the longitudinal and transverse sectors.We will start by looking at the latter whose equation of motion is simpler.Choosing the wavenumber k x along the x direction without loss of generality, the transverse fluctuations δπ y obey This is a simple diffusion equation with the shear diffusion constant D ⊥ .We can now use the (Fourier-)Laplace transform 7 such that the transverse equation of motion becomes The solution is formally given in terms of the retarded correlator for the transverse momentum which is defined as Using that ĜR π y π y (ω = 0, k) = χ ππ,0 in the hydrodynamic long wavelength limit (i.e., we only keep the leading term in an expansion in k; see footnote 6.), we have ĜR The correlator exhibits a pole on the imaginary axis at ω = −i D ⊥ k 2 indicative of a purely diffusive mode.
We can now carry the same analysis in the longitudinal sector where the dynamical equations are coupled.Denoting δφ a (t, k) = (δε(t, k), δn(t, k), δπ x (t, k)), the dynamical equations can be written succinctly as We can once again use a (Fourier)-Laplace transform to rewrite this system of equations as with the dynamical matrix where we defined d χ = χ εε,0 χ nn,0 − (χ nε,0 ) 2 the determinant of the susceptibility matrix in the ε, n sector.The poles of the Green's functions associated to this system are the frequencies for which det K = 0.The roots of this polynomial in the long wavelength limit are a diffusion mode (originating in charge diffusion) and two propagating sound modes where the speed of sound and the two diffusion constants are defined as Note that for n 0 = 0 the speed of sound reduces to the familiar expression c 2 s = δP/δε (using Eq. ( 13) and the inverse susceptibility matrix).Similarly to what we did in the transverse sector, we can compute the retarded Green's functions by inverting the dynamical system [16] ĜR where the middle equation enforces the condition that the static ω = 0 part reduces to the longitudinal part of the thermodynamic susceptibility matrix χ L,0 .The various correlators can then be obtained ĜR with the normalized determinant of the dynamical matrix The other correlators can be obtained via the Ward identities ĜR as well as the Onsager reciprocal relations

Hydrodynamic fluctuations in a lattice background
We shall now redo the fluctuation analysis in a lattice background.This lattice will be sourced by a periodically modulated external chemical potential such that the fluid is still at rest and in local equilibrium, but all its constituents will now slowly vary in space.In particular, this last assumption of local equilibrium means that the scale of spatial fluctuations of µ ext and other local quantities must be larger than the local equilibration scale.Therefore, we must have G ≲ µ, T . 8his lattice background manifestly breaks translation invariance.Momentum is therefore no longer a strictly conserved quantity.However, as the breaking is sourced through a hydrodynamic variable and as we assume it is weakly broken, we can still use hydrodynamic analysis [5,6,17,18].The spectral function of the associated operator to this deformationthe charge density J t = n -, evaluated in the homogeneous background, can be used to compute the momentum relaxation rate.This is known as the memory matrix formalism and was thoroughly detailed in e.g.[5,19].The momentum relaxation rate induced by an operator O sourced at wavenumber G with strength g takes the form [20] For a cosine ionic lattice Eq. ( 27), g = µ 0 A/2, and we have two deformation sources, one copy each at ±G -noting that the expression ( 28) is parity invariant in G. Therefore, the memory matrix relaxation rate for an ionic lattice is It will prove useful to separate the terms according to their scaling with G in this expression as Γ ionic,mem.= Γ η + Γ d with Using the Einstein relations Eq. (22a), together with χ ππ,0 = ε 0 + P 0 , χ πn,0 = n 0 these are a convective shear drag term Γ η and an intrinsic diffusive term [6,20,21] We will recover this same expression for the momentum relaxation time from our Bloch wave analysis.This analysis improves on the memory matrix technique by understanding how all the hydrodynamic fluctuations behave.In a periodically modulated background, every background quantity in local thermal equilibrium X (x) = X ( μ(x), T 0 ) now admits Fourier series expansions In order to apply the same method as in the previous section, we must first know how to relate perturbations of sources and responses in this new background.Because the background is static, the susceptibilities will also be static.However, because the thermodynamic quantities are position dependent and have non-vanishing Bloch modes, the susceptibilities will now also be position/momentum dependent.In principle, they depend on two Bloch momenta.However, in the slowly varying hydrodynamic background we may approximate them as local functions χ(x) (see also footnote 6) that follow the expansion (32). 9The relation between 9 One can analyze the general behavior of two-point functions under lattice symmetries of the background [22].

Given a two-point function G(x, y), one can pick a center of mass point
k+2nπ/L)r .For hydrodynamic susceptibilities χ = G J t J t we assume that they are local, i.e., we can restrict to δ = 0 to leading order.In a strictly periodic background there is no structure beyond the lattice scale and hence only the χ (n) (k = 0, δ = 0) modes are non-vanishing.perturbations in the sources and responses is then The breaking of isometry by the lattice means there is no longer a decoupling between a longitudinal and transverse sector, i.e., φ A , λ A collectively denote the responses {δε, δn, δπ x , δπ y } and the sources {δλ ε , δλ n , δv x , δv y }.Both perturbations are likewise expanded on Bloch modes matching the discrete lattice symmetry for X ∈ {δφ A , δλ A }.As a result of the spatial dependence in the background different Bloch modes of the perturbations cross couple The dynamical equation can then be written, after Laplace transform, as The indices n, m indicate the Brillouin zones while each block K (n,m) is a 4 × 4 matrix.The diagonal blocks K(n,n) correspond to the couplings between the responses in the same Brillouin zone while the off-diagonal blocks will account for coupling between different zones.These are due to the presence of the lattice and will vanish in the limit where the lattice amplitude goes to zero A → 0. We will be interested in a weak lattice where the lattice amplitude A is very small, and keep only terms up to order A 2 . 10The coupling between two modes with momenta k + nG and k + mG for m > n will be of order A m−n .Moreover, within perturbation theory, terms of order A in the off-diagonal blocks will contribute to the same order as terms of order A 2 in the diagonal blocks; we can therefore drop terms of order A 2 and higher in the off-diagonal blocks.This also means we can consider "nearest-neighbor" interactions onlyby which we mean off-diagonal terms with m = n ± 1.In the long wavelength approximation we therefore can narrow our study to the three momenta k + nG with n ∈ {−1, 0, 1}, i.e., the first three Brillouin zones.It is important to note that the diagonal terms even in the n = 0 Brillouin zone can still have non-trivial higher order corrections in A. A similar setup was already considered in [23].
We will discuss this momentarily.We shall, however, first make one more simplification.It will prove more useful to use the equations in terms of the sources δλ A with where we can relate the two matrices using the susceptibility matrix χ by K = K • χ.In this language, the A = 0 dynamical matrix (20) takes the form This choice seems to a priori obfuscate the relationship between modes more than (20) due to the off-diagonal frequency dependency.However, because χ is a static matrix, the determinants of K and K have the same poles in the complex frequency plane, and in the lattice case where the inverse susceptibilities present in (20) are more complicated, this form will prove clearer.
In the next few sections, we will determine this matrix K in a lattice background with lattice vector G for both finite k momentum fluctuations and k = 0 momentum fluctuations to order A 2 in the lattice amplitude.As standard, the zeroes of its determinants will indicate the position of the dynamical modes of this system.We will then compute the conductivity as an example of how the various correlators are modified by the presence of the lattice.

Finite momentum aligned fluctuation spectrum
For a generic fluctuation with momentum k, even in the long wavelength limit, the fluctuation matrix truncated to nearest neighbor cross-coupling sufficient for the leading order in A correction will be a 12 × 12 matrix.This is because there is no decoupling into transverse and longitudinal sectors for a generic momentum.However, if one chooses the fluctuation momentum k to align with the lattice wavevector, a decoupling does occur.Choosing k along a lattice vector defined to be in the x-direction, a parity symmetry in the y-direction remains.The even and odd sectors decouple into the longitudinal and transverse parts: x , δλ (1)  ε , δλ (1)  n , δv (1)   x , (39a) The dynamical matrix is then diagonal in a 9 × 9 and a 3 × 3 block.

Transverse sector
Starting with the transverse sector, the associated dynamical matrix KT is In the hydrodynamic approximation the local static susceptibility χππ (x) = ∂ π x ∂ v x (µ(x)) (see footnote 6 & 9) now also depends on the lattice amplitude as can be seen from its Bloch components In the previous expression, we have introduced the expansion for a given Bloch mode . Note that by definition, X (0) 0 = X 0 which we will keep this way.The poles of the transverse fluctuation matrix can now be found easily, and we have ππ,1 χ ππ,1 χ ππ,1 χ The poles remain purely diffusive, and we see that the only effect of the lattice on the transverse sector is to renormalize the shear diffusion constants D ⊥ at order O(A 2 ).We do see an Umklapp-like pole in the dispersion relation at the edges of the Brillouin zones k = ± G 2 .Formally, this value of k is outside of the regime of validity of the expansion in small A. One has to resum the perturbative expansion and then one finds level repulsion, as is well known; see also the discussion at the beginning of Sec.3.2 and footnote 11.It is distinct from conventional Umklapp, however, in that it is not level-repulsion in the dispersion (the real part of the pole in the complex frequency plane), but in the width of the fluctuation.At the edge of the Brillouin zone the width narrows and vanishes at exactly k = ± G 2 .

Longitudinal sector
The longitudinal sector is characterized by a 9 × 9 dynamical matrix KL of the form of 3 × 3 blocks The K (OD) L (ω, k, p) block with k < p belongs to the Bloch sector n = −1, and the one with k > p to the Bloch sector n = 1.In a cosine lattice, however, all background quantities are parity-invariant X (−n) = X (n) , and so from here on out, we will only use the n > 0 expressions.The 3 × 3 blocks To leading order in A 2 , the n = 0 Bloch momenta X (0) still have a dependency in A just as in the previous section.
While difficult, it is possible to find the poles associated to this 9 × 9 matrix generically.For very small momentum k = O(ϵ 2 ) and G = O(ϵ), they take the form ) with n ∈ {−1, 1} and ". .." indicate corrections of order O(A 2 k) and higher.The relaxation rates Γ d , Γ ionic,mem.are of order A 2 and equal to the memory matrix expressions given in Eqs.(31).
For large k the expressions are not easy to express.However, we can use the mixing with Umklapped Bloch waves analysis to understand numerical simulations.In the longitudinal sound sector we do observe genuine level repulsion in the modified dispersion relation at the edges of the Brillouin zone.For a visualization in an explicit example later, see Fig. 8.There is thus a true sound "band gap".Sound modes with frequencies ω = ±c s G/2 do not exist in this latticized medium.Or more precisely put, sound with wavelengths λ = 2πk ≪ G propagate normally with essentially unaltered speed of sound c s = dω dk = c 0 s .As the wavelength approaches the edges of the Brillouin zone, sound slows down, and right at the edge of Brillouin zone for λ = (2π) G 2 , they cease to propagate as the group velocity c s = dω dk | k=G/2 = 0.The medium is opaque to sound at these wavelengths.Considering possible applications to condensed matter physics, we note for completeness that all these results are of course derived assuming a fixed infinitely stiff external lattice.Lattice vibrations/phonons are not taken into account.Were one to include these in the analysis, this will likely make the level repulsion and opaqueness to sound less sharp.

The k = 0 zero momentum perturbation
The k = 0 zero momentum is special and asks for a separate discussion.This is for three reasons.Again in the context of condensed matter physics, the k = 0 fluctuation describes the homogeneous responses of the system to outside probes.These are the observed macroscopic thermal and electrical conductivities, and warrant being singled out.Secondly, we shall see that in the limit of k → 0 several modes becomes degenerate.One must always be careful with accidental degeneracies.This is also the case here.The degeneracy is lifted in the presence of the lattice deformation.However, since we only consider the lattice perturbatively, this implicitly means we consider AV int ≪ k where V int is a characteristic scale denoting the strength of the interactions between the Bloch modes.The degeneracy limit and the small lattice amplitude limit do not commute.We shall illustrate this in more detail below.We can still do a perturbation analysis in A, but this must be done from the k = 0 starting point separately. 11 Finally, mathematically, the k = 0 fluctuation is special in that at vanishing momentum, parity in the x-direction (G ↔ −G) is restored.In the 1D lattice we consider -with lattice vector in the x-direction -the longitudinal and transverse fluctuations at k = 0 therefore break up 11 A simple example that illustrates the point is the toy model fluctuation matrix This has poles at E = ± k 2 + A 2 V 2 int signaling level repulsion at k = 0. Expanding these poles in A gives into odd and even superselection sectors under G ↔ −G δY L− = δλ (1)  ε − δλ (−1) ε 2i , δλ (1)  n − δλ (−1) n 2i , δv (1)  x + δv (−1) x 2 , δv (0) x , (47a) δλ (1)  ε + δλ (−1) ε 2 , δλ (1)  n + δλ (−1) n 2 , δv (1)  x − δv (−1) δv (1)  y + δv (−1) For the sake of brevity, as k = 0 we have suppressed all k arguments in the dynamical expressions δ In this basis, the overall dynamical matrix K′ = U KU −1 is diagonal by block and the dynamical equations take the form where U = U L 0 0 U T is the matrix that reorders the fields from the basis in Eq. ( 39) to Eq. ( 47)

Transverse sector
Let us again consider the transverse sector first.The dynamical matrices KT− and KT+ are These have the following diffusive poles ) There are several aspects to note: firstly, as mentioned above these poles do not correspond to the k → 0 limit of the finite k fluctuations in Eq. ( 42).The indicated emergent degeneracy at k → 0 between ω T,−1 and ω T,1 is obvious in Eq. ( 42).The lattice perturbation A lifts this degeneracy and therefore the limits k → 0 and A → 0 do not commute.Secondly, there is a zero mode in the T +-sector.This is the standard transverse k = 0 excitation, that corresponds to a change of the static homogeneous transverse pressure background and as a zero mode should not be considered in the fluctuation spectrum.

Longitudinal sector
Consider the G-parity odd longitudinal sector first.Its dynamical matrix KL− is given by with The lines are there to highlight the coupling between two sub-sectors.The top-left block is equivalent to the coupling matrix (20) in the homogeneous system at momentum G encoding the G-parity odd part of two sound modes and a charge diffusion mode (see Eq. ( 21)), while the lower-right block reflects the conservation of momentum in the homogeneous case.In the presence of the lattice deformation the momentum conservation mode now couples with the Umklapped finite G sound-, and charge diffusion modes through the off-diagonal terms.We can find the modes of this dynamical matrix in the same way we did before, and we find ω We therefore see that the poles in (54) are those of the homogeneous system (21) at momentum G with corrected diffusion constants due to the effects of the lattice.The procedure to compute the explicit form of the corrections s is detailed in Appendix C. For a generic relativistic fluid these are quite involved; in the explicit case of a fluid with conformal invariance they simplify greatly and we give the expressions below in Eq. (61).
The most noteworthy part is the Drude pole ω Drude .As previewed at the beginning of this Sec.3 the lattice breaks translational symmetry and this induces a momentum decay rate.The more detailed Bloch wave analysis which gives us all fluctuations at finite k, ω in the hydrodynamic regime beautifully recovers the finite ω memory matrix result (29), as it should.In the condensed matter context, it is this sector specifically that governs the k = 0 thermoelectric conductivities, where the presence of the second diffusive mode ω D (originating in Umklapped charge diffusion) in addition to the Drude mode has significant observable consequences as expounded in [10].
For the G-parity even sector the dynamical matrix is given by Once again, the lines show how the various sub-sectors couple through the off-diagonal Adependent terms.The L + sector is again a pair of sound modes and a charge diffusion mode, but now the part that is even under G-parity.The n = 0 charge and energy conservation modes are reflected in the upper-left blocks.The modes of this system are We see again that the latter three poles in (57) are those of the homogeneous system (21) at momentum G with corrected diffusion constants due to the effects of the lattice.The shifts differ, however, from the L− sector.The explicit form of the corrections D (L+),2 ρ , c (L+),2 s , D (L+),2 s can again be derived through the method detailed in Appendix C, and tractable expressions for the special case of a conformal fluid are given in Eq. ( 61).The first two poles are the ones encoding charge conservation and energy conservation; they remain unshifted at this order in perturbation theory.As we explain in Appendix C, one of these modes is an exact conservation mode and by rotating the system back to K instead of K, it is easy to see that this mode corresponds to overall charge conservation -associated here to δλ (0) n .The other mode is merely unshifted at this order in perturbation theory.
Having computed the corrections to the k = 0 modes, one clearly sees that these are not equal to the limit k → 0 of (45).In that limit, at leading order, ω D,±1 = −i D 0 ρ G 2 + i 2 Γ d and ω S,±,0 = − i 2 Γ ionic,mem., whereas the explicit k = 0 computation has ω (L−) Drude = −Γ ionic,mem., and an additional (to leading order in A 2 ) zero mode.This difference is due to the noncommutativity of the k → 0 and A → 0 limits.In the next section, where we illustrate the emergence of these hydrodynamic modes in an explicitly computed example, we will show precisely how these poles are related in the k → 0 limit.

Bloch wave hydrodynamics emerging from holographic models: A comparison
We will now validate the understanding of charged (relativistic) hydrodynamics in a periodic potential by comparing it with the low energy physics of holographic models.Holographic models describe the strong coupling regime of quantum field theories in a manifestly real time formalism.Uniquely so, this includes the emergence of hydrodynamics at low frequencies and long wavelengths ω, k ≪ T, µ [24][25][26][27].This last part is also known as fluid-gravity duality [28].By considering a strongly coupled quantum field theory in a spatially periodic chemical potential background, i.e., an ionic lattice, 12 described holographically in terms of its dual gravitational description, we will see that the Bloch wave hydrodynamics described above emerges.There is one simplifying feature in the two holographic models we choose here.Both describe a conformally invariant system for which the equation of state takes the scaling form which directly implies P = ε/d.Furthermore, due to the conformal symmetry such a system must also have a vanishing bulk viscosity, i.e., ζ = 0.
The two specific models we consider are the strongly coupled theories holographically dual to the Reissner-Nordström (RN) black hole [29][30][31][32][33] as well as the Gubser-Rocha (GR) black hole [34].We will solve the fluctuations in these systems modulated by a finite chemical potential numerically and compare to the predictions from Bloch wave hydrodynamics as presented in the previous sections.These two systems are chosen as their ground states are possible candidates to explain the mysterious strange metal physics underlying high T c superconductors.There is reason to believe that this physics is indeed that of strongly coupled electrons in an ionic lattice.The possible relevance of Bloch wave hydrodynamics in the context of strange metal physics is described in a companion article [10].
A brief description of the numerical holographic set-up is provided in appendix D.2; more details can be found in [10].The crucial aspect of relevance here to the comparison of the numerics with our Bloch hydrodynamic analysis is the equation of state of the two models.The 2 + 1 dimensional (finite temperature) field theories dual to AdS 4 RN and GR black holes are conformal charged fluids 13 with equation of state where we have defined T ≡ T /µ.A direct consequence of this conformal equation of state is that χ εε,0 = d χ ππ,0 and χ nε,0 = d n 0 .As a consequence many previous expressions simplify.Specifically the order A 2 corrections to the poles for k = 0 are now given by the tractable expressions: The explicit expressions for the thermodynamic quantities in the grand canonical ensemble for the RN and GR black holes can be found in Appendix D.1.We will use the longitudinal optical conductivity )〉 as a probe.Generically this current will receive contributions from all hydrodynamic fluctuations; these essentially determine the low frequency long-wavelength response.At finite k, this means we should see all 9 modes described in (45).At k = 0, however, the current is part of the L− sector and we will only see the first sector with its 4 modes.Fig. 1 gives a schematic sketch of what the spectrum of the current-current correlator -and therefore the optical conductivity -should look like in the complex frequency plane based on our hydrodynamic predictions.
Precisely this expectation is reproduced by the numerical results in holographic duals to RN and GR black holes where hydrodynamics is emergent.Fig. 2 plots the density of the absolute value and argument of the optical conductivity for small values of the real part of the frequencies, i.e., zoomed in near the imaginary frequency axis.Each picture is at a different value of k/µ ∈ {0, 0.001, 0.005}.We see that for k = 0 there are only two purely diffusive poles, as predicted, that split into two propagating and two diffusive poles at finite k.For finite k, there should also be a third diffusive pole very close to the real axis.Its weight is very low, however, but it can be unveiled by zooming in carefully.This was plotted in  A priori, a fifth pole ω D,0 at k > 0 also couples to the electrical current, but it is not visible on the range plotted.A more refined computation, does reveal it (Fig. 3).Fig. 3 for k/µ ∈ {0.006, 0.008, 0.01} (this choice of momenta proved more convenient to display).In all cases the location of these poles can be compared with the predictions from our hydrodynamical analysis after substituting in the relevant equation of state.The match is perfect for both k = 0 and k finite as denoted by the white circles in Fig. 2 and the triangles in Fig. 3. Similarly, we plotted in Fig. 4 the argument of the optical conductivity near the sound at momenta G ± k.The weight of these poles is very small and they are therefore difficult to identify in |σ|.The argument of σ, on the other hand, displays a jump at the poles.A similar analysis holds for the conjugate pair of poles at −G ± k.These sound poles give rise to a characteristic peak in the real conductivity, first noted in [36][37][38].
To illustrate in more detail the hydrodynamical origin of all these poles and their full explanation in terms of thermodynamic quantities, one can track the location of the poles as a function of temperature.Focusing only on the purely diffusive poles (two for k = 0, and three for k ̸ = 0), as they are more easily extracted numerically by scanning carefully over the negative imaginary frequency axis, we also find here a perfect match between numerics and hydrodynamic prediction, Eqs. ( 54) and ( 45) respectively, but now as function of T /µ; see Fig. 5 and Fig. 6.
With our computational RN and GR examples we can also illustrate the subtle nature of the k → 0 limit.As we saw in the previous section, the naive extrapolation to k → 0 of the two n = 0 sound modes ω S,0,± = − i 2 Γ ionic,mem.+ O(k) does not correspond with the physical = 0. To emphasize this once more, the origin of this difference comes from the non-commutativity of the k → 0 and A → 0 limits.In Fig. 7, we have carefully analyzed the low k regime of the GR black hole.For k/µ = 10 −4 , the two diffusive poles are close to their k = 0 values (54) and (57).As we increase k, they get closer and collide, leading to the two sound modes of ( 45).This diffusionto-sound crossover happens when k G ∼ A 2 illustrating the non-commuting limits k → 0, A → 0 which means we can estimate the characteristic length scale of the interactions to be V int ∼ AG (see footnote 11).
Finally, to re-emphasize the underlying Bloch wave Umklapp physics, Fig. 8 shows the real part of the momentum-dependent optical conductivity σ(ω, k) in the ω, k plane.The right-hand plot is a zoomed-in version of the left-hand plot near the edge of the Brillouin zone k = G 2 .The gray dots are numerically obtained solutions of det KL = 0 for the same parameters, showing that the hydrodynamic description of the matrix (43) at order O(A 2 ) matches the data over the entire Brillouin zone.At low frequency, we see the expected sound mode ω ∼ c s k dominating the low frequency regime, but we also can see its interaction with the sound mode ω ∼ c s (G − k).They meet at the edge of the Brillouin zone k = G 2 and in the right-hand plot, we see the traditional level repulsion of Umklapp and the opening of a gap in the sound mode spectrum.

Conclusion
The crucial message of this paper is that hydrodynamic fluctuations in a periodically modulated background should be understood based on a Bloch wave analysis instead of simple plane waves.If the typical length scale of this modulation is sufficiently large and the amplitude sufficiently small, we can still use hydrodynamics to study the long-time response of the conserved charges.This is an old observation in neutral hydrodynamics, but deserves restudy for charged hydrodynamics given the novel experimental progress of observed hydrodynamic flow in electronic condensed matter systems [7][8][9].This is particularly so in the presence of a charged fluid, which introduces an additional intrinsic diffusive mode.The presence of a spatial periodic modulation introduces Brillouin zone copies, also for this additional mode, and due to Umklapp at the Brillouin zone boundary this higher Bloch mode mixes with the    54)).The corrections to the diffusion constants are smaller than our numerical accuracy for our choice of parameters, so we can simply ignore them here.long distance late time k = 0 sound modes.We showed how one can compute the explicit pattern and strengths of these mixings from the underlying hydrodynamics.As is standard but ever so useful in hydrodynamics is that the behavior of both the patterns and the strengths can be expressed in underlying thermodynamic quantities, notably the susceptibilities, combined with the transport coefficients.
An important feature of a periodic modulation -well known in the condensed matter physics context -is that it breaks translational symmetry.For a perturbatively small lattice the correction to the momentum pole can be interpreted as the momentum relaxation rate and our result agrees with the relaxation rate obtained through the memory matrix formalism, as it should.It is important to emphasize once more that even though there is one relaxation rate, this relaxation rate has two contributions Γ d and Γ η corresponding to the two longitudinal diffusive processes.A priori these can have different scaling in temperature. 14They also exhibit different scaling in the lattice wavevector G. Due to this, in systems with charge disorder parametrized as an averaging over many independent lattices, one of these terms will dominate.It is rather the other aspect of the periodic modulation -the presence of Bloch modes in higher Brillouin zones that we wish emphasize here.At finite density this includes an Umklapped charge diffusion mode.As we analyzed in a companion paper, this mode may be of relevance in condensed matter physics [10].The strange metal phase of high T c superconductors shows the development of a mysterious mid IR peak in the optical conductivity at temperatures T ≃ 300K; see e.g.[39,40].The phenomenology of this peak is almost exactly reproduced by a collision between the Drude pole and the Umklapped charge diffusion pole in a holographic model of the strange metal dual to the Gubser Rocha black hole [10].If it can be experimentally verified that charge transport in the strange metal is in fact hydrodynamical, this will be the explanation of that phenomenon.
Finally, we verified our results by numerically computing response functions in strongly coupled systems holographically dual to Reissner-Nordström and Gubser-Rocha black holes.The important feature is that hydrodynamics emerges naturally in holographic systems and is not an input.In the computed optical conductivities, we found precisely the poles matching those predicted by our hydrodynamics computation.As as function of varying parameters such as momentum and lattice strength, these poles show complicated behavior including pole collisions and level repulsion denoting various regime changes.
We conclude with emphasizing that the hydrodynamics description of those holographic systems remains valid throughout these collisions and level repulsions.This contrasts with recent studies on the validity of hydrodynamics postulated as a pole collision/level repulsion with a first UV (gapped) pole [41][42][43].Our result here shows that this identification has to be done with care.The Umklapped modes are also a priori gapped modes in the zero momentum limit k → 0. However, they remain modes of the conserved charges, can be fully captured a hydrodynamic description and play a different role from non-hydrodynamic UV modes.
The analysis carried in this paper crucially relied on a static background charge distribution to mimick the effects of a frozen ionic lattice.This ignores the effect of lattice vibrations.Including phonon modes would require a different setup.Moreover, the assumption of local thermal equilibrium rather strongly constrains the hierarchy of scales as ω, k ≪ G ≪ T .While the results we have achieved are rather general and only rely on the presence of global symmetries and periodicity -which would seem to imply this is valid for a wide range of metallic systems -one must remain cautious as to whether such hierarchy of scales is realized within physical systems.Dutch Research Council/Ministry of Education.The numerical computations were carried out on the Dutch national Cartesius and Snellius national supercomputing facilities with the support of the SURF Cooperative as well as on the ALICE-cluster of Leiden University.We are grateful for their help.

A Thermodynamics and susceptibilities
In this section of the supplementary material, we will briefly review some key thermodynamic identities related to the static susceptibilities.We will be interested in the conserved charges {δε, δn} and their associated sources {δλ ε , δλ n } = {δT /T, δµ− µ T δT }.Since we will be focusing on thermodynamics, we will only be interested in the equilibrium solution and therefore we will drop the X notation for background thermodynamics quantities.As a reminder, the susceptibility matrix in the (ε, n) sector is defined as while the momentum susceptibility is χ ππ = ε + P. Furthermore, we have the thermodynamic identity Using this relation and the first law dε = T ds + µdn, we have Looking back at (A.1), this means that χ nε = χ εn .When considering conformal matter in Sec. 4, we used that P = ε/d which is directly implied by the equation of state (58).In that particular case,

B Onsager relations
One of the important checks we must make that our dynamical system is well-defined is that it respects Onsager's relations.These can be derived by considering how the system behaves under time-reversal invariance.Given the anti-unitary operator T such that [H, T ] = 0, we can classify each of the operators associated to our hydrodynamical variables by their representation under this operator.For a given operator δX a , we will have T δX a (t, x)T −1 =η a δX a (−t, x) with η a = ±1.Denoting the retarded Green's function associated to a dynamical matrix K ab by G R a b , we have We can then see that, due to the anti-unitarity of T , In Fourier space, this means ĜR a b (ω, p, p ′ ) = η a ĜR ba (ω, −p ′ , −p)η b and specifically for our periodic background, we can write the Green's function as [22] The elements of N are simply the coefficients of the equations ( 6) written in terms of the sources and expanded in the basis (34).By using that ĜR = (1 + iωK −1 ) • χ, the (B.4) can then be written as Let us apply this to the conservation equations ( 12) for an in-going momentum p and an outgoing momentum p ′ N εn (p, p ′ ) = σ Q dq (−pq μ(q))δ(p + q − p ′ ) , (B.7a) N nε (p, p ′ ) = σ Q dq (pq μ(q) + q 2 μ(q))δ(p + q − p ′ ) , (B.7b) N επ (p, p ′ ) = − dq (pχ ππ,0 (q) + q ε(q))δ(p + q − p ′ ) , (B.7c) N πε (p, p ′ ) = dq (−p ′ χ ππ,0 (q) + q ε(q))δ(p + q − p ′ ) , (B.7d) We now want to check the Onsager condition for N using that η ε = η n = −η π = 1.This can be done as follows for the (ε, n) sub-sector while for the momentum-charge sector, we have N nπ (p ′ , p) = dq (−pn(q))δ(p ′ + q − p) , (B.9a) Finally, we only have to check the energy-momentum sector We see therefore that the Onsager reciprocal relations are obeyed by our equations (12).

C Second order corrections in lattice strength
Let us consider here a dynamical matrix K of size N × N .The modes of this matrix are given by the solutions to the polynomial equation P(ω) ≡ det K = 0, where we can write Suppose the coefficients a n = p a n,p A p have a power series expansion in a parameter A. We are now interested in perturbative solutions ω = ω + A 2 ω 2 , around a given A = 0 solution P( ω)| A=0 = 0. To do so we can define auxiliary polynomials P p (ω) = N n=0 a n,p ω n such that We can now expand the equation P(ω) = 0 in A at leading and subleading orders, and we find the following two conditions The first equation is simply the leading order of the mode when there is no lattice while the second equation gives us the subleading correction.Finally, all the coefficients a n,p are themselves polynomials in G, k which can be further expanded in order to get the corrections at higher order in momentum. 15o far, it was implicitly assumed that the modes have no degeneracy when A → 0 as that would imply that P ′ 0 ( ω) = 0.When that happens, the correction is given by higher order terms with This is the case for KL+ with ω = 0.However, it turns out that for this matrix, P 4 ( ω) = 0 and P ′ 2 (ω) = 0, so the two degenerate poles remain degenerate at this order in perturbation theory, with ω 2 = 0.A simpler way to see this is also to notice that the first line of KL+ is proportional to iω and therefore so will be det KL+ .Consequently, this sector admits an exact conservation mode and we can use our non-degenerate method on the 4 × 4 lower-right subblock of this matrix where there is no degeneracy left.We would then see that the other ω = 0 pole also remains unshifted ω 2 = 0.
All that is therefore needed to compute the corrections in (54) and ( 57) is to know the coefficients a n,p for a given matrix K.In the case of KL− , we have On the other hand, for the 4×4 lower-right sub-block of the matrix KL+ , we find the following coefficients for the determinant The coefficients for the determinant of the full longitudinal matrix (45) are too involved to be written down here but can be obtained in the exact same way.All these corrections were derived using Mathematica.

D Numerical computations in strongly coupled field theories dual to Reissner-Nordström and Gubser-Rocha AdS black holes: Setup D.1 Thermodynamics
In Sec. 4, we focused on the specific conformal hydrodynamics that emerges at long wavelength and low frequencies from the holographic dynamics of the RN and GR black holes.In equilibrium the thermodynamic equation of state of each is given by ( 59) and (60) respectively.From these, we can determine the charge density n 0 = ∂ P ∂ µ T as well as the entropy density s 0 = ∂ P ∂ T µ while the energy density just follows from conformal invariance and is given by ε 0 = 2P 0 .One can further compute the various susceptibilities χ ab,m appearing in the hydrodynamics expressions of Sec. 3. As a reminder from Appendix A, the conformal equation of state also imposes χ εε,0 = 2χ ππ,0 and χ nε,0 = χ εn,0 = 2n 0 .For the RN black hole, the various susceptibilities and thermodynamic quantities are Lastly, we need to know some information on the transport coefficients η and σ Q to compute the hydrodynamic response.
These can be determined in the momentum-dependent homogeneous systems through η = lim ω→0 1 ω Im G T x y T x y (ω, k = 0) and σ Q = lim ω→0 1 ω Im G J x J x (ω, k = 0).In the case of conformal-to-AdS 2 solutions like the RN and GR black holes, these expressions can be solved analytically for the two transport coefficients.The shear viscosity η saturates the minimal viscosity bound η = s 0 4π [44] while σ Q was computed for a wide class of scaling black hole solutions [45] and here is given by , for RN, (D.3a) , for GR.(D.3b)

D.2 Numerics
We briefly review here how we compute the optical conductivity in the 2+1 dimensional strongly coupled conformal field theory holographically dual to the RN and GR black holes in the presence of a lattice.More details about the numerical methods used to compute these is the emblackening factor and A t (z) = µ(1 − z) a U(1) gauge field.In the above expressions, z is the radial coordinate ranging from the AdS boundary at z = 0 to the horizon of the black hole at z = 1.The temperature of this black hole is T = 12−µ 2 16π . 16he GR black hole is similarly obtained by extremizing the Einstein-Maxwell-Dilaton action with the potentials Z(φ) = e φ/ 3 and V (φ) = 6 cosh φ/ 3 Its metric is dz 2 + g(z)(dx 2 + d y 2 ) , (D.7a)This model is similar to the Einstein-Maxwell model with the addition of a neutral scalar field φ which controls the strength of the U(1) charge through the potential Z(φ).A consequence of this is the ability to discharge some of the black hole charge near the horizon such that the extremal T = 0 solution of this GR black hole will have a vanishing horizon and therefore vanishing entropy S T =0 = 0.The distance from extremality is controlled by the parameter Q: it is related to the chemical potential through µ = 3Q(1 + Q) and the temperature of the non-extremal black hole is given by T = 3 1+Q 4π .To obtain backgrounds with an explicit lattice, we will allow for a more general ansatz for GR with h(z) and g(z) unchanged from (D.7), but every field Q i j , a t , ϕ is now a priori a function of x and z.We require these fields to be regular near the horizon 17 and that their UV behaviour at z = 0 recovers AdS asymptotics. 18Moreover, to encode the modulation of the chemical potential, we must impose the following boundary condition on the gauge field a t (z = 0) = 1 + Acos(G x) .(D.11) The system is then solved numerically for the unknown functions Q i j , a t , ϕ.
To compute the optical conductivity in the holographically dual field theory, we must consider small fluctuations on top of this spatially modulated background.We linearize the Einstein equations around our lattice background .In the language of our hydrodynamic setup in Sec. 2, this is akin to simply turning on an external electric field δE x with momentum k and frequency ω.This response is also solved for numerically.
The numerical solutions to these equations were obtained using a publicly available custom package [46] and computed on the Dutch national Cartesius and Snellius supercomputers with the support of SURF Cooperative.

Figure 1 :
Figure 1: Drawing of the expected position of the poles in the current-current correlator at finite k (left) and at k = 0 (right) based on the hydrodynamical predictions in terms of Bloch waves in Sec. 3. Generically there are 9 poles: the standard two sound modes plus a (charge) diffusive mode of charged hydrodynamics cross coupled with the n = 1 and n = −1 Umklapp copies of each.At k = 0, there is an emergent symmetry due to which the longitudinal current-current correlator only probes the first sector KL− that contains the 4 modes that are odd under inverting the lattice momentum G ↔ −G.

Figure 2 :
Figure 2: Density plot of (the logarithm of) the GR longitudinal conductivity log |σ(ω)| (top) and its argument Arg σ(ω) (bottom) in the complex frequency plane close to the imaginary axis for T /µ = 0.1, A = 0.05 and G/µ = 0.1, for four values k/µ ∈ {0, 0.001, 0.005, 0.01}.At k = 0, we see the poles ω (L−) Drude and ω (L−) D while at k > 0, the Drude pole splits into the two sound modes ω S,±,0 (denoted by a white ⋄) and the diffusion pole splits into ω D,±1 (denoted by a white •).The markers indicate the analytical position of these poles prescribed by our hydrodynamical derivation.A priori, a fifth pole ω D,0 at k > 0 also couples to the electrical current, but it is not visible on the range plotted.A more refined computation, does reveal it (Fig.3).

Figure 3 :
Figure 3: GR conductivity log |σ(ω)| in the complex plane close to the imaginary axis for T /µ = 0.1, A = 0.05 and G/µ = 0.1, varying k/µ ∈ {0.006, 0.008, 0.01}.We see a purely diffusive pole on the imaginary axis which matches the hydrodynamic diffusion pole ω D,0 (denoted by △).The area plotted is zoomed on the origin compared to Fig.2.There the diffusive pole was too small to be visible.

Figure 4 :
Figure 4: Argument of the GR conductivity arg σ(ω) in the complex plane for T /µ = 0.1, A = 0.05 and G/µ = 0.1, varying k/µ ∈ {0, 0.001, 0.005, 0.01}.At k = 0, we see the sound pole ω (L−) + .Its real part is precisely at c s G with c s = 1/ 2 in a d = 2 conformal fluid.At k > 0, it splits into the two sound modes ω S,±,1 (denoted by ).The markers indicate the analytical position of these poles prescribed by our hydrodynamical derivation.These poles are more difficult to observe in |σ| than those on the imaginary axis and are easier to see as jumps in the argument of the complex function.

Figure 5 :
Figure 5: Comparison between the position of the poles on the imaginary axis (points) and the analytical hydrodynamical formula (54) at k = 0, as a function of T /µ, and for A = 0.05 and G/µ = 0.1.This is done for GR on the left and RN on the right.The blue data is the Drude pole ω (L−) Drude and the red data corresponds to the Umklapped diffusion pole ω (L−) D (Eqs.(54)).The corrections to the diffusion constants are smaller than our numerical accuracy for our choice of parameters, so we can simply ignore them here.

Figure 6 :
Figure 6: Comparison between the position of the poles ω D,0 and ω D,±1 on the imaginary axis (points) and the analytical hydrodynamical expressions (45) at k/µ = 0.01, as a function of T /µ, and for A = 0.05 and G/µ = 0.1.This is done for GR on the left and RN on the right.The corrections to the diffusion constants are smaller than our numerical accuracy for our choice of parameters, so we can simply ignore them here.

Figure 8 :
Figure 8: (Left) Re σ(ω, k) plotted in the (k/G, ω/µ) plane for A = 0.05, T /µ = 0.1 and G/µ = 0.1.(Right) Zoom on the Brillouin zone boundary at k = G/2 (region indicated by a black frame from on the left-hand plot) showing the level repulsion and the gapped sound mode at the edge of the zone.The gray dots are the hydrodynamic prediction given by numerically finding the roots of the determinant of Eq. (43).

g
µν = ḡµν + δh µν e −iωt+ikx , (D.12a)A µ = Āµ + δb µ e −iωt+ikx , (D.12b) ϕ = φ + δψe −iωt+ikx , (D.12c)and solve for these fluctuations with infalling boundary conditions, corresponding to choosing the response sourced through the retarded Green's function.The response in the radial electric field F z x in answer to an oscillating source in the potential δ∂ t A x (ω) ≡ δb x keeping the other components sourceless19 evaluated in the limit z → 0 then translates through the holographic AdS/CFT correspondence into the longitudinal optical conductivity σ = lim z→0F z x ∂ t A x = δJ x δE x