Holographic entanglement entropy for relativistic hydrodynamic flows

We study the behaviour of holographic entanglement entropy (HEE) in near equilibrium thermal states which are macroscopically described by conformal relativistic hydrodynamic flows dual to dynamical black brane geometries. We compute HEE for strip-shaped subsystems in boundary dimensions d=2,3,4, which provides us with general qualitative inferences on the interplay between fluid flows and entanglement dynamics. At first, we consider the zeroth order in hydrodynamic derivative expansion, holographically described by stationary boosted black branes. Working non-perturbatively in fluid velocity, we find that, as the fluid velocity approaches its relativistic upper limit, the UV regulated HEE exhibits a divergence, at arbitrary temperature. Also, the holographic mutual information between two relatively close subsystems vanishes at some critical fluid velocity and remains zero beyond it. We then compute HEE in an excited state of the fluid in the presence of a sound mode, considering first order in derivative expansion. As a simplified setup, we first work with non-dissipative dynamics in $d=2$, where the time evolution of HEE is studied in the presence of the sound mode and a propagating pressure pulse. In d= 4, we find that dissipative sound modes produce an additional dynamical UV divergence which is subleading compared to `area law divergence'. No such divergence is observed for dissipative sound mode in d=3.


Introduction
The translationally invariant macroscopic thermal dynamics of a quantum field theoretic system is believed to be captured by a set of effective hydrodynamic degrees of freedom. The theory governing these hydrodynamic variables is primarily based on the underlying symmetries, which makes its applicability significantly universal. In the hydrodynamic description, the relevant microscopic information is packaged into a few parameters -the transport coefficients. The fluid variables capture the near equilibrium dynamics very efficiently, while the corresponding microscopic dynamics is extremely complex. The intricate interplay between the macro and micro dynamics usually provides us with enourmous insights into a physical system, but usually is a difficult question to address. The difficulty mainly arises due to the intractability of the microscopic degrees of freedom in such a situation. The fluid-gravity correspondence [1] provides us with a very special setting where a complete picture is provided and the transport coefficients are determined from the dynamics of a dual black brane. In this picture, the various length scales of the hydrodynamic system are encoded in a holographic direction.
The fluid-gravity correspondence maps the dynamics of a d dimensional conformal relativistic fluid to that of a AdS (d+1) black brane. The microscopic description of this holographic fluid is provided by a strongly coupled conformal gauge theory. The black brane has a complete knowledge of this strongly coupled microscopics. This is why, all the hydrodynamic transport coefficients can be computed from the black brane geometry, in a regime where a direct field theoretic calculation is intractable. Similarly, the holographic geometry can be used to extract other microscopic information about the near-equilibrium thermal states. Particularly, the entanglement structure of these underlying quantum field theoretic states, which are evolving towards equilibrium, is significantly interesting. An excellent measure of this entanglement is provided by the holographic entanglement entropy (HEE) [2][3][4][5][6][7]. In this paper, our objective is to apply the holographic dictionary to study the dynamics of HEE for the states whose macrodynamics corresponds to specific hydrodynamic flows.
In order to compute entanglement entropy in a continuum field theoretic set up, we first choose a subsystem A, usually taken to be a space-like sub-region on a chosen timeslice. The entanglement entropy of a given state, is then given by the Von Neumann entropy of the reduced density matrix obtained by integrating out the degrees of freedom outside A. The knowledge of this entanglement entropy for various subsystems provides us with enourmous quantum information about the underlying state. Holographically, this entanglement entropy is computed by the celebrated Ryu-Takayanagi formula [2,3], which was generalized to a more general covariant form in [4]. In this prescription, the HEE for A is given by where A Γ is the area of the co-dimension 2 surface homologous to A which has been extremized on the bulk geometry dual to the the field theoretic boundary state. The subscript Γ refers to this extremal surface. In our work here, using the fluid-gravity map, we will explicitly write down the geometries dual to specific hydrodynamic configurations. Subsequently, we shall compute the area of the extremal surface in these geometries to obtain HEE (S A ) for the hydrodynamic states. For the subsystem A, we choose a stripshaped sub-region on a canonical time slice of the AdS boundary. Our choice of subsystem is time-independent, while the underlying state is a fluctuating fluid configuration (see fig.1).
In this paper, we shall primarily focus on two varieties of fluid configurations. In the first case, we shall consider steady state flows where the fluid is moving with a constant velocity (see section 2). The computation of HEE for such thermal systems has been attempted by several authors in the recent past [8][9][10][11][12][13][14][15][16]. However, in all these earlier works, some simplifying assumption have been made, and therefore, a thorough analysis, particularly in the regime where the fluid velocity approaches its relativistic upper bound, seems to be missing so far. For this steady state case, apart from HEE, we also compute the holographic mutual information between two non-overlapping strip-shaped subsystems. Among various other results, we find that there is a critical fluid velocity beyond which mutual information vanishes due to an exchange of dominance of extremal surfaces. This is similar to the 'disentangling transition' previously observed in [8,17], where the abrupt vanishing of mutual information was observed as the separation between the two subsystems in question exceeded a certain critical value 1 . In this paper, we observe such a 'disentangling transition' taking place as a function of the fluid velocity. This we believe, is a new and important result, which was not possible to obtain without the non-perturbative approach of our paper.
In another case, we consider fluid configurations where a propagating sound wave have been set up with some initial conditions. We then proceed to study the change in HEE produced due to such fluid fluctuations in the linearized approximation, assuming a small amplitude of the excitation (see section 3). Since, sound waves constitute one of the basic fluid excitations, it is interesting to understand the HEE dynamics associated with it. In d = 2 the sound mode is non-dissipative, while in higher dimensions these fluctuations dissipate away to equilibrium. In the latter case, we aim to capture the time-evolution of HEE as the fluid settles back to equilibrium. In d = 2 where the sound mode is nondissipative, we can superpose them to constitute a pressure or energy density pulse. In this case, we study how HEE evolves with time as the pulse passes across the subsystem. Our investigation here is reminiscent of [18] where HEE has been computed in the AdS plane wave geometry with an energy flux flowing through the subsystem. We are also inspired by the body of work on the study of HEE following quantum quenches (see for instance [19][20][21]). In fact, we may consider our fluid dynamical set-up to be a thermal state which has been quenched by a macroscopic observable associated with the energy-momentum tensor, and we study the subsequent late time evolution of HEE.
In our work here, we employ a mixture of analytical and numerical techniques. In particular, most of our results in d = 2 is analytical. We use this analytical result to optimize our numerical procedure, which are then applied to d = 3, 4. The qualitative inferences which we draw from our results in d = 2, 3, 4 may be straightforwardly generalized to higher dimensions.

Summary of results
Before proceeding to the details, let us summarize the important new results of our paper below. In this paper we investigate the behaviour of entanglement entropy in a class of conformal fluid states, using holographic techniques, for a strip like subsystem.
• We compute HEE for boosted black holes in d + 1 (d = 2, 3, 4), viewed as holographic geometries dual to fluids moving with a constant velocity. Our subsystem lies on a 'canonical time slice' and its size is fixed as we vary the fluid parameters. Our calculation is performed without making any assumption on temperature or velocity parameter. All earlier studies of this system had made some simplifying assumption about these parameters and lifting those assumptions is a new aspect of our study. In particular, we curiously observe that the regulated HEE diverges as the fluid velocity approaches its relativisitc upper bound (v → 1) (see fig.4 and fig.6). See section 2 for more details.
• For the constant velocity fluid flows, we also compute holographic mutual information. We find that this quantity undergoes an abrupt 'disentangling' transition for some critical velocity for any separation between the subsystems (see fig.5 and fig.7). This critical velocity tends to 1 as the separation between the subsystem tends to zero. The existence of this critical velocity is a new and interesting result. See section 2.2.2 and section 2.3 for more details.
• For holographic fluids in d=2, we have studied the time dependence of HEE as a non-dissipative propagating pressure pulse passes through the demarcated subsystem.
The geometric dual to this pressure (or energy) pulse is an exact solution to the bulk Einstein equations. We find that HEE has a curious dip when the pulse enters the subsystem, followed by an expected rise and plateau while the pulse resides inside the subsystem (see fig.10). We report these results in section 3.2.
• For holographic fluids in d = 3, 4, we study the behaviour of HEE in the presence of a dissipative sound mode. In this case HEE has the usual expected decay closely following the damping of the sound mode. However, surprisingly, we find that in d=4, there is an additional UV divergence which is subdominant compared to usual area law (see (3.44)). We present the details of these results in section 3.3.
2 Entanglement in steady state relativistic fluid flows

General Procedure
Before proceeding with the hydrodynamic excitations, we would like to consider the simple steady state case, when the fluid is moving with a constant velocity. More precisely, we consider a relativistic conformal fluid, with a constant velocity, in a flat Minkowski background with η µν = Diag{−1, 1, 1, . . . }.
Holographically, this fluid configuration is described by a bulk metric of a boosted black brane, which in d + 1 dimensional Schwarzschild coordinates 2 , takes the form 3 Here, 0 ≤ v ≤ 1 is the boost parameter which is identified with the fluid velocity. Also, the AdS boundary is reached in the limit z → 0. We denote the horizon radius with ζ, which is related to the temperature of the black brane as This is also identified with the steady state temperature of the fluid. The equations of hydrodynamics in steady state is provided by the conservation of the energy-momentum tensor, which has the following constitutive relations where P µν = η µν +u µ u ν is the projector orthogonal to the fluid velocity u µ , E(T ) and P(T ) are respectively the energy and pressure densities, while O(∂) represents higher derivative corrections to this perfect fluid form. In steady state, the temperature and fluid velocity are constants over space-time, implying that all the higher derivative terms vanish. Also, for such a steady state it is immediately clear that the fluid variables satisfy the conservation equations ∂ µ T µν = 0. The form of the energy momentum tensor (2.3) follows from the boosted AdS black brane (2.1), via the AdS-CFT dictionary. The duality also determines the precise functional forms In this section, we use the Schwarzschild coordinates since this has been predominantly used in the literature on HEE. The derivative expansion in the fluid-gravity correspondence [1] necessitate the use of more regular Eddington-Finklestein (EF) coordinates, so as to formulate a well defined derivative expansion. Therefore, while discussing hydrodynamics modes with first order derivative corrections (in section 3) we will use EF coordinates and in appendix A we have explicated the transformation between the two coordinates. 3 We shall set the AdS radius R AdS = 1 throughout this paper.
On the boundary of the boosted black brane, we choose a strip like subsystem (A) on a canonical time-slice 4 and we wish to quantify the entanglement of A with the region outside the strip. The underlying fluid is moving with a constant velocity in some given direction (see fig.1 for a schematic of the situation).
For zero boost (v = 0), the HEE for black branes in arbitrary dimensions have been computed in [8][9][10]. The extension of this analysis to boosted black branes have been attempted in several papers earlier [11][12][13][14][15][16]. However, in all these papers some simplifying approximations have been employed, such as small boost in [12][13][14], or small temperature [11]. To our knowledge, a complete computation of HEE for arbitrary boost and arbitrary temperature seems to be missing in the literature. This subsection aims to perform a thorough analysis of this situation, particularly focusing on the situation when the fluid velocity approaches its relativistic upper bound (v → 1) at an arbitrary background temperature.
As pointed out in [11], for a non-perturbative HEE computation in this scenario, the extremal surface is not expected to lie on a bulk time-slice. Therefore, even if we are using Schwarzschild coordinates (2.1), for the boosted case, we must employ the covariant protocol of [4].
Here, we shall consider a strip-shaped subsystem in a d dimensional space-time. The symmetries of such a strip-like subsystem simplify some of the calculations significantly. From the boundary point of view, there are two special spatial directions in this set-up. The short edge of the strip is directed along some spatial direction, while the boosted fluid may be flowing in any other independent direction. The result of HEE is expected to depend on the angle between these two vectors. But again for the sake of simplicity, in this paper, we will assume these directions to be collinear, i.e. the short edge of the strip is extended along the same direction as the boost. With this in mind, the embedding ansatz for the extremal surface, which is a co-dimension 2 surface homologous to our strip like subsystem on the boundary (see fig.1), is taken to be In other words, the world-volume of our co-dimension 2 extremal surface is parameterized by φ andỹ i . The φ-coordinate on the world-volume takes values between −π/2 ≤ φ ≤ π/2, and z (< ζ) is the turning point of the extremal surface, which will be determined in terms of the width of the strip on the boundary. Further, on the boundary, to regulate the infinite volume arising due to the long sides of the strip, we will consider theỹ i coordinates to be bounded Now, since the boost is along the x-direction (along the short-edge of the strip) (2.1), we retain translation invariance along theỹ i -directions. Thus, we need to solve for two 4 Our boundary time-coordinate is chosen to be the boundary limit of bulk Eddington-Finkelstein-time or Schwarzschild-time, both being identical on the boundary of AdS. Constant time-slices on the boundary, with the time coordinate defined in this way, will be referred to as 'canonical time-slice'. Note that, since the Eddington-Finkelstein coordinates are indispensable for the fluid-gravity bulk metric [1], the boundary time coordinate defined in this way, is somewhat special for the holographic hydrodynamics.
functions of one variable, x(φ) and t(φ). The area functional of Γ, to be extremized, is given by (2.6) where x and t denotes derivatives of these functions with respect to φ. Here, the factor of L d−2 arises due to the integrals over theỹ i coordinates on the world-volume of the extremal surface. Also, in (2.6) g = g(z cos φ), and f = f (z cos φ) and their functional forms are given in (2.1). The equation of motion for x(φ) and t(φ) which follows from the extremization of A Γ is given by We now have to solve these equations (2.7), with the following boundary conditions Since, we are dealing with a stationary scenario, we expect the HEE to be independent of time. Hence without any loss of generality we have chosen the boundary time slice at t b = 0. The equations (2.7), (2.8) is difficult to solve analytically. Therefore, we have to resort to numerical analysis for computing the area of the extremal surface.
In the absence of boost (v = 0), t(φ) = constant is a solution to the equations (2.7). This implies that, in the unboosted case, if we use Schwarzschild coordinates, the extremal surface always lies on a constant time-slice throughout the bulk. This, however, is not the case in the presence of boost, even if we work in Schwarzschild coordinates.

Extremal surface reparameterization: φ → z
It is often practically useful to consider the embedding of the extremal surface such that the AdS-radial coordiate z is identified with one of the world-volume coordinates. With such a choice, the embedding equations for the strip-subsystem would read with the equations of motion (2.7) being minimally modified to ∂ ∂z (2.10) Now x and t denote derivatives with respect to z, while f = f (z) and g = g(z) as given in (2.1). In this case, the boundary conditions (2.8), are implemented with a minor difference. The solution to (2.10) has two distinct branches (x ± (z), t ± (z)) which are valid in the domain z ∈ [z , 0], with z being the turning point. The boundary conditions for the two roots are implemented in the following way The two branches are smoothly stiched together at the turning point z (see fig.3 for the extremal line in d = 2). The area functional, after this reparameterization, is given by (2.12) This area functional evaluated on the extremal surface gives the HEE by the formula (1.1). We shall regulate the UV divergence by removing the HEE corresponding to pure AdS The ∆S E , as defined above, will be useful for quoting our numerical results.

Fluid flows in d = 2 dual to boosted BTZ black brane
In 1 + 1D boundary, an exact expression of HEE for the steady state fluid flow dual to boosted BTZ black brane, may be obtained by a simple argument using a result available in [22]. After presenting this argument in section 2.2.1, we proceed directly to solve the equations (2.7) numerically and verify our argument. Further, in d = 2 we can also use our analytical result to compute mutual information for configurations with constant fluid velocity, which we present in section 2.2.2. Our numerical procedure in d = 2, having been verified by an exact analytical result, is then straightforwardly carried over to higher dimensions in section 2.3.

Holographic Entanglement entropy Analytical exact result
Let us first consider the unboosted BTZ black brane in Schwarzschild coordinates where, f (z) is given by the same function in (2.1) for d = 2. This metric is dual to a static fluid configuration in d = 2 with fluid velocity u µ = {−1, 0} and temperature given by (2.2). Now, we wish to define our boundary subsystem (A) as the region between two space like points, say P and Q. Without the loss of generality, we can consider the boundary points to be P ≡ (0, 0) and Q ≡ (τ, ξ). The HEE for this situation was worked out in [4,22] and the answer is given by where uv is the ultraviolet cut-off 5 . Let us boost this set-up with velocity v, so that the fluid velocity becomes u µ = γ(−1, v). Under this boost, the bulk BTZ metric (2.14) goes over to (2.1) (for d = 2). Besides affecting the underlying fluid state, the boost also transforms the end points of our subsystem P and Q. Consequently, our subsystem A will now lie on a different time-slice compared to the time-slice before boosting. We wish to compute HEE for the boosted BTZ when our subsystem A lies on a canonical time-slice (see footnote 4).

Figure 2:
Schematic representation of the argument leading to the computation of HEE for the steady state of a 1 + 1D conformal fluid, moving with a constant velocity. We start with the static fluid and find HEE for a subsystem defined by the points P and Q as shown in the left figure. This situation is boosted so that the system, now defined by the pointsP andQ, returns back on to a canonical time slice as shown in the right figure (see footnote 4). The boosting ensures that the background fluid is moving with a constant velocity, equal to the value of the boost parameter. This procedure enables us to write down an exact expression for HEE in this case, as given by (2.19).
The boosted coordinates would be related to the old coordinates by the transformatioñ where P µν = η µν + u µ u ν is the projector orthogonal to u µ . Now, prior to boosting if the end points of A are P and Q, then after boosting they must transform toP = (0, 0) and Q = (0, ), which will ensure thatP andQ both lie on the canonical time-slice. Thus, we get back to a canonical time-slice only after boosting. This would ensure that the underlying state is that of a fluid moving with a constant velocity, while we are measuring entanglement on a canonical time slice (see fig.2). To accomplish this we choose (τ, ξ) so where {tQ,xQ} are the coordinates of the pointQ which is the boost transformed version of point Q = {τ, ξ}. This fixes τ and ξ to following Plugging this back into (2.15), we obtain This is the final expression of HEE for the 1 + 1D hydrodynamics state, where the fluid is moving with a constant velocity. In this next subsection, we shall reproduce this numerically directly from the equations (2.7), (2.8).
Note that the result for this boosted case (2.19) is significantly different from the unboosted thermal state for the same subsystem lying on a canonical time slice. The later is given by setting α = 1 in (2.19) Although trivial, it is also interesting to apply the logic of this subsection for the vacuum state. The boundary vacuum state is dual to pure AdS, for which the analogue of the formula (2.15) applied to an identical subsystem A with boundary points P and Q, is given by which is identical to HEE for unboosted vacuum. We find this match simply because the vacuum is boost invariant. On the other hand, the boundary state that corresponds to a black brane is dual to static fluid (deconfined plasma). Such a state is not boost invariant. Viewed as hydrodynamic states, the fluid at rest is different compared to the fluid moving with a constant velocity and their structure of entanglement is also distinct.
It is instructive to explicate the large velocity limit of the exact analytical formula (2.19).
If we first perform a small temperature (large β) expansion of (2.19) we get This result is exact in v and matches with [11] for d = 2. Note that, from (2.23) it appears that the regulated entanglement entropy ∆S A := S A − S (vac) A diverges as γ 2 when v → 1. This conclusion, however, is slightly misleading, as far as generic large velocity (v → 1) behaviour is concerned. This is because, if we work out the v → 1 limit of (2.19) at arbitrary temperature, we find which shows that the leading order divergence in ∆S goes as γ for large velocity. Note that the leading order behaviour in (2.24) is same as that obtained in the large temperature limit, although no such assumption was made to arrive at (2.24). This shows that the small temperature result reported for arbitrary dimensions in [11], although valid for arbitrary boost, may not capture the generic large boost behaviour correctly 6 .

Numerical verification
We can adapt the procedure outlined in section 2.1 to d=2, and find a numerical solution to the equations for the extremal surface (2.7). In fig.3, we provide a representative plot of the extremal surface (which is a line for d = 2). In our solutions we clearly observe that t(φ) has a non-trivial dependence on the world-volume coordinate φ, implying that the extremal surface, even in the Schwarzschild coordinates (2.1), does not lie on the constant time-slice extended into the bulk. With this extremal surface determined numerically, we can directly evaluate the regulated HEE for the boosted fluid (bf) Here, the pre-factor c/6 = R AdS /4G N has been introduced as a normalization which is convenient for presenting numerical results 3 . Our numerical results are exact in boost parameter v and valid for arbitrary temperature. In fig.4, we numerically plot ∆S E as a function of v at fixed T and we find an excellent match with the analytical answer in (2.19).
In particular, the divergence structure (2.24) is correctly reproduced by the numerics. The fixed temperature for this plot is T = 1/2π and subsystem size is = 2. Note that since we have regulated our HEE with the vacuum HEE, the v → 0 limit is non-zero and provides the regulated HEE for the static fluid at fixed temperature. The answer diverges as v → 1 following (2.24).

Holographic Mutual information
Another useful measure of quantum entanglement for a given state is quantum mutual information. For non-overlapping subsystems A and B, mutual information can be expressed in terms of entanglement entropy as follows where S A represents the entanglement entropy between the subsystem A and its compliment. Given the property of strong subadditivity of entanglement entropy, this quantity is guaranteed to be positive (I(A, B) ≥ 0). In addition, it is UV finite. With the help of (2.26) the holographic prescription for entanglement entropy can be adapted to compute quantum mutual information (see [5,8,26]).
For the strip-shaped subsystems, let us consider both A and B to have the same size , while the distance between the two is denoted by x. For the fluid flows with a constant velocity v, I(A, B) depends on three dimensionless parameters. Apart from the fluid velocity v, the other two parameters are the ones which may be constituted with dimensionless ratios between , x and temperature T . Since, we have observed a peculiar diverging behaviour of entanglement entropy as the fluid velocity approaches its relativistic upper bound, we would now like to understand the behaviour of mutual information for large fluid velocities. For d = 2, the existence of the analytical result (2.19) immediately implies an analytical formula for mutual information I(A, B) = I( , x, β, v), which is given by (2.27) where α(v) is given in (2.18) and β is the inverse temperature (2.2). When the subsystems A and B are very far away from each other I(A, B) vanishes generically. In fact, I(A, B) continuously reduces to zero at some finite critical distance (x c ) and remains zero for larger distances between A and B. This fact is well understood in the holographic context [8]. The term S A∪B in (2.26) has two extremal surfaces homologous to A ∪ B. At the critical distance, dominance flips between these two extremal surfaces. This critical distance (x c ) depends on the properties of the underlying state, for example, it is determined by the temperature for a thermal state. In [8], the functional dependence of this critical distance on temperature was computed for v = 0. Given the exact result (2.27) in d = 2, we can immediately generalize their result for a non-zero fluid velocity.
In fig.5a, we have plotted I(A, B) versus the fluid velocity v at fixed x/ and T . We observe that, for a given separation x, there exists a critical velocity v c where mutual information vanishes. The value of v c increases as we decrease x, approaching its relativistic upper bound as x → 0. Once the mutual information reaches zero, it remains so beyond v c . At v c , there is an exchange of dominance between the two extremal surfaces for S A∪B . This is exactly the same phenomenon which gives rise to a critical distance x c between A and B, as we have described above. Here, we find that the deformation of the extremal surfaces, produced due to the fluid velocity (see fig.3), forces this exchange of dominance between the extremal surfaces, even when the two subsystems are close to each other (x < x c ). In fig.5b we demarcate the region in the {x/ , v} plane, where mutual information is non-zero for a fixed temperature. As expected from the results of [8], this region shrinks with the increase in temperature.

Higher dimensional fluids
In higher dimensions, we do not have any analytical results. But, it is straightforward to generalize the verified numerical procedure outlined for d = 2 in section 2.2.1. The general procedure has already been outlined in section 2.1 for general dimensions. Therefore, skipping the details, we directly present the results of our numerical computations in this section.
The variation of the regulated HEE (2.25) versus the fluid velocity v, keeping T fixed, has been plotted in fig.6a for d = 3, 4. We observe that just like in d = 2, the   regulated HEE diverges as v → 1. In fig.6b, we demonstrate that the divergence occurs as ∆S E ∼ (1 − v) −1/2 as v → 1. Again, this behaviour is identical to that observed for d = 2. With this observation, it is extremely tempting to conjecture that, in this setup, the divergence in ∆S E occurs due to a factor of γ. Perhaps this is universal across dimensions for a strip-shaped subsystems, with the fluid moving parallel to its short edge. It would be extremely interesting to understand the nature of this divergence for more general subsystem geometries.
Once we are able to compute the entanglement entropy with reliable numerics, it is straightforward to compute the mutual information (2.26) in higher dimensions. In fig.7, we plot the mutual information versus the fluid velocity in d = 3. Even here, we notice a behaviour which is qualitatively identical to that in d = 2 (see fig.5a). In particular, we observe that there is a critical velocity beyond which the mutual information vanishes. Again, this feature is also perhaps true across all dimensions.

Outline of the general procedure
After considering the steady state fluid configurations with a constant velocity, in this section, we proceed to analyse the entanglement structure of fluid states which are slightly out of equilibrium. More specifically, we shall focus on the propagating sound mode, which is ubiquitous in all systems admitting a macroscopic fluid description. Our set-up here, will again be holographic, which technically limits our observations to conformal relativistic fluids. However, our qualitative observations in this section would perhaps be more universal than the holographic context.
The sound mode is a fluctuating small amplitude perturbation over a stationary background 7 . We will compute the HEE for the sound mode using the fluid-gravity correspondence [1]. We shall now provide a brief overview of the general procedure adopted for this HEE computation.
The dual black brane solution providing a holographic description of the fluctuating states of a relativistic conformal fluid in d = 4, with a controlled derivative expansion, was 7 For simplicity, throughout our analysis in this section, we shall consider the sound mode over a static background. Our analysis can be straightforwardly generalized to the case with a background fluid velocity (see section 4 for a discussion about this point). first provided in [1], which was subsequently generalized to other dimensions in [27][28][29]. Using these constructions, it is simple to read off the black brane geometries dual to the sound modes. Since these bulk geometries are constructed in the derivative expansion, they are valid when the wavelength (λ = 2π/k) of the sound mode is large compared to the temperature scale (λT 1). Beside this small parameter k/T , the amplitude of the sound mode (say, ) also must be small 8 , since it is a linearized solution to the fluid equations. We shall consider the amplitude ( ) as the smallest of the two parameters.
Due to the linearized nature of the sound mode, the dual geometry of the fluctuating black brane has the following structure Here g (1) ab is the -order correction to the background metric g ab has a further expansion in powers of k due to the derivative expansion corrections to the bulk metric.
Similar to the previous section, we would like to compute HEE for a strip-like subsystem over the fluid configuration with a propagating sound mode (see fig.1). The short-edge of the strip is of length , and again, to keep the computations simple, we shall restrict ourselves to the scenario where the sound mode is propagating along the direction of the short edge 9 . Apart from k/T and , the dimension of the subsystem introduces a third independent dimensionless parameter in the problem, which is given by T . While there are no a priori restrictions on this parameter, in the holographic set-up, if we are to capture the thermal effects adequately, we must ensure that T ∼ O(1) or greater.
The HEE computation is accomplished by computing the area of the extremal surface as prescribed in (1.1) [4]. Unlike the stationary scenario discussed in section 2, we now expect a non-trivial time dependence of HEE. If we parameterize our co-dimension 2 extremal surface (Γ) with the coordinates {ξ}, the induced metric on Γ takes the form However, due to the decomposition (3.1), the induced metric also decomposes as ab .
(3.4) 8 As we will see in the next subsection, is the dimensionless parameter which may be taken as the ratio δT /T0, where δT is the amplitude of temperature fluctuations induced by the sound mode, while T0 is the background temperature. 9 This can be easily generalized to the case where there is an angle between the direction of the short edge and that of sound propagation. Naively, it may be expected that the component of the sound mode along the long-edge of the strip would not contribute to ∆SE.

The area of the co-dimension 2 surface Γ is immediately given by the induced metric
where G is the determinant of the induced metric (3.2). However, due to (3.3) this area also has an expansion in amplitude Here A (0) Γ has the contribution from the background static geometry and constitutes the leading order piece in HEE. It is given by By the prescription (1.1), this area when evaluated over the extremal surface provides the entanglement of the strip subsystem for the fluid in static equilibrium. On the other hand the corrections due to the sound mode is captured by A Γ , which is given by Note that, at first order in , the functional in (3.8) is to be evaluated over the same surface Γ that was obtained by extremizing A (0) in (3.7). The change in HEE produced due to the fluid sound mode is therefore given by As in the previous section (see (2.13)), for our numerical plots, it will be convenient to work with the quantity ∆S E = A

Sound mode in d = 2
In d = 2, conformal relativistic fluid dynamics is trivial and there are no non-trivial transport coefficients, in the uncharged sector 10 . Therefore, the energy-momentum tensor is always of the perfect fluid form The key difference of (3.10) with that of in higher dimensions (2.3), is that, (3.10) does not have any higher derivative corrections, even if the fluid variables (T and u µ ) have arbitrary space-time dependence. Here, the energy density E and pressure P must be equal (see (2.4)) so as to ensure tracelessness of the energy-momentum tensor -a requirement of conformality. All solutions to this system are non-dissipative. The two fluid variables -T and the spatial component of normalized u µ -are to be determined by the fluid equations which are given by the conservation of the energy-momentum tensor Here, (3.11) gives us two independent equations, providing solutions for the two fluid variables. These solutions can have a non-trivial dependence on space-time coordinates, without requiring higher derivative corrections. Thus, unlike higher dimensions, the restrictions due to derivative expansion on the sound wave momentum k (i.e. k/T 1 ) does not apply for the d = 2 fluid.
For any arbitrary solution of (3.11), an explicit holographic dual metric was constructed in [28] The derivative corrections in the second line of (3.12) ensures that the metric is an exact solution to the Einstein equations, whenever the fluid parameters u µ (x) and T (x) = 1/ (2πζ(x)), as functions of the boundary coordinates x µ , solves (3.11). The explicit form of the fluid equations (3.11) is given by The metric (3.12) is essentially diffeomorphic to the BTZ black hole [28], which on the gravity side, reflects the triviality of conformal hydrodynamics in d = 2. The dynamics induced by these non-trivial diffeomorphisms serves as a toy model for hydrodynamics with a non-dissipative sound mode. Note that, (3.12) is similar to a fluid-gravity metric in higher dimensions (see (3.19), (3.31)) expanded upto first order in derivatives, but (3.12) is an exact solution. This absence of the requirement of a derivative expansion makes this setting particularly suitable to initiate a study of the HEE for the sound mode.

Linearized solutions
Let us denote the boundary coordinates with x µ = {t, x} and express the two fluid functions in the following way Then (3.13) admits the following linearized sound wave solution where ζ 0 determines the background temperature.
Exploiting linearization, it is straightforward to superimpose these sinusoidal solutions with a Gaussian weight to constitute a propagating wave-packet solution of the form where σ is the width of the pulse. Clearly, (3.16) is also a linearized solution of the fluid equations (3.13). Note that, physically this solution corresponds to a energy density or pressure pulse propagating in the positive x-direction. The amplitude of this pulse is proportional to , which in general, denotes a small positive or negative correction to their background values. In our subsequent numerical computations, we shall consider in (3.16) to be positive, corresponding to a positive correction to the background energy density. We can plug in the solutions (3.15) and (3.16) into (3.12) to obtain the corresponding bulk metric, which will be clearly of the form (3.1). We can then use the procedure outlined in section 3.1 to compute corrections to HEE corresponding to these fluid solutions.

Extremal surface in the background
In order to implement the algorithm of section 3.1, we need to find the extremal surface in the background metric (see (3.1)) This is simply the static BTZ black brane written in Eddington-Finkelstein (EF) coordinates and is obtained by plugging in the O( 0 ) terms of (3.15) and (3.16) into (3.12).
If we now consider a boundary subsystem (A) of length centered at the origin (x ∈ [− /2, /2]). Here, the extremal surface in EF coordinates anchored on and homologous to A can be written down analytically. It is simply a coordinate transformation (see appendix A) of the corresponding extremal surface in Schwarzschild coordinates [10]. The explicit form of the embedding equations for this extremal surface in EF coordinates is given by where z = ζ 0 κ, and κ = tanh 2ζ 0 .

(3.18)
Here t b denotes the boundary canonical time-slice which contains the subsystem A. The dependence of HEE on t b provides its time evolution.

Numerical results for HEE
In order to find ∆S E , we must evaluate the functional (3.8) over the extremal surface (3.18). The metric correction g (1) ab in (3.1), and hence the corrections to the induced metric G (1) αβ is obtained by computing the O( ) term from (3.12) after substituting the linearized solutions (3.15) and (3.16) 11 . We will now present our results for ∆S E in terms of numerical plots.

Sound-waves:
We now compute ∆S E for the scenario where a non-dissipative sound-wave (3.15) with momentum k constantly passes over our subsystem A. Once we assume a small value of the parameter , there are two more independent dimensionless parameters in the problem -T 0 , k/T 0 , T 0 being the background temperature.
The variation of ∆S E with boundary-time (t b ) has been shown in fig.8. We find a periodic variation of HEE, which was expected due to the sinusoidal nature of the fluid configuration. The period of the oscillations are given by 2π/k and is independent of temperature. This was expected from the fact that the time-dependence only arises through the term cos kt appearing in (3.15). In the embedding equations (3.18), t goes linearly with t b , which gives rise to a term proportional to cos kt b in final linearized answer.
The amplitude of these oscillations is determined both by the temperature and k, for a fixed . From fig.8 we see that for fixed k, this amplitude increases with increase in temperature. Recall that there exists no restriction on the value of k here, due to the lack of any derivative expansion for the d = 2 conformal fluid. So, we can plot this amplitude for arbitrary k at a fixed temperature. In fig.9a, we study the variation of ∆S E with respect to k on a fixed boundary-time slice t b = 0. It is very interesting to note that for large values of k, ∆S E behave sinusoidally with k. In particular, there exists special values of k for which ∆S E vanishes. The distance Ω between the first two zeros of ∆S E (t b = 0), in 11 Note that the expression for G (1) αβ is significantly long. Therefore, to avoid clutter we refrain from reporting it in this paper. It may be straightforwardly computed with any software supporting symbolic manipulation ( 'Mathematica' has been used in our work). The final answer for ∆SE is obtained after integrating the functional in (3.8), which we were only able to perform numerically.  fig.9a, as well as the period of oscillation for large k , is determined by the temperature.
In fig.9b, we have shown that the value of Ω decreases with the increase in temperature. Propagating pulse: We now consider a traveling wave-packet -a pressure or energy density pulse (3.16). We record the change in HEE as this pulse traverses the subsystem A. As expected, when the pulse is completely outside the subsystem A, travelling towards or away from it, ∆S E is zero. But, while the pulse moves across the subsystem, ∆S E exhibits an interesting time dependence, as shown in fig.10.
After fixing to a small value, there are three length scales in this problem as well. Apart from and T , there is σ which is the width of the pulse. When the pulse of positive energy density moves completely inside the subsystem A, it perhaps drags in some 'matter' within A from outside. From this naive picture we can expect that there will be an increase in ∆S E , which would eventually go down to zero, when the pulse passes out. In fig.10a we see that this expectation is almost correct, except that there is a additional feature. When the pulse enters (or exits) the subsystem, we find that there is an interesting dip in ∆S E before the expected rise (or fall).
In fig.10b, we demonstrate that if the width of the pulse σ is comparable to or larger than the size of the subsystem , there is a unique peak of ∆S E , which corresponds to the time when the maxima of the peak coincides with the central point of the subsystem. On the other hand, when σ is significantly smaller than , the ∆S E saturates to a constant value for the duration the pulse spends inside the subsystem. In both the case, the dip in ∆S E , during entry or exit, persists.
The extent of the dip (h 2 ), and the maximum value of ∆S E during the process (h 1 ), both are determined by the background temperature T 0 , for a fixed width of the pulse. In fig.11, we plot the variation of h 1 and h 2 with respect to background temperature. We note that h 2 decreases with temperature, while h 1 increases. It appears from fig.11b that h 2 would perhaps tend to a finite value as T 0 → 0. This suggests that the physical origin of this dip may be due to some non-thermal effect, which comes into play when propagating fluctuations of energy-momentum tensor tries to enter or exit the subsystem.  fig.10) as a function of temperature.

Computation of ∆S E for sound mode
The change in the area of the extremal surface Γ in (3.23) due to linearized sound wave fluctuations, which captures the corresponding change in entanglement entropy, is given by (see (3.8)) where, I(z) has the following form where, We have expanded I I and I R in z, about the asymptotic boundary, to demonstrate the following interesting feature of the integrand I(z).
Let us recall that since ∆S E captures the change in HEE compared to the thermal background value, we do not expect any leading order divergence in ∆S E , but there may be subleading divergences. Also, since we are in d=3, a subleading divergence would mean a logarithmic divergence. Clearly, in (3.29) we see that there is a term in I I which is O (1/z). This term may lead to a potential UV logarithmic divergence upon performing the integration in (3.27). However, the term is proportional to (a − b) and we know for the fluid metric to satisfy the Einstein's equation we must have a = 1 and b = 1. Thus consistency with the right fluid-gravity metric necessitates the vanishing of this potentially divergent term. We should note that this cancellation occurs only at the linear order in k. This is the order up to which our calculation is reliable since we have considered only first order in derivative expansion. It is possible that at higher order this subleading logarithmic divergence exists with some universal coefficient. As we shall see below, this magical cancellation does not happen in d = 4, where a subleading UV divergence is observed.
From the structure of (3.28) and the embedding equation (3.23), it is immediately clear that ∆S E has a straightforward dependence on t b , which factors out of the integral. In fact upon performing the integration in (3.27), we should have the following schematic form (3.30) The real and imaginary parts of ∆S E (t b ) has been plotted in fig.12a, which is clearly damped harmonic due to form (3.30). The phase of S, which we denote by Φ S , represent This variation is trivially damped harmonic as predicted by (3.30). In (b) we plot the phase of the amplitude Φ S as a function of background temperature T 0 .
a phase shift in ∆S E compared to the original sound wave perturbations (3.21). The variation of this phase with respect the T , for a fixed k, has been shown in fig.12b. We find that this phase has a weak dependence on temperature, and exhibits a growing trend towards high temperature.
to the 'area law' divergence [5], which for d = 4 is proportional to 1/ 2 uv 13 . We should emphasize that our observation is limited to O(k), and it is possible that there are other subleading divergences at higher order in the derivative expansion. The integral in ∆S E can be performed numerically only after extracting out the divergence. The integrated answer has the final schematic form (3.44) Here, the function S(k, , T ) is UV finite and has the same qualitative form as in d = 3 (see fig.12). Also, recall that we have rescaled HEE by a factor of L 2 . Hence, this UV divergent piece in the full HEE is proportional to (kL 2 / uv )e −λ t b exp i k t b √ 3 , which is subdominant compared to the area law divergence proportional to L 2 / 2 uv . Clearly, the subleading divergence decays with time and disappears in the final equilibrium state.

Discussions
In this paper, we have studied holographic entanglement entropy for two varieties of fluid configurations. At first, we have considered a stationary configuration where the fluid moves with a constant velocity. Here we find that HEE generically increases with an increase in fluid velocity when the other parameters are held fixed. Also, as the fluid approaches its relativistic upper bound v → 1, the regulated HEE diverges. From our analysis, it appears that this divergence is proportional to a factor of γ in d = 2, 3. This is perhaps a universal behaviour for strip-shaped subsystems which we consider here. Admittedly, we do not have a clear understanding of this additional divergence. However, we would like to point out that a similar divergence was also present in the earlier perturbative results of [11], although the structure of the divergence is slightly different in our non-perturbative treatment.
In this context, we also study holographic mutual information between two nonoverlapping subsystems of the same size, separated by a distance x. When the fluid velocity is zero (v = 0), it is known from previous work [8] that this mutual information vanishes at large x, while it is non-zero at small x. There is a 'phase-transition'-like effect for mutual information at some critical distance x c . For our steady state fluid flows, we find that even when the subsystems are close by (x < x c ), there exists a critical value of fluid velocity v c above which mutual information vanishes and remains zero for higher velocities. This occurs due to the flip of dominance between extremal surfaces, similar to what is observed when at v = 0, the separation between the subsystems is increased. This critical velocity exists for all values of x < x c . We find that v c increases with a decrease in the separation between the two subsystems, and v c → 1 as the separation approaches zero. This result 13 The presence of this UV divergence in ∆SE definitely puts a question mark on our perturbative approach. But, since this divergence is subdominant compared to the leading 'area law', our perturbation theory is perhaps not completely invalidated.
follows from an analytical formula in d = 2 (see (2.27)), while it has been verified numerically in d = 3 (see fig.7). Again, we think, the qualitative aspects of this observation straightforwardly generalize to higher dimensions.
Our computation of entanglement measures for the steady-state fluid flow has a few immediate generalizations. In our work, we have always assumed that the fluid velocity is along the short edge of the strip. This may be generalized to a general angle, giving rise to a richer structure of HEE. It is tempting to speculate that, as we vary this angle, there may be a phase-transition-like effect as observed in [18]. Also, along this direction, generalizing our work to other shapes for the subsystem should also be interesting.
We then proceed to study HEE for the sound mode excitation of the fluid. Firstly, in d = 2 the sound mode is non-dissipative. The variation of HEE with respect to boundary time is, in a way trivial since it has a predictable sinusoidal behaviour mimicking that of the fluid variables in the linearized approximation. However, it is observed that the amplitude of HEE has a very interesting dependence on the momentum k of the sound waves. This behaviour is quasi-periodic, with HEE vanishing at specific values of k, which are determined by the background temperature (see fig.9a). The sound modes may be superimposed to form a traveling Gaussian wave-packet, which is a positive pressure or energy density pulse. While this pulse moves across the subsystem of interest, HEE generically experiences a jump. However, we find that when the pulse enters or exits the subsystem, there is a curious dip in HEE (see fig.10). The pulse may also be thought of as a 'semilocal' quench with the energy-momentum tensor, and similar effects on HEE may also be visible for a non-thermal background as well. It would be very interesting to investigate the physical origin of this dip in future work.
Finally, we have considered the sound mode in d = 3, 4, where it is dissipative 14 . The fluid fluctuations are damped due to the shear viscosity, and at late times these configurations settle down to equilibrium. Due to our linearized approximation, the time dependence of HEE follows this damped oscillatory behaviour of the fluid variables (see fig.12). But, the amplitude of these oscillations is significantly interesting due to a phase which is physically relevant. The most noteworthy feature in this calculation is the presence of an additional UV divergence in d = 4, which is subleading compared to the 'area law divergence'. This divergence is dynamical in the sense that its coefficient is exponentially damped at late times and therefore disappears in the final equilibrium state. In d = 3, potentially, there was the possibility of such a subleading logarithmic divergence, but it is absent due to a miraculous cancellation at first order in derivative expansion. The possibility of the existence of this divergence in d = 3 when higher order derivatives are taken into account remains open. The physical origin of this divergence in d = 4 is not very clear. Dissipation is clearly not the reason, as it is absent for the d = 3 sound mode 15 . It is possible that this is some artifact of the derivative expansion. In future investigations, it would definitely be very interesting to develop a better understanding of this mysterious 14 There is a significant recent interest in the entanglement structure of non-unitary systems (see [30] for a recent work). This part of our work is definitely in line with such developments. 15 The divergence also appears to be absent for diffusive effects in d = 4 which can be seen from the asymptotic fall off of the metric functions in [31,32]. divergence.
In our analysis of the sound mode, we have considered the fluctuations as a linearized excitation over a static background. It would be very interesting to turn on a background fluid velocity for this question. The interplay between the directions of the background fluid velocity and that of the sound mode would definitely have important consequences for holographic entanglement.
The fluid states are generally very complex macroscopic states whose underlying entanglement structure would be extremely difficult to estimate using quantum field theoretic techniques. For conformal fluids, the holographic set up provides us with a unique opportunity to compute such quantities relatively easily. Our HEE results here, are primarily of theoretical interest, and we believe, it adds to our understanding of the interplay between quantum information measures and macroscopic dynamics. With this in view, extending our line of investigation to holographic superfluids would be particularly interesting. Since our calculations rely heavily on holographic methods, this is perhaps its main limitation. Technically our results are only valid for CFTs admitting a holographic dual. However, we would expect some of the qualitative features of our result to be valid beyond the purview of holography. Also, the fluid-gravity solutions we have used are approximate solutions written in a derivative expansion. Hence our results of HEE for the dissipative sound mode in d = 3, 4 is limited by this approximation. In all our investigations we have used a strip-like subsystems, and generalization of this assumption may also lead to more interesting results in this context. Extending our entanglement entropy investigation to more complicated fluid flows beyond the sound mode (such as turbulent flows) would also be very interesting future work.