Holography and hydrodynamics with weakly broken symmetries

Hydrodynamics is a theory of long-range excitations controlled by equations of motion that encode the conservation of a set of currents (energy, momentum, charge, etc.) associated with explicitly realized global symmetries. If a system possesses additional weakly broken symmetries, the low-energy hydrodynamic degrees of freedom also couple to a few other"approximately conserved"quantities with parametrically long relaxation times. It is often useful to consider such approximately conserved operators and corresponding new massive modes within the low-energy effective theory, which we refer to as quasihydrodynamics. Examples of quasihydrodynamics are numerous, with the most transparent among them hydrodynamics with weakly broken translational symmetry. Here, we show how a number of other theories, normally not thought of in this context, can also be understood within a broader framework of quasihydrodynamics: in particular, the M\"uller-Israel-Stewart theory and magnetohydrodynamics coupled to dynamical electric fields. While historical formulations of quasihydrodynamic theories were typically highly phenomenological, here, we develop a holographic formalism to systematically derive such theories from a (microscopic) dual gravitational description. Beyond laying out a general holographic algorithm, we show how the M\"uller-Israel-Stewart theory can be understood from a dual higher-derivative gravity theory and magnetohydrodynamics from a dual theory with two-form bulk fields. In the latter example, this allows us to unambiguously demonstrate the existence of dynamical photons in the holographic description of magnetohydrodynamics.


I. INTRODUCTION
The past decade has seen a resurgence of interest in developing a systematic understanding of hydrodynamics as an effective field theory, describing the relaxation of locally conserved quantities towards global equilibrium in terms of long-lived (low-energy) degrees of freedom [1][2][3][4][5][6][7][8][9]. While the formulation of a dissipative hydrodynamic theory from an action principle was a long-outstanding problem, the rapid development of its reformulation in terms of effective field theory, along with other formal approaches to its classification [10][11][12][13][14], were largely ignited and accelerated by the advent of gauge-string duality (holography) [15], in particular, its ability to describe hydrodynamics of strongly interacting states [16].
The lines of research that have sprung from these developments have led to a number of important applications, ranging from a vastly improved understanding of thermalization and hydrodynamization of the quark-gluon plasma resulting from heavy-ion collisions [17][18][19][20], a comprehensive formulation of the theory of magnetohydrodynamics [21,22], to an understanding of the dynamics of electrons in exotic 'strange' metals [23][24][25][26]. For recent overviews see [27,28]. Many of the applications listed here pertain to old problems in physics. As a result, various phenomenological hydrodynamic approaches to their resolution have been known for decades. Unfortunately, these phenomenological approaches often lack rigor. In particular, a strategy common to many of these attempts is a rather ad hoc coupling of fluid degrees of freedom (particle density, momentum density, velocity, etc.) to non-hydrodynamic degrees of freedom (magnetic field, chemical reactant concentration, etc.). Another is an explicit breaking of various conservation laws (energy, momentum, charge, etc.).
A classic example, which we will study in detail in this work, is the textbook formulation of magnetohydrodynamics (MHD) (see e.g. [29,30]). MHD combines a theory of fluid degrees of freedom that obey the continuity equation and the forced Euler (or Navier-Stokes) equation, while the dynamics of electromagnetic fields obeys Ampere's law (neglecting displacement current), Faraday's law and magnetic Gauss's law. Momentum conservation of the fluid sector is explicitly broken (forced) by the addition of an external Lorentz force and the electric field is commonly expressed in terms of the magnetic field via the boosted ideal Ohm's law in the limit of infinite conductivity (see [31]). Standard formulation of MHD can be seen as lacking a systematic coupling between the separated fluid and electromagnetic degrees of freedom, as well an understanding of (global) symmetries by which one can organize a theory of long-range excitations in plasmas.
These questions were addressed recently in [21] where MHD was reformulated and extended by using the language of higher-form symmetries [32]. As a result of the above-mentioned issues with the historical approach to MHD, the standard formulation of MHD explicitly breaks the gradient expansion; this is particularly acute when MHD is coupled to dynamical electric fields.
Nevertheless, while from a formal effective field theory point of view such a theory should be viewed with suspicion, standard MHD and its simple phenomenological extensions make a number of extremely successful predictions about the dynamics of complicated astrophysical plasmas and processes in fusion reactors.
Another classic example of a system with an explicitly broken conservation law is fluid dynamics with weak momentum relaxation [28,33,34]. For example, such systems are known to describe both the drag of the Earth's surface on the atmospheric fluid [35], and the hydrodynamics of electrons in clean metals [26].
In light of these and other examples, the main goal of this work is to elucidate a formal approach to studying hydrodynamic theories with weakly broken symmetries, both from the point of view of field theory and holography. We will show that a large number of existing phenomenological theories can be precisely understood within this framework.
Hydrodynamics is a theory which is valid on length and time scales that are long compared to the mean free time τ mft and mean free path mfp [36]. Since hydrodynamics is a gradient expanded theory, it is convenient to express these statements formally as τ mft ∂ t 1, mfp ∂ i 1, (1.1) expressing the fact that gradients of relevant fields need to be small compared to two inverse time and length scales. The hydrodynamic equations of motion take the schematic form where ρ A correspond to the conserved densities, and J i A to conserved currents which can be expanded in a gradient expansion series of higher-order hydrodynamics (see e.g. [13,37]). The expectation values ρ A are computed in the equilibrium state of the system, e.g. the thermal equilibrium ρ A = Tr[ρ A e −βH ], with β = 1/T . Now, suppose that there exists a single non-conserved operator P in the theory, which has a relaxation time (decay rate) τ 1 , so that P (t)P (0) ∼ e −t/τ 1 , while all other 'orthogonal' operators have a much smaller relaxation time {τ 2 , τ 3 , . . .} τ 1 . 1 Formally speaking, the hydrodynamic 1 The notion of operator orthogonality here can be formally understood in terms of thermodynamic susceptibilities, see e.g. [28].
degrees of freedom do not include P , and τ mft τ 1 . However, with the addition of a single degree of freedom associated to P , we may replace the equation of motion (1.2) with ∂ t ρ A + ∂ i J i A [ ρ A , ∂ ρ A , · · · , P , ∂ j P , · · · ] = 0, (1.3a) ∂ t P + ∂ i J i P [ ρ A , ∂ ρ A , · · · , P , ∂ j P , · · · ] = − P τ 1 . (1.3b) The resulting equations of motion may be extremely accurate so long as {τ 2 ∂ t , τ 3 ∂ t , . . .} 1, which is a parametric improvement over the hydrodynamic expansion without P as a degree of freedom. As will be discussed in the main body of this work, there is a well-developed theory for explicitly computing τ 1 , either from field theory or from holography, as the answers are known to exactly match [38].
In many cases, τ 1 is a divergent function of a single parameter, e.g. the coupling constant λ, or a density of impurities that break translational invariance: τ 1 ∼ λ −2 . At strong coupling (λ 1), typically, τ 1 ∼ τ 2 τ mft . Thus, hydrodynamics is a very good description of low-energy (late time) dynamics. On the other hand, at similar length scales, but at weak coupling (λ 1), τ 1 τ 2 , so the dynamics of P can become important when τ mft ∼ τ 1 . As we will review, in the regime of λ 1, there is a controlled prescription for computing τ 1 so long as P is one of a finite number of long-lived quantities. Obtaining (1.3) using this prescription, hydrodynamics can be improved to systematically account for weakly broken symmetries with approximate conservation laws. It is straightforward to generalize (1. 3) to the case where a finite list of operators P α is long-lived.
For lack of a better phrase, we will call the theory in (1. 3), in which an approximately conserved quantity is treated on the same footing as exactly conserved quantities, a quasihydrodynamic theory.
The purpose of this work is to make a three-fold contribution to an already existing theory of quasihydrodynamics, which should help in unifying a number of past results under a common language as well as establish a systematic way to study such theories in the future. Firstly, we point out that-at least within linear response-a large number of well-known phenomenological theories are quasihydrodynamic: these include not only the momentum-relaxing fluid, but also magnetohydrodynamics and plasma physics with dynamical photons, simple models of viscoelasticity, (at least in some cases) the Müller-Israel-Stewart (MIS) theory of relativistic hydrodynamics, and (quantum) kinetic theory (in some respects). In every case, we identify both the approximately conserved quantities and the perturbatively small parameters which govern their decay rates. Some of these identifications are novel and may lead to new insights into old phenomenolgoical theories.
Secondly, we note that many quasihydrodynamic theories exhibit a universal "semicircle" law in which a hydrodynamic diffusion pole collides with a quasihydrodynamic pole to create a propagat-ing wave. A well-understood example is the formation of sound waves in a momentum-relaxing fluid at frequencies τ −1 1 . Examples that (to the best of our knowledge) have never yet been understood as quasihydrodynamic include transverse sound waves in elastic solids and electromagnetic waves in plasma physics. Thirdly, we study strongly coupled quantum theories with holographic duals [28] and describe how quasihydrodynamics arises from the bulk perspective. Within linear response, we show that the quasihydrodynamic regime exists and that quasihydrodynamics can be resummed to all orders in ωτ 1 by a controlled bulk computation of the quasihydrodynamic correlation functions and equations of motion. This allows us to-for the first time-present an unambiguous identification of a photon in a holographic plasma, thereby justifying claims made in Refs. [31,39]. We also demonstrate how in specific holographic model, MIS phenomenology and quasihydrodynamics with an approximately conserved stress tensor can arise. Earlier calculations along similar lines can be found in [40].
The outline of this paper is as follows. In Section II, we discuss the general quasihydrodynamic framework, and explicitly show that theories including MIS theory, magnetohydrodynamics and viscoelasticity are quasihydrodynamic in certain controllable limits. This part of the paper will serve as a more detailed summary of our results. Section III summarizes our holographic algorithm for analytically computing quasihydrodynamic equations of motion for the boundary theory, which we expect will find broad applicability to similar holographic problems. Sections IV and V deal with holographic theories of magnetohydrodynamics and higher-derivative Einstein-Gauss-Bonnet gravity, respectively. In each case, we show analytically how the quasihydrodynamic limit arises.

II. HYDRODYNAMICS WITH WEAKLY BROKEN SYMMETRIES
In this section, we outline how our framework of quasihydrodynamics for systems with approximately conserved quantities can be used to understand a large number of phenomenological theories known from past literature. In particular, we begin by giving precise formulas for the phenomenological τ 1 introduced in (1.3). We will then describe multiple examples of quasihydrodynamic theories and discuss the consequences of the weakly broken symmetry on the quasihydrodynamic modes.

A. Linear response
Let us first summarize what is known about the equations of motion of a quasihydrodynamic theory, within linear response. Suppose that the many-body Hamiltonian H 0 admits a number of local conservation laws associated with charge densities ρ A and ρ a : One of the ρ A is always the energy density. The operators associated with ρ a are also conserved, i.e. their charges commute with the Hamiltonian H 0 . Let us now perturb the Hamiltonian as so that ρ a no longer commute with the full Hamiltonian H: The ρ a quantities are our approximately conserved quantities. The quasihydrodynamic expansion is a derivative expansion in which-as we will see-∼ ∂ α will scale with derivatives.
Let µ A and µ a denote the thermodynamic conjugates (i.e., generalized chemical potentials) to ρ A and ρ a respectively, and let denote the susceptibility matrix of the hydrodynamic operators. The susceptibility matrices χ Ab and χ ab are defined in a similar manner. Suppose that when = 0, the hydrodynamic equations where ρ A,a and J i A,a are implicitly functions of µ A,a and its derivatives. One can show that when = 0 [41], , which has been rescaled out). Note that N ab = −N ba [41]. All dissipative contributions to quasihydrodynamics are contained in the matrix M ab , proportional to the spectral weight ofρ a .
The left-hand-side of (2.6) can also be expanded as a power series in , but as we will see, such terms will always be subleading compared to the orders in the derivative expansion in which we are interested.
If N ab = 0, then we take ∼ ∂ in the gradient expansion. An example of such a system is applying a small, non-dynamical external magnetic field B to a charged fluid, which breaks momentum conservation in two spatial directions. In this case ∝ B, and the approximately conserved operators a are two spatial momenta: e.g. P x and P y . The M ab corrections are subleading in the derivative expansion and contribute at the same order as viscosity.
If N ab = 0, which is the case we will focus on in this paper, then we take 2 ∼ ∂ in the gradient expansion. A similar observation was made in [42]. In this case, M ab is treated as a first-order term in the gradient expansion.
B. Diffusion-to-sound crossover A classic and possibly simplest quasihydrodynamic model is an example of the diffusion-to-sound crossover that interpolates from Fick's law of diffusion at low frequencies ω(k) to a propagating, sound-like linear (in k) waves at high frequencies ω(k). This example is illustrative, and we will see how it arises in a diverse set of physical systems in later sections.
One historical motivation for this model is as follows [43]: consider a diffusion equation governing, e.g., the magnetization in a spin chain with a locally conserved spin S = σ z : Here, D is the spin diffusion constant; the hydrodynamic spin current is J = −D∂ x S + · · · , where · · · denotes higher-derivative corrections. One can compute hydrodynamic Green's functions that follow from (2.8) [44], which accurately describe two-point functions in the quantum theory at small k and ω. These Green's functions have a pole with the dispersion relation ω(k) given to first order by Such poles are often called hydrodynamic quasinormal modes. The leading-order correction to (2.9) arrises from third-order hydrodynamics with a term proportional to O(k 4 ) [13].
There are several reasons why a dispersion relation of the form of (2.9) is problematic. Its group velocity ∂ω/∂k is proportional to k, which means that at large k, the propagation of diffusive modes is superluminal and thus acausal. A related problem is that the hydrodynamic Green's functions fail to obey microscopic sum rules, thus breaking unitarity [43]. A sum rule typically The dispersion relation for a single diffusive mode (2.9) is then replaced by a pair of modes that follow from the quadratic polynomial equation for ω, 11) and have dispersion relations At small momentum k, (2.12) gives a diffusive mode (2.9) and an additional massive gapped mode with the gap controlled by the inverse relaxation time 1/τ : The existence of an additional mode is a direct consequence of introducing a new dynamical field J with ∂ t J in (2.10), which makes the final (determinant) equation (2.11) quadratic in ω. At high k, the dispersion relation becomes linear which is why it is often said, somewhat imprecisely, that such modes become sound modes at large k. This dispersion relation is causal so long as its group velocity is The full dispersion relations from Eq. (2.12) exhibits the simplest signature of quasihydrodynamics: the collision of the two poles on the complex ω axis as a function of increasing momentum k. This collision occurs at which introduces a length scale into the problem above, well below which the strict hydrodynamics (2.8) applies; otherwise the quasihydrodynamics (2.10) applies. For the theory presently under consideration, the collision is plotted on the dimensionless complex ωτ plane in Figure 1. We also plot the real and imaginary parts of dimensionless ωτ as a function of kτ in Figure 2. Note that the behavior of the imaginary part of ωτ displays the "semicircle law" mentioned in the Introduction.
In terms of the quasihydrodynamic framework of this work, the physical motivation behind Eq. (2.10) can be understood as promoting the spin current operator itself to being long-lived.
There is no generic reason for this to occur, but when it does, then the present framework becomes applicable. As summarized in the Introduction, what we claim is that then, this type of a pole collision should be seen as a signature of the presence of an approximately conserved quantity. The existence of such an approximate symmetry can then either arise from microscopic dynamics or exist accidentally at certain special points in the parameter space of couplings and other tuneable parameters.
By following an argument of Kovtun [45], we end this section by asking the question of whether the gapped mode that causes the collision is real or an artifact of the expansion? Imagine an where we assume that |∂ t S| ∼ |∂ x S| ∼ 1. When treated ad verbum, this equation gives rise to exaclty the two modes of Eq. (2.12). However, one may perform a field redefinition Since, in hydrodynamics, S has no microscopic definition, it is just as good of a field as S . Physical spectrum should thus be invariant under the field redefinition (2.18). We can then choose α = τ , so that to O(k 2 ), there is only a single unambiguous diffusive mode, ω = −iDk 2 + O(k 3 ), which is invariant under field redefinitions.
While one may conclude that the gapped mode ω + (cf. Eq. (2.13)) is unphysical within hydrodynamics, the conclusion that one must reach in the quasihydrodynamic theory of Eq. (2.10) is different. In Eq. (2.10), S and J are independent fields, equations contain only first-order derivatives and most importantly, the theory can no longer be expressed as a gradient expansion.
Hence, there is no reason to expect that field redefinitions of the type S = S + ∆S and J = J + ∆J , where ∆S ∼ ∆J ∼ should keep the physics invariant. In the absence of a large scale separation between τ mf t and τ , such field redefinitions cannot be used to remove the gapped mode ω + and therefore, we conclude that gapped modes such as ω + are physical when interpreted within the context of quasihydrodynamics, which contains more degrees of freedom than hydrodynamics.

C. Müller-Israel-Stewart theory
As for the second example, we turn our attention to the Müller-Israel-Stewart (MIS) theory of relativistic hydrodynamics [46][47][48], in particular, to its complete form, which is a direct extension of second-order hydrodynamics [10] (the BRSSS theory). In the language of Schwinger-Keldysh field theory, the study of MIS theory was initiated in Ref. [6]. In order to place the MIS theory within the context of quasihydrodynamics and theories with approximately conserved currents, let us consider the four-dimensional, neutral and conformal stress-energy tensor expanded to second order in derivative expansion [10]: where D = u λ ∇ λ is the longitudinal derivative, ∆ µν = u µ u ν + g µν , σ µν = 2∇ µ u ν is the shear stress tensor, Ω µν = ∇ µ u ν − ∇ ν u µ the vorticity, R ρµνσ the Riemann tensor of the manifold the fluid occupies, and the bracket A µν denotes the transverse, symmetric and traceless part of A µν : The coefficients ε and p are the thermodynamic energy density and pressure, with ε = 3p, η is the shear viscosity and ητ Π , κ, λ 1 , λ 2 and λ 3 are the second-order transport coefficients.
The conservation of the stress-energy tensor T µν , generates a coupled set of four hydrodynamic partial differential equations for the three independent components of u µ and any scalar function, e.g. ε, or alternatively, the near-equilibrium temperature field T (x), which in general have acausal solutions, just as the simple diffusive case in Section II B. 2 These equations correspond to the purely hydrodynamic set of equations of the (1.2) type, with only strictly conserved energy and momentum. Hereon, we will only consider linearized solutions of (2.22)-the quasinormal modes-in flat space.
Just as in Section II B, the full set of diffusive and propagating sound modes suffers from acausal behavior. Similarly, this problem can be cured by treating Π µν = −ησ µν as an independent set of five degrees of freedom (Π µν is transverse, symmetric and traceless) and rewriting T µν in terms of 2 For some recent progress on formal discussions of causality in fluid dynamics, see [49].
Π µν as T µν = εu µ u ν + p∆ µν + Π µν , (2.23) and introducing an independent equation of motion for Π µν . In this work, we will only consider the linearized stress-energy tensor. What we find is the following set of nine quasihydrodynamic differential equations: which have the form of (1.3), or in terms of the diffusive example from Section II B, Eq. (2.10).
The coefficient τ Π plays a role of the relaxation time for the extra Π µν degrees of freedom.
We now look for the dispersion relations ω(k) of the quasinormal modes by solving the system of equations (2.24a)-(2.24b) in the transverse (shear) and longitudinal (sound) channels. The quadratic and cubic polynomial equations for shear and sound channels, respectively, are shear: where c s = 1/ √ 3 is the speed of sound of a conformal fluid in four dimensions. Notice that the equation (2.25a) in the shear channel is identical to Eq. (2.11) from Section II B, up to the identification of τ = τ Π and the diffusion constant with the usual ratio where s is the entropy density. The solutions of (2.25a) are thus again with the same small-k and large-k characteristics as the pair of modes in Section II B-the collision of the diffusive and gapped poles at critical momentum (2.28) After the collision, at large k, the speed of propagation becomes (cf. Eq. We refer the reader to Figures 1 and 2

. 3
In the sound channel, the cubic (in ω) polynomial equation (2.25b) gives three modes with dispersion relations that can be found explicitly, but we do not state them here in full. The modes consist of a pair of hydrodynamic modes with the sound dispersion relation at small k and a gapped mode. In the small-k expansion, they are We note that because of a single relaxation time τ Π , the gap is the same in both channels. As a result, at k = 0, when the SO(3) rotational invariance is restored, both shear and sound channels exhibit only a single gapped mode each with ω = −i/τ Π . At large k, We plot the full dispersion relations in Figures 3 and 4. The precise motion of the quasihydrodynamic modes in the complex plane depends on D and τ Π [51].
As in the previous subsection, it is tempting to associate this theory with quasihydrodynamics, yet, a priori, there is no reason why Π µν would correspond to an approximately conserved quantity.
Moreover, since the "UV completion" of the MIS theory is highly phenomenological, it is difficult to understand the precise microscopic origin of Π µν or how the theory should be correctly extended.
For this reason, it would be highly desirable to have a well-defined microscopic way of deriving the MIS theory. As we will show in Section V, the structure of the MIS theory and its equations can in fact be systematically derived from holography, using dual higher-derivative theories, such as the Einstein-Gauss-Bonnet theory [51][52][53]. This will be done to first order in fluctuations of the equilibrium fields. Indeed, the fact that the structure of MIS should be reproduced by higher-derivative theories could be anticipated from the studies of coupling-dependent properties of thermal spectra in [51,53]. It was shown there that the coupling between hydrodynamic and new, purely relaxing modes at intermediate coupling gives rise to the polynomial equations (2.25) and thus dispersion relations (2.27) and (2.30). The relaxation time τ is controlled by the coupling constant, and in the regime of large higher-derivative terms-weak field theory coupling-the new modes become long-lived with τ T 1. In Section V, we will explicity demonstrate how Π µν arises from dual gravitational perturbations, and derive the (linearized) structures through which it couples to T µν in equations of motion-i.e., the MIS equations (2.24).

D. Magnetohydrodynamics
We now turn our attention to the second main example of a quasihydrodynamic theory studied in this work: magnetohydrodynamics of a plasma with dynamical photons. Plasma is an ionized gas, which is charge-neutral at long distances-the long-range electric force is Debye screened and in the plasma phase, photons are massive. Nevertheless, the constituents of the plasma interact electromagnetically. While its long-distance charge neutrality implies that the equilibrium electric field E = 0, the equilibrium magnetic field B can be arbitrarily strong. In its standard formulation [29,30] (see also Ref.
where v is the velocity, ρ is the mass density, p the pressure and µ 0 the vacuum permeability. In Faraday's law, E is completely fixed by the assumption of the idea Ohm's law E = −v × B and the leading-order dissipative correction to this relation are suppressed by 1/σ where σ is the conductivity. Gauss's law plays no role in ideal MHD. For details, see Ref. [31].
The theory of MHD was recently reformulated by using the language of higher-form (generalized) global symmetries in [21]. 4 All known MHD equations and their systematic dissipative extensions follow from the realization that plasma possesses a global U (1) one-form symmetry associated with the conservation of the number of magnetic flux lines. As a result, the theory has a two-form conserved current J µν so that the full set of equations of motion is [21] In terms of the microscopic photon field A µ , the two-form current is The holographic dual of MHD, which we will study in Section IV, was constructed in Ref. [31]. 6 Unlike standard MHD, a symmetry-based classification of the stress-energy tensor T µν and the two-form J µν permits a systematic gradient expansion with all possible transport coefficients. 4 The relation between MHD formulated in the language of higher-form symmetries and MHD written directly in terms of electromagnetic fields was discussed in [22]. A recent work [54] claimed that writing down a consistent hydrodynamic partition function for MHD (on a thermal circle) requires the addition of an extra scalar degree of freedom. For a hydrodynamic theory with such a mode, the relation between the two different formalisms requires additional considerations, which were worked out in [54]. 5 For the construction of a theory of fluids with q-form symmetries, see Ref. [55]. 6 For a study of different aspects of the same holographic theory dual to a state with a one-form symmetry, see [39].
Moreover, the equation of state of the plasma can be an arbitrary function of temperature and the strength of the magnetic field, parametrized by the chemical potential µ conjugate to the number density of the magnetic flux lines ρ. In particular, the pressure can be an arbitrary function p(T, µ).
Beyond µ, the other hydrodynamic fields are the usual T , u µ and the spacelike vector field h µ .
The vectors are normalized in the following way: u µ u µ = −1, h µ h µ = 1, u µ h µ = 0. The relevant thermodynamic relations are ε + p = sT + µρ and dp = s dT + ρ dµ. Explicitly, to first order in the gradient expansion (see Ref. [21]), with η ⊥ , η , ζ ⊥ , ζ , ζ × , r ⊥ and r the seven transport coefficient in a charge conjugation symmetric state (see also Ref. [22]). The projector transverse to u µ and h µ is As in standard MHD, the electric field is not a dynamical variable. In the formalism of [21], using the microscopic relation between J µν and F µν , we can write the (relativistic) magnetic and electric fields in the fluid's comoving frame as where the electric field E µ is proportional to one-derivative tensors m µ and s µν and thus, the two resistivities r ⊥ and r -a generalized (inverse) Ohm's law. Note that in the above relativistic definition, the electric field E µ is defined with respect to the rest frame of the fluid. For this reason, in equilibrium, E µ = 0. This should be contrasted with electric and magnetic fields measured in equilibrium of the plasma, i.e. when u µ = (1, 0). Adopting the definitions in (2.35), we write When perturbed away from the equilibrium, i.e. for u µ = (1, v) with |v| 2 1, then the constitutive relations for J µν together with the definitions (2.36) give rise to the ideal Ohm's law used above, E = −v × B, and that the conservation law ∂ µ J µν = 0 automatically gives Faraday's and magnetic Gauss's laws. Also, note that the identifications (2.35a) and (2.35b) require the state to be perturbatively connected to the vacuum so that the "concept" of a microscopic photon field A µ exists.
It is clear from the present discussion that the effective theory of MHD is only valid so long as the electric field in a plasma is small-of the order of momentum. Therefore, for example, an initial state with a separation of positive and negative charges, and thus a large initial electric field that may lead to plasma oscillations cannot be described by MHD. Moreover, such scenarios would necessarily break the gradient expansion. Additional degrees of freedom are required in the effective theory. These degrees of freedom are gapped massive photons. Consider for example a thermal plasma with √ B/T 1 in quantum electrodynamics (QED). Then, to leading order, the Debye mass is m 2 D = 1 3 e 2 T 2 . In the (extreme) low-energy limit, k/m D ∼ k/T 1, so the spectrum of the effective theory contains no gapped photons. This is MHD. However, when either the electric field or k become comparable to T , the inclusion of massive photons is needed.
One may attempt to use the MIS-type procedure from Section II C to perform a simple phenomenological extension of the MHD equations (2.33a)-(2.33b) to include additional massive degrees of freedom. Since dynamical electric field arises from the charge sector, we can extend the two-form as where J µν is an anti-symmetric two-form, which we treat as encoding massive degrees of freedom that are independent of u µ , h µ , T and µ. In order to keep the definition of electric and magnetic fields consistent with Eqs. (2.35b) and (2.35a), we introduce the following constraint which enforces Eq. (2.35a). Hence, J µν contains three independent degrees of freedom and governs the dynamical electric field (cf. Eq. (2.35b)): The simplest quasihydrodynamic equations that we can write down are then with J tz = ρ = B z setting the strength of the magnetic field in equilibrium. Using Eq. (2.36), it is then straightforward to show that at the leading order in derivatives, Eq. (2.40c) becomes the advertized Ampere's law in a plasma: for perturbations that only depends on t and z. The right-hand-side is nothing but the boosted Ohm's law in the presence of a magnetic field. It turns out that the holographic model of [31,39] precisely encodes the above equation for any strength of the magnetic field, including B = 0. This will be shown in Section IV.
To investigate the simplest prediction of this theory, we study the corrections from new modes to the transverse Alfvén waves (see [21]). To this end, we perturb to first order the (quasi)hydrodynamic fields with the Fourier decomposition e −iωt+ikx sin θ+ikz cos θ and compute the spectrum in this channel. Note that J µν has a zero equilibrium value (E = 0 in equilibrium).
Again, the magnetic field is chosen to point in the z-direction. θ denotes the angle between the background magnetic field and momentum; without loss of generality we set k y = 0. In the (transverse) Alfvén channel, the only non-zero fluctuations are δu y , δh y , δJ xy and δJ yz . Instead of a pair of Alfvén modes, we now obtain a cubic polynomial for the quasinormal modes, which has the form where the anisotropic combinations of viscous and resistive transport coefficients are defined as In the small k expansion, the three modes take the form The pair of modes ω 1,2 is the pair of (non-perturbative) Alfvén waves from [21] and ω 3 is the new gapped mode. The form of the dispersion relations is completely analogous to the sound modes in MIS theory, cf. Eqs. (2.30a)-(2.30b), with the gapped mode ω 3 having the k 2 term be twice the dissipative O(k 2 ) contribution to ω 1,2 . Note that because we only extended J µν and not T µν , this contribution is purely resistive. Because of the lack of UV completion of T µν , which would introduce a second relaxation time, the large k limit is still acausal.
The conclusion that can be drawn is precisely the same as in the MIS theory. While the hydrodynamic part of the theory-MHD-can be constructed unambiguously, the construction of its MIS-type quasihydrodynamic extension, while consistent with expectations, is highly phenomenological. In an analogous manner, [22] extended MHD by adding the simplest E 2 terms to the partition function and treating E µ as a dynamical field. This allowed them to find Langmuir waves in a plasma with non-zero charge density. However, as we will see in Section IV, the theory of MHD and its systematic extensions can be derived directly from holography, in this case by using the holographic dual of MHD constructed in [31].
Now that we have discussed the need to extend MHD to a quasihydrodynamic theory, the question that remains is what is the approximately conserved current when τ ∼ τ mf t ? The answer is simple. In such systems, there exists an approximate one-form symmetry associated with an approximate conservation of a number of electric flux lines crossing a co-dimension two-surface. In the limit of τ → ∞, in which the mass of the new mode becomes zero-i.e., the photon becomes massless-we arrive at an electromagnetic vacuum state. In terms of Eqs. (2.40a)-(2.40c), the two terms in (2.40b) combine into the full dual electromagnetic field strengthF µν = 1 2 µνρσ F µν (u [µ h ν] = u µν is the magnetic part ofF µν [21]) with both electric and magnetic components, we can write F = dA, so that both dF = 0 and d F = 0.
E. Theory of elasticity, higher-form symmetries and topological aspects of quasihydrodynamic currents It is instructive to look at another well-known example of a quasihydrodynamic theory, which is in fact similar to magnetohydrodynamics of Section II D. This is the theory of elasticity, a historically extensively explored subject. While not normally phased in this language, the recent reformulation of elasticity from the point of view of higher-form symmetries [56] makes the connection to MHD (through the symmetry structure) and quasihydrodynamics rather transparent.
Moreover, it will enable us to further elaborate on the meaning of some approximately conserved currents, which can be interpreted from a purely symmetry-based perspective.
The theory of elasticity (see e.g. Ref. [57]) describes the dynamics of an elastic medium in a d-dimensional Euclidean space, embedded into the (d + 1)-dimensional spacetime x µ = (t, x i ). In its most common incarnation, the theory uses a dynamical variable, which is the displacement field φ I , where I = {1, . . . , d}. One can think of φ I as describing the coordinates of a given piece of the solid. The dynamics of φ I is governed by the conservation of momentum where P 0 I = ρ∂ t φ I and P i I = C ij IJ ∂ i φ J . The tensor C ij IJ is known as the elastic tensor (see e.g. Refs. [58][59][60] for reviews). One may also think of this theory as a theory of massless Goldstone boson that arise from spontaneously broken translational symmetry due to the presence of a lattice. Elasticity theory is then the low-energy limit in which the inhomogeneity of the state can be neglected [61].
For example, when all translational symmetries are broken by the lattice with some corresponding spatial group L, then the Goldstones φ I parametrize the quotient space R d /L U (1) d .
In the solid free of crystalline defects, such as dislocations and disclinations, there exists an additional conserved current when φ I are single-valued. In particular, This property plays the role of the topological Bianchi identity. As in the case of electromagnetism and MHD discussed in Section II D, where the Bianchi identity should be thought of as encoding conservation of the number of magnetic flux lines [21,32], Ref. [56] argued that the topological current J I should also be treated as a genuine conserved global current, which can facilitate the dynamics of an effective field theory of elasticity.
In fact, the periodicity of the Goldstones, i.e. φ I (x) ∼ φ I (x + ), where is the lattice spacing, gives a straightforward interpretation to J µ 1 ...µ d I as the current corresponding to a conserved number of domain walls of φ I in a direction perpendicular to the (µ 1 , . . . , µ d ) plane. Their conserved number density (the charge) is S 1 J I , which in conventional language counts the number density of the lattice sites of the elastic medium. The conservation of the higher-form current implies that there can be no crystalline defects in the system, which is consistent with the interpretation of (2.47) noted above, which stems from analytic properties of φ I . By focusing on the theory in 2 + 1 dimensions, one may work with a dualized description of the momentum operator in which In this picture, the conservation of momentum is equivalent to the conservation of magnetic flux lines in electromagnetism (point-like objects in 2+1 dimension) and the conservation of J I play a role of the conservation of electric flux lines [59,60].
For a theory with a non-vanishing equilibrium domain wall density J I = 0, it is clear that the current J I and the momentum P I are overlapped (in the terminology of the memory matrix formalism). 7 Among other things, this gives rise to a propagating mode in the transverse channel captured by linearized conservation laws. To be more explicit, consider the fluctuation of the and where C −1 is the inverse of elastic tensor. The conservation equation for J ⊥ then couples to the conservation of momentum in the following way: which can be combined into a wave equation for T ty upon decomposition of C −1 (see e.g. [60]).
This is in contrast with the spectrum of transverse fluctuations of fluids, which are controlled only by the conservation of momentum and, as a result, exhibit only (non-propagating) diffusive modes (see also [56] and [62,63]).
The higher-form current J I can become approximately conserved in the presence of small nondynamical dislocations. This scenario can occur in the melting of two-dimensional crystals [64][65][66][67].
In this situation, one can write down the broken Ward identity for J I as With these definitions in hand, we can write the system of quasihydrodynamic equations (2.10) as where the operator Q = j counts the number density of the domain walls in the system.
The number is conserved when τ → ∞.
• Higher-dimensional systems and domain walls: One can straightforwardly generalize the (1 + 1)-dimensional discussion to a (d + 1)-dimensional system. Imagine a theory with a conserved one-form j. Then, we define Q µ 1 ...µ d ≡ µ 1 ...µ d νj ν , which can be conserved or approximately conserved, as in We will discuss a class of holographic models which exhibit this property in Section III.
• MIS theory: In the shear channel of the MIS theory, the linearized equation of motion is identical to the theory with an approximately conserved d-form that appeared above in the context of the theory with elasticity. The higher-form current is J ⊥ = P ⊥ (there is no appearance of an analogue of the elastic tensor C µν IJ ). Similarly to the spin system discussed above, the antisymmetric structure of J ⊥ also implies that the spectrum in the scalar channel does not depend on the wave vector. In particular, in 3 + 1 dimensions, we will have becomes approximately conserved and the IR limit of the plasma quasihydrodynamic.

F. Kinetic theory
A more subtle appearance of the quasihydrodynamic formalism arises in (quantum) kinetic theory, within linear response [68][69][70][71][72][73][74][75][76]. Let f p (x) be the one-particle distribution function of a kinetic theory for particles of momentum p at position x, and let δf p (x) denote perturbations away from thermal equilibrium. The linearized kinetic equations read where W pq is the linearized collision integral, whose precise form we will not write here. As noted in [77] for a (quantum) kinetic theory of fermions, W pq is related to the spectral weight of the quasiparticle density operators. Thus, (2.55) exactly reproduces (2.7), suggesting that kinetic theory also fits into our broader framework of hydrodynamics with weakly non-conserved operators.
Here, the set of conserved and approximately conserved operators includes all quasiparticle density operators: e.g. c † p c p , where c † p and c p are quasiparticle creation/annihilation operators at momentum p.
There is an important physical distinction between (general) kinetic theory described here and other examples of quasihydrodynamic theories studied above. As emphasized in the Introduction, quasihydrodynamic theories typically have a small number of approximately conserved quantities, along with a "soup" of many degrees of freedom with much shorter lifetimes. In the present example, however, all of the interesting degrees of freedom obtain lifetimes that are comparable to the eigenvalues of W pq .
To the extent that kinetic theory is quasihydrodynamic, we anticipate that the quasihydrodynamic analogy will also persist to other integrable models in 1 + 1 dimensions, following the recent advances in "generalized hydrodynamics" [78][79][80]. We further expect that the quasihydrodynamic formalism may allow for a proper description of dynamics in such theories when they are weakly perturbed by integrability breaking deformations.
One common approximation in kinetic theory is that the relaxation to equilibrium is dominated by only a few non-conserved modes [10,81]. In other words, W pq will have a hierarchy of nontrivial eigenvalues and we focus on the smallest ones. In the context of the linearized MIS theory, the stress tensor itself is an approximately conserved quantity. This implies that can be approximately written in terms of only u µ , T and Π µν . In this expression, one conventionally approximates the collision integral in RTA [10,81]. Let us emphasize that in a generic kinetic theory, there is no reason to expect that X µνρ (T, u, Π) is not just as long lived as Π. What is special about the holographic models we describe below is that for certain parameter regimes and due to special microscopic reasons, Π is (alone) a long-lived quantity.
An alternative approximation in kinetic theory is that all non-zero eigenvalues of the collision integral are identical (see e.g. Ref. [77]). Such models will not generically recover the MIS phenomenology. To the extent that such models are quasihydrodynamic, every single mode δf p is quasihydrodynamic.

G. Other examples
There are a few further examples of quasihydrodynamics that we will not address in this paper, but which we note here for reference: the "zero-sound" physics of low temperature holographic probe brane models [40,[82][83][84], and quantum-fluctuating superconductivity [85]. Other systems of interest that may be re-examined in the future from the point of view of quasihydrodynamics include theories such as the one studied in [86].
As mentioned in the Introduction, recent experiments [23][24][25][26] have observed evidence for hydrodynamic flow of electrons. Unfortunately, this hydrodynamic flow is not exact due to the presence of impurities, phonons and Umklapp scattering, all of which can relax momentum away from the electronic fluid. As derived in [28,33], in the presence of momentum relaxation, the linearized hydrodynamic equations of motion for momentum density g x and energy density in a charge neutral fluid become exactly analogous to (2.10), with S replaced by and J replaced by g x . The conventional sound wave of a fluid is split into a diffusive mode and a "gapped" mode, as depicted previously in Figure 2.

III. OUTLINE AND SUMMARY OF THE HOLOGRAPHIC METHOD
The purpose of this section is to introduce a simple holographic calculation of quasihydrodynamics. The physical system that we study is a transition from diffusion to ballistic physics, analogous to Section II B. So our focus here is not on the physics, but on the technical strategy that we use to analyze the holographic model. The methods that we develop below are a streamlined version of the algorithm of [38,40]. For the remainder of the paper, we assume that the reader is familiar with holographic theories; we will not try to explain the basics of the correspondence.
Our goal is to study linear response in a large-N (matrix) field theory in d + 1 spacetime dimensions. We postulate that this theory is holographically dual to classical gravity in d + 2 spacetime dimensions. Single-trace operators which are rank-s tensors in the field theory are dual to rank-s fields in the bulk. In this section, we will be studying a gauge field A a in the bulk with The dual operator of this gauge field is a conserved U (1) current associated to a conventional (zero-form) symmetry. Suppose that the electromagnetic part of the action is where the tensor X abcd is anti-symmetric under a ↔ b and c ↔ d and symmetric under ab ↔ cd.
In the four-dimensional bulk, there are only six components which are These tensors X are functions of the background bulk fields. We will assume that the background geometry is isotropic and translationally invariant in the boundary spacetime directions: We also assume that on the background geometry, A a = 0. These assumptions suggest that X should only be functions of the background metric, i.e. of a(r) and b(r). Isotropy implies that X 1 = X 2 and X 5 = X 6 everywhere in the bulk. Furthermore, the regularity of the background solution implies a relation between t and r components of X abcd at the horizon: X 1 (r h ) = X 5 (r h ) and X 2 (r h ) = X 6 (r h ). Holographic electromagnetic bulk actions that can be cast in this form include the probe brane theory [40] as well as theories with higher-derivative couplings for F ab [53,[87][88][89]. Our convention will be that the boundary theory, which lives in the UV, is at r = 0; the geometry ends in the IR at a planar black hole horizon at r = r h . From general principles, near the horizon, while b(r) is finite. If the UV theory is conformal, then a(r) ∼ b(r) ∼ r −2 as r → 0. For simplicity we can assume these UV scalings in the discussion that follows, though the general algorithm we describe is not sensitive to this assumption.
Since A a couples quadratically in S EM , we conclude that the linearized equations of motion for bulk fields depend only on fluctuations δA a . Using boundary spacetime translational invariance, we may write and need only solve ordinary differential equations for δA a (r). For convenience, we will assume that k points in the x-direction. Let us first consider the equation of motion for δA t and δA x , which couple and determine the one-point functions δj t and δj x . The equation of motion for these modes are In general, even for the simplest geometries, these equations cannot be solved exactly.
However, we do not need the exact solution. We are only interested in the quasihydrodynamic regime, which (at least in holography) necessitates ω T. (3.7) Roughly speaking, in holographic models, the radial coordinate r corresponds to the "energy scale" at which physics takes place. At energy scales large compared to temperature, a field obeying (3.7) appears approximately static. Indeed, if we look at (3.6), as r → 0, all of the ω and k dependence is negligible. So we could perturbatively construct a solution in ω and k near the boundary. In contrast, very close to the horizon, (3.4) dominates the equations of motion. The near boundary expansion fails and instead the solutions must obey suitable infalling boundary conditions (things fall into black holes, not out of them). Indeed, (3.6) again becomes a simpler differential equation to solve in the limit r → r h .
Since we have two separate regimes in which the solution to (3.6) appears simple, we may attempt to perform a matched asymptotic expansion, as sketched in Figure 5. We will first determine the complete set of solutions to the bulk equations of motion-with ω = k = 0-in the outer region, which begins at r = 0 and extends in the direction towards the horizon. Then, we will construct a solution to (3.6) in the inner region, which extends from the horizon r = r h outwards, usually by a distance ∼ T −1 . This solution will be non-perturbative in ω, but also assume k = 0. Our claim is that the outer and inner regions interlap in an intermediate region: We do not have a proof that this is universally an intermediate region in which to match solutions, but this inequality will hold in every example studied in this paper. This intermediate region is close to the horizon, and so the matching coefficients between inner and outer solutions will generally depend on the near-horizon geometry.
Nevertheless, there are a few subtle differences between our implementation and in those that exist in the literature, which we will comment on as we proceed. The matching algorithm that we describe generalizes straightforwardly to much more sophisticated problems than have previously been tractable in the literature; indeed, our discussion of holographic magnetohydrodynamics in Section IV would be a highly non-trivial calculation using other previous holographic matching methods.
a. Outer region: In the region away from the horizon, we expand the solution for δA t and δA x as ω, k → 0 as: The radial functions Φ 1,2 (r) can be obtained by setting ω = k = 0 in (3.6), as mentioned previously: ds 1 a(s)X 5 (s) . (3.9) Both of these functions vanish at the boundary while, as one get closer to the horizon, Φ 1 (r h ) is finite but Φ 2 (r h ) diverges logarithmically. For future purposes, it is convenient to denote the finite part of the functions Φ 1,2 by φ 1,2 , and subtract out the logarithmic divergence as follows: ln a(r). (3.10) Hence, To determine the value of r at which this outer region ends, we can consider the first order in the ω, k expansion. This can be done either in the scheme where ω ∼ k ∼ or ω ∼ k 2 ∼ , where the small parameter 1. In either scheme, we can write   x is identical to those at zeroth order in , namely, This implies that δA [1] x ∝ Φ 2 (r). The fact that Φ 2 (r) diverges deeper in the bulk implies that the expansion breaks down when a (r h )X 5 (r h ) ln a(r) ≈ 1. (3.14) In other words, in the low frequency limit ω/T , the regime of the validity of the outer region solution can be expanded up to at least r ≈ r h − T −1 e −4πT /ω , as depicted in Fig. 5.
Using (3.6c) along with the results of the previous paragraph, we immediately see that (after performing an inverse Fourier transformation from ω and k): This is the exact Ward identity associated with charge conservation. It is one of the two equations of motion in the quasihydrodynamic regime. The other, "approximate" conservation law cannot be fixed without matching the solution into the near-horizon regime.
A few points are in order. Firstly, we have not used gauge-invariant "master field" variables.
While this approach is easy to use here, it becomes increasingly difficult in more complicated systems, where gauge-invariant variables may become difficult to relate to physical quantities of the boundary theory. Moreover, the equations of motion of the fluctuations in the radial gauge can be easily written down, even in a more complicated problem (as we discuss in Section IV). Secondly, our aim is to find the explicitly broken Ward identity and not simply the dispersion relation of the quasinormal modes. Thus, it is easier to simply work directly with variables whose boundary conditions relate to the boundary observables of interest. This also helps us avoid ambiguities which arise at quasihydrodynamic pole collisions. These features will be explicitly visible in the discussion below.
b. Inner region: In the region close to the horizon, the solution to (3.6) can be written as where the functions A  (1) a (r) is a pure gauge solution: it is absent in computations with gauge-invariant variables, but it is nevertheless useful to keep track of it. For example, we can determine both the source and the response when interpolating to the boundary in the radial gauge computation [97]. One may argue that since A (1) a (r) is pure gauge, it has to satisfy a first-derivative constraint and have no effect on physical quantities. However, on the operational level, one may also just substitute the

ansatz (3.16) into the equation of motion (3.6) near the horizon and solve for A
(1) a (r) and A (2) a (r). In order to do this, one can start by demanding that, upon the above substitution, the coefficients of divergent terms, such as {a(r) −1 , a(r) −1−iω/4πT , . . .}, have to vanish. This turns out to give all constraints on the near-horizon solutions at the order in ω and k to which we are working. In this case, the relevant constraints are t (r h , x µ ) = 0. (3.17) Observe that the gauge-invariant part of A (1) a vanishes when evaluated at the horizon, as regularity demands. Note also that we will be using the shorthand notation A  We first match the logarithmic pieces in the solutions in both regions. In the inner region, this is the O(ω) contribution arising from the infalling contribution to (3.16). We obtain where we used the fact that a(r) = 4πT (r h − r) + O(r h − r) 2 in our coordinates. Similarly, the matching condition for the finite part implies that These matching conditions together with the regularity condition imply that, up to first order in derivatives, we have One can recast this relation in the form of a weakly broken conservation law as in Eq.(2.10), with This computation is almost parallel to the charge diffusion in Einstein-Maxwell theory presented in [90] where one finds that τ T ∼ 1 and thus 1/τ ω in the hydrodynamic limit ω/T 1. In this regime, it simply means that the relaxation time ofĵ x is very small at large temperature and j t is not a long-lived operator anymore. Thus it is sensible to drop the time derivative term in (3.21) and obtain the constitutive relationĵ x = D∂ xĵt . However, in non-Maxwell theories as considered such as DBI or higher-derivative theories, one can show that the relaxation time can be tuned such that (τ T ) −1 1 [40,53,87,88]. In this case, one has to keep all the terms in (3.21) to properly capture the quasihydrodynamic behavior of the system.

d. Results:
The quasihydrodynamic equations of motion for this theory consist of (3.15) and (3.21): where χ j = 1/φ 2 (r h ) is a thermodynamic susceptibility which determines the coupling to the background gauge fieldâ.
e. Scalar channel: Similar analysis can also be carried out for the scalar channel involving δA y , which results in a finite life-time of the operatorĵ y . The relevant equation of motion is which yields the following solution in the outer region δA y =â y +ĵ y Φ 2 (r). (3.25) Here, Φ 2 (r) is the same function found in the previous (longitudinal) sector (3.6). Since here

A
(1) y = 0, we find that the matching conditions imply Note the absence of spatial derivatives. The parameter τ is the same as those in (3.22). We can also wrap the two relations in both channels, (3.21) and (3.26), into one approximately conserved two-form current Q µν = j µν , wherej µ = (ĵ t /v,ĵ x ,ĵ y ). Combining with the conservation of j µ , we finally arrive at the quasihydrodynamic equations Let us end this summary by quickly discussing how quasihydrodynamics arises in the presence of double trace deformation. Such a deformation implies that the definition of the physical source is a linear combination of the radially dependent and independent pieces (analogues ofâ andĵ in (3.21) and (3.26)). Depending on the double trace coupling, which determine the linear combination, it is possible to for τ to be large (τ T 1); in this case, the current becomes approximately conserved.
We will discuss an approximately conserved current of this type in Section IV.

MAGNETOHYDRODYNAMICS WITH DYNAMICAL PHOTONS
In this section, we apply the holographic algorithm of Section III to derive the quasihydrodynamic, low-energy excitations in a field theory with a conserved two-form current. The holographic dual of a strongly interacting theory with a one-form symmetry was proposed in [31,39]. In particular, [31] argued that this holographic setup furnishes a description of a strongly interacting field theory with a matter sector gauged under dynamical U (1) electromagnetism, thereby providing a holographic dual description of a plasma described by MHD in the limit of low energy. As a result of electromagnetism, the theory was argued to contain dynamical photons [31,39]. Hence, one should be able to use this holographic setup to not only describe MHD but also its (quasihydrodynamic) extension to a theory of plasma with dynamical, and potentially strong, electric fields.
At the level of linear response, we will explicitly demonstrate that this is achieved in holography.
Before delving into the details of the calculation, we will review a few essential features of this model. The bulk theory is composed of the Einstein-Hilbert gravity action, a negative cosmological constant and the two-form gauge field B ab , written as where H = dB, Λ is the UV cut-off and n a is the unit vector normal to the boundary. S bnd combines the boundary Gibbons-Hawking term and counter-terms for the gravitational part of the theory. A detailed exposition of the model, including the holographic renormalization, can be found in [31].
The diffeomorphism invariance of the metric and the gauge symmetry of the two-form bulk fields imply that the boundary duals possesses a conserved stress-energy tensor T µν and a conserved two-form current J µν . In particular, in the presence of an external two-form source b µν , where H, here, is the three-form field strength H = db of an external two-form gauge field. In terms of the generating functional, For the source to be physical, one has to formally impose that it be cut-off independent, i.e. that ∂b µν /∂Λ = 0, which results in the logarithmic running of the double-trace coupling κ with a Landau pole. In other words, one finds that 1 κ(Λ) = finite − ln(Λ/L), (4.5) where the finite part that sets the renormalized electromagnetic coupling needs to be imposed as the renormalization condition at a new scale L, or equivalently, by the choice of the Landau pole scale. In practice, the value can be chosen by hand. For details, see [31,56]. As we will see below, it is the (finite) value of the electromagnetic coupling that governs the relaxation time of the operator J zi , or the electric flux density, and sets the regime of validity of the quasihydrodynamic approximation with an approximately conserved electric flux density.
We take the background metric ansatz and two-form gauge field to have following equilibrium forms (cf. Eq. (3.3)): where h(r) parametrizes the density of the two-form conserved current, which equals the Hodge dual of the magnetic flux density. The functions a(r), b(r) and c(r) in the metric all asymptote to one at the boundary, r → 0, so that the geometry is asymptotically AdS 5 . At the horizon, b(r) and c(r) are regular and approach non-vanishing constants b(r h ) and c(r h ). The regularity conditions impose a stronger constraint on the "emblackening factor" a(r) and the gauge field h(r), which are forced to behave as close to the horizon at r = r h .
Finally, we note that in this work, we will study the transverse fluctuations where the index i denotes the (x, y) plane coordinates, with the plane being perpendicular to the z-direction of the background magnetic field lines.
The two subsection, which contain the bulk of the derivation of aspects of quasihydrodynamics in a holographic plasma can be summarized as follows: • First, in Subsection IV A, we consider the case with a zero background magnetic flux density.
At small momentum, we find that the system exhibits a diffusive mode, which will be referred to as magnetic diffusion, along with a massive, non-hydrodynamic mode with an imaginary gap set by the renormalized electromagnetic coupling. The system will then be shown to exhibit the collision of the diffusive and gapped poles described in Section II. We will find explicit analytic expressions (in terms of background bulk quantities and the U (1) coupling) for the relaxation time of the approximately conserved operator J zi and the asymptotic speed of the pair modes (after the collision) at large k. In the large-k limit and for small electromagnetic coupling, this speed will tend to the speed of light. We thus explicitly show the presence of photons in the plasma. We begin by considering a setup with h(r) = 0. In this limit, the function c(r) = 1. Moreover, the metric and the two-form field fluctuations are decoupled and we can focus on the equations of motion for two-form gauge field. In the radial gauge, they are δB ti + kra(r)δB zi = 0. a. Outer region: In the low frequency limit, we can show that the solution is where the functions Ψ 1,2 (r) satisfy the following first-order differential equations: Ψ 1 (r) = 1, ra(r)Ψ 2 (r) = 1, (4.11) with the UV boundary conditions at r = 1/Λ, with Λ/L 1, giving rise to logarithmically divergent behavior (the Landau pole), Ψ 1,2 (r = 1/Λ) = − ln(ΛL) + · · · . which is finite, and ln (a(r)) + ψ 2 (r), (4.15) where Both ψ 1 and ψ 2 are constructed so that they are finite everywhere in the bulk.
As before, (4.9c) guarantees that This is the exact conservation law for magnetic flux lines. The quasihydrodynamic mode for electric flux lines will arise by matching to the near-horizon region.
b. Inner region: Near the horizon, the general solution for the two-form field as a sum of a regular pure-gauge and infalling parts: where δB µi are regular function at r = r h . We proceed by substituting the above ansatz into the equations of motion near the horizon and demanding that the coefficients of terms containing a(r) −1 , and a(r) −n−iω/4πT for n > 0 vanish. This procedure implies several constraints on B ti and B zi . Firstly, demanding that the coefficients of a(r) −1 in both (4.9a) and (4.9b) vanish, gives (1) where we denote B µν (r h , x µ ) as simply B µν . Similarly, by demanding that the coefficient of a −2−iω/4πT in (4.9a) vanishes, we find that On the other hand, the function B (2) zi (r h , x µ ) is non-vanishing and satisfies a more complicated constraint. However, there is no need for us to solve explicitly for B (2) zi . It is only required to be finite at the horizon.
c. Intermediate region: In the intermediate region, we match the solutions from inner and outer regions. Staring from the inner region, we isolate the logarithmically divergent piece in δB xz by writing ln a(r) zi ln a(r), where we take Φ 2 (r) ≈ Φ 2 (r h ). The matching between solutions in the two regions implies that (1) ti = δB ti + δĴ ti (ψ 1 (r h ) + ln(r h /L)) .  To interpret this relation in the dual field theory language, we recall that the physical source is defined by the following linear combination of the coefficients of the r-independent term and the logarithmic term, as outlined in Eq. (4.4) (see also [31,39]): where M ≡ Λ exp(κ(Λ)) is an RG-invariant scale. The constraint (4.24) then becomes By turning off the source, we finally find the Ward identity for the approximately conserved current: The meaning of this equation is clear. This is a quasihydrodynamic relaxation equation for δJ yz , which corresponds to the electric field in the thermal plasma. The electric field is therefore not only induced through fluctuating magnetic fields, as in MHD, but is a genuine dynamical degree of freedom, with a relaxation time τ . A direct consequence is that the system indeed contains dynamical photons, as claimed in [31,39]. At large τ , the approximately conserved quantity is the conservation of the number of electric flux lines crossing a co-dimension-two surface. We note that in this expression, temperature (or r h ) sets the characteristic energy scale of the system. Mr h is the parameter that controls the renormalized electromagnetic coupling: For the case at hand, the AdS 5 -Schwarzschild black hole, we find that ψ 1 (r h ) = ψ 2 (r h ) = 0.
Hence, the speed of the propagation of photons is actually the speed of light in vacuum, which is independent of the electric charge: v = 1. As e −2 r 1, the time scale τ obeys τ T 1, and so the electric field is indeed a quasihydrodynamic mode. For this geometry, (4.29) This may be compared with the numerical estimate of this lifetime found in [39], which is about 3.5% larger.

d. Results
: We now collect our results. Using (4.17), (4.26) and (4.27), and defining δE x = δJ yz , δB y = δJ t y , we conclude that ∂ t δB y + ∂ z δE x = 0, (4.30a) These are precisely the dynamical Maxwell's equations in the presence of an Ohmic current J ∝ E, as discussed in Section II D. The conductivity of the plasma is proportional to 1/τ ∝ e 2 r and is small in the quasihydrodynamic limit. Indeed, as the electromagnetic coupling vanishes, the photons decouple from matter and there is no Debye screening. Hence, as T → 0 and e r → 0, the quasihydrodynamic mode-the photon-becomes arbitrarily long lived.

B. Alfvén waves and photons
We now consider the case with a non-zero background magnetic field, i.e. with h(r) = 0. The background solution is the magnetic black brane solution [98], which is the same background one considers in the case of external, non-dynamical magnetic field (see Ref. [31] for a discussion on the relation between holography with and without dynamical magnetic fields). The thermodynamic quantities, transport coefficients and the low-energy spectrum as a function of temperature and background magnetic field were computed numerically in [31,99] (see also [100,101]). Here, we show that the form of the quasihydrodynamic low-energy spectrum can in fact be obtained analytically.
The equations of motion for the background metric, i.e. for {a, b, c}, are not particularly illuminating and we will not write them here. However, it is important to note that the background equation for the background two-form gauge field is In other words, we can write the radial derivative of h(r) as where h 0 is a constant. We choose a gauge with h(r h ) = 0 and, as a result, it is convenient to write down the solution to (4.32) as This form is particularly convenient for analyzing the solution near the boundary, r ≈ 1/Λ 1.
The result is Similarly, there are two additional radially conserved quantities, Q 1 and Q 2 , such that dQ 1 /dr = dQ 2 /dr = 0, which arise from linear combinations of the xx-and tt-components, and yy-and xx-components, respectively. Namely, Using h(r h ) = 0, we find that 37a) where s is the entropy density of the system, as measured by the horizon area.
These relations are crucial to simplify the equations of motion for the metric h µν and the gauge field fluctuation δB µν , which consist of four second-order ODEs: and two first-order (constraint) equations: a. Outer region: The equations in the outer region can be organized into two decoupled sets: where {π ti ,π zi ,Ĵ ti ,Ĵ zi } are integration constants. The first-order equations (4.39) imply that ∂ tπ ti + ∂ zπ zi = 0, ∂ tδ J ti + ∂ z δĴ zi = 0, (4.42) which are nothing but the exact conservation laws for transverse momentum and the two-form current (magnetic flux).
Let us first look at the first pair of equations, i.e. Eqs. (4.40a)-(4.40b). Solving for h ti , we find This implies that are two integration constants of the equation (4.43). Similarly, using (4.40a)-(4.41b), we find that 1 (r), (4.45) where Ψ 1 and Φ 1 satisfy √ c 46) which are solved by This solution can be used to construct the other linearly independent solution to (4.46) via the Wronskian trick. Applying this method to Ψ 1 , we find that It is convenient to use the relation (4.46) to solve for Φ (2) 1 in terms of Ψ (2) 1 in order to avoid the unnecessary appearance of a logarithmic divergence. We find 1 (s), (4.49) where the value of Φ 1 (r = 0) is fixed by the first derivative of Ψ 1 at r → 0 through the first relation in Eq. (4.46). In summary, the solutions for h ti and δB zi in the outer region are 1 (x µ ), (4.50a) 1 , 1 (r h ) c(r h ) r h a (s)ψ (4.51b) To better understand the meaning of variablesĈ (1) 1 andĈ (2) 1 , it is helpful to rewrite them in terms of the zero magnetic field case. First, we need to consider how the fluctuations h ti and δB zi behave close to the boundary and defineĥ Upon substituting the above definitions into the outer region solutions and expanding them near the boundary, we find h ti =ĥ ti + 1 4 r 4 π ti + 4h 0 δB zi + O(r 5 ), (4.53a) where, in terms of the radial conserved quantity (4.36), Note also that the above near-boundary expansion reduces to the one in the zero background magnetic field case of Section IV A as we take h 0 = 0.
An analogous sequence of manipulations can be applied to the second pair of fluctuations, h zi and δB ti from Eqs. (4.41a)-(4.41b). The solutions are then 2 (r). (4.55b) The functions Ψ 2 and Φ 2 satisfy 2 . (4.57d) 1 and Φ (2) 2 are finite at the horizon while Ψ 2 is not. The latter solution can again be split into a finite and a divergent part as and Ψ (1) 2 . (4.60) Similarly, as before, we can relate the functionsĈ 2 to their h 0 = 0 counterparts aŝ Hence, this gives the following two near-boundary expansions: b. Inner region: The solutions in the near-horizon region can be written in the following form: where H and a(r) −2−iω/4πT ) vanish. Among many constraints that emerge from this procedure, demanding that the coefficients of a(r) −1 vanish implies that (1) Similarly, requiring that the coefficients of a(r) −2−iω/4πT in (4.38) vanish implies that where we have used a shorthand notation H µν = H µν (r h , x µ ). These are the essential ingredients required to derive the approximate conservation law of interest to this section.
c. Intermediate region: We proceed by matching the two sets of solutions in the intermediate region.
The matching conditions imply that the coefficient of the terms that scale as ln a have to match. In particular, Similarly, we find the following relations for the matching of the finite part of the solutions: 1 , (4.67a) 2 , (4.67b)

B
(1) 1 . Henceforth, for simplicity, we will also turn off the sources for bothπ µν and δĴ µν . The first set of relations that follows from regularity conditions implies that where . (4.69) . (4.71) d. Summary: We now provide the physical interpretation of the above holographic results.
For simplicity, we will focus on the regime h 0 T 2 , where we will see how quasihydrodynamics emerges. The analysis is similar for more generic systems. In this regime, |Q 1 | 4h 2 0 ln(Mr h ), and we find that where α is a dimensionless coefficient which can be computed.
It is instructive to change variables to those used in Section II D. Namely, the 'velocity field' along with δB y = δĴ ty and δE x = δĴ yz ; recall the definition of e r in (4.28). Here, χ denotes the susceptibility of the transverse momentum and is proportional to (ε + p) for the system described in Section II D. Eq. (4.42) then implies that χ∂ t δu y + ∂ z π zy = 0, (4.74a) The first equation is the conservation of momentum. The second equation is Faraday's law.
Finally, we turn our attention to the remaining two equations of motion. Firstly, in the limit described above, we find that (4.68) becomes Note that we have not included the term with ∂ t π zy as it is subleading relative to the constant term in π zy when ω T (the regime of validity of our analysis). Hence, π zy is not a quasihydrodynamic degree of freedom. Nevertheless, it was convenient to carry it through the calculation as if it was quasihydrodynamic. Indeed, (4.75) reduces to a first-order constitutive relation within hydrodynamics, and as h 0 → 0, we recover the classic holographic result that the shear viscosity η obeys η = s/4π. At first order in h 0 , we observe a correction to the stress tensor arising from the dynamical electric field, which would not have been included in conventional MHD. The ellipsis in (4.75) denotes higher-order corrections in h 0 , which we have neglected.
The second equation of motion follows from (4.70): where 1/τ is given by (4.29). This is precisely Ampere's law in a plasma. Since the Ohmic current must be evaluated in the rest frame of the fluid, we obtain a velocity-dependent term in (4.76), exactly analogous to (2.32c). Again, · · · denotes corrections in h 0 which we have neglected for simplicity.
Combining (4.74), (4.75) and (4.76) we thus arrived at a holographically derived theory of quasihydrodynamic MHD with dynamical electric fields in the transverse Alfvén channel, which is exactly analogous to our discussion in Section II D. In holography, we may carry out the computation to higher orders in h 0 in principle, although we will not systematically discuss the higher-order contributions here. We expect that at sufficiently low T , the electric field will remain a quasihydrodynamic mode at higher orders in h 0 , but have not found a simple analytic proof of this.
As we have shown in this section, a holographic analysis of the bulk theory (4.1) with a dynamical metric and a two-form gauge field of [31,39] opens the door to a systematic derivation of not only MHD but also its various extensions that pertain to the physics of plasmas. shear: δT µν denotes first-order perturbation of T µν and in the sound channel, and Π zz is defined as Π zz ≡ δT zz − δp, where δp is the perturbation of the pressure.
To show how the MIS equations arise from a dual gravitational description, we focus on the Einstein-Gauss-Bonnet theory [52,53,102,103], which is a useful holographic toy model for analyzing thermal field theories at intermediate coupling strengths. The crucial feature of the theory, and of other higher-derivative theories, is the emergence of quasinormal modes with a purely relaxing, imaginary dispersion relation ω(k) [51,53,104]. As was demonstrated in [51,53], these modes reproduce many expected spectral properties of thermal theories in transition from strong, towards weak coupling. They exhibit the emergence of quasiparticle-like excitations in the spectrum, formation of branch cuts [51,53,105] and the quasiparticle transport peak in the spectral function at k = 0 [106]. Moreover, the existence of purely relaxing quasinormal modes leads to the breakdown of hydrodynamics, formally defined at the critical momentum k c at which the collision occurs, given some λ GB . In the shear channel, this is a results of the collision of diffusive and gapped poles. The breakdown of hydrodynamics also leads to a longer hydrodynamization time [18,107], relevant in heavy-ion collisions [17,20,108], and a longer isotropization time [109,110].
The action of the Einstein-Gauss-Bonnet theory in five dimensions is where the negative cosmological constant is Λ = −6/L 2 . We set the Anti-de Sitter (AdS) radius to L = 1. An especially useful feature of this theory lies in the fact that its equations of motion involve only first and second derivatives. Thus, in principle, the coupling λ GB ∈ (−∞, 1/4] can be treated non-perturbatively. The background solution for this theory can be written in the form where The symbol r h denotes the radial position of the horizon. The AdS boundary is at r = 0. We set N 2 GB f (0) = 1 to ensure that the boundary speed of light equals to one. Moreover, the Hawking temperature is given by T = N GB /(πr h ). For further details regarding this theory, see [53].
Increasing the negative value of the higher-derivative coupling constant λ GB , i.e. increasing −λ GB , can be though of as tuning the dual coupling constant away from infinity. In the regime of large −λ GB , the relaxational quasinormal mode comes into the regime of |ω|/T 1 and has a non-trivial interplay with the hydrodynamic modes [51,53]. Its gap can be defined as .
This is the mode that plays the role of Π µν in the MIS theory. As we show below, its associated relaxation time τ can be found analytically when τ (λ GB )T 1 . Moreover, we will also find the speed of propagation of modes at large k, v, which was discussed in Sections II B and II C, cf. Eqs.
(2.15) and (2.29). In particular, v 2 = D/τ . The results derived from the action (5.2) are where ψ 1 and ψ 2 are function of the background metric: Performing the integral ψ 2 explicitly, we obtain the relaxation time where γ GB = √ 1 − 4λ GB . The expression agrees precisely with the (zero-momentum) frequency of the gap of the non-hydrodynamic modes in all three channels of Einstein-Gauss-Bonnet gravity obtained from the analytic quasinormal mode calculations of Refs. [51,53]: (5.9) In particular, see Eq. (2.44) of Ref. [53]. Importantly, we note that unlike in MIS theory, cf. Sec. II C, the relaxation time τ that is derived from the gravity bulk, which encodes all higher-order corrections to hydrodynamics, does not equal the second-order transport coefficient τ Π [53,102,103]: However, in the large −λ GB limit, as λ GB → −∞, Finally, the speed v is given by (5.12) Note that the reason that v is superluminal stems from the fact that it is computed from taking the limit of k → ∞-the limit in which the Einstein-Gauss-Bonnet theory suffers from various known UV problems and instabilities [53,[111][112][113]. The theory should always be thought of as one with a UV cut-off on ω and k for any non-zero λ GB . For this reason, acausal UV behavior is irrelevant for the ω T and k T regime of interest to this work in which the theory of quasihydrodynamics is applicable.
The remaining of this section is devoted to demonstrating how the equations (5.1) emerge from the matching conditions in the bulk.

A. Scalar channel
We begin by deriving the scalar channel equation (5.1a). To do this, we study the decoupled δg xy fluctuation of the black brane metric (5.3).
δ(ds 2 ) = δg xy (r, t, z)dxdy ≡ 2 r 2 h xy (r, t, z)dxdy. (5.13) The bulk equation of motion for this mode, in the Fourier basis h xy (r, t, z) ∼ h xy (r)e −iωt+ikz , can be written as where the coefficient A is The matching procedure then proceeds in the manner outlined in Section III.
b. Inner region: The inner, near-horizon region is defined by r h − r 1/T . There, we make the following ansatz for h xy : where H (i) µν (r, x µ ) are regular function at the horizon. We then substitute the above ansatz into the equation of motion (5.14) and demand that the coefficients of all divergent pieces vanish in the near-horizon expansion. For the coefficient of f (r) −1 , we find that

B. Shear channel
In this shear channel, we turn on the metric fluctuations where we have used the radial gauge h rµ = 0. The relevant set of coupled dynamical equations of motions is and a first-order constraint equation where the coefficient A was defined in (5.15). As before, it is the first-order constraint (5.26) that implies the conservation of transverse momentum T ti . This is not the case for the approximate conservation of δT iz , which instead arises from the matching condition in the bulk. Note that the outer, the inner and the intermediate matching regions are defined in the same way as in the scalar channel.

VI. CONCLUSION
In this paper, we have developed a general framework for discussing linearized hydrodynamic theories with additional approximately conserved currents, which we called quasihydrodynamics.
As shown, this framework can be used to understand a large number of phenomenological theories discussed in previous literature, including our main two examples: magnetohydrodynamics coupled to dynamical electric fields (which thus includes dynamical photons), and the relativistic Müller-Israel-Stewart theory. Since a generic feature of such theories is the presence of not only massless hydrodynamic modes, but also of "massive" long-lived modes, a systematic constructions of effective quasihydrodynamic theories is a formidable task, as the formal gradient expansion is not directly applicable. An even more difficult task is a derivation of such effective theories from their underlying quantum microscopic description-for some recent progress in the hydrodynamic setting (without approximately conserved quantities) see [114].
With a view towards a long-term goal of building systematic quasihydrodynamic theories, we studied the emergence of quasihydrodynamic theories in strongly coupled theories with holographic duals. As we have shown, holography can be used as a tool to systematically derive a low-energy description of systems with long-lived excitations. In particular, we have developed an explicit holographic algorithm for analytically computing the linearized quasihydrodynamic equations in a generic system, extending and simplifying earlier developments of [38,40]. We first showed in detail how to carry out this procedure to unambiguosly demonstrate the existence of dynamical photons in a holographic dual of magnetohydrodynamics [31]. This explicitly confirms the claims made previously in [31,39] that the holographic dual of a theory with a one-form symmetry encodes dynamical electromagnetism on the boundary. Moreover, it provides a systematic method for deriving extensions of MHD, which we anticipate will find a wealth of applications in plasma physics. In our second example, we uncovered the equations of MIS theory from a holographic higher-derivative Einstein-Gauss-Bonnet gravity theory, which was previously shown to be able to exhibit long-lived massive excitations within the low-energy (hydrodynamic) regime [51,53]. Beyond the ability to derive dynamical equations of motion and effective field content of quasihydrodynamic theories, we expect that our methods will be useful also for analytically recovering hydrodynamic dispersion relations in the holographic duals of ordinary fluids, superfluids, solids and more. The abstract nature of the method allows one to carry out the entire calculation in terms of gravitational background fields until the very end, be they analytically or numerically known. One is at most required to perform a set of simple final numerical integrals to obtain those dispersion relations.