Dissipation triggers dynamical two-stream instability

Two coupled, interpenetrating fluids suffer instabilities beyond certain critical counterflows. For ideal fluids, an energetic instability occurs at the point where a sound mode inverts its direction due to the counterflow, while dynamical instabilities only occur at larger relative velocities. Here we discuss two relativistic fluids, one of which is dissipative. Using linearized hydrodynamics, we show that in this case the energetic instability turns dynamical, i.e., there is an exponentially growing mode, and this exponential growth only occurs in the presence of dissipation. This result is general and does not rely on an underlying microscopic theory. It can be applied to various two-fluid systems for instance in the interior of neutron stars. We also point out that under certain circumstances the two-fluid system exhibits a mode analogous to the r-mode in neutron stars that can become unstable for arbitrarily small values of the counterflow.


I. INTRODUCTION AND CONCLUSIONS
Two-fluid (or multi-fluid) systems are realized in various settings. They can be created in the laboratory, for instance through Bose-Fermi mixtures of cold atomic gases [1][2][3]. They may also appear in neutron star cores, where neutrons, protons, and electrons, possibly together with hyperons, can be described as a multi-fluid or even a multisuperfluid system [4][5][6][7]. Besides the core of a neutron star, the inner crust hosts a complicated system of superfluid neutrons (and electrons) immersed in a lattice of ions [8,9]. This is somewhat comparable to an atomic superfluid in an optical lattice, which shows instabilities analogous to the two-fluid systems discussed here [10]. Also, any single superfluid at finite temperature is a two-fluid system, in which the condensate and the collective excitations act as independent fluids [11][12][13]. Two-fluid models have also been proposed in the context of the quark-gluon plasma to account for strongly and weakly coupled sectors [14] and for dark matter being gravitationally coupled to normal matter in neutron stars [15,16] and quark stars [17].
Here we are interested in two miscible relativistic fluids which are coupled to each other but have independent velocity fields. We shall work in a hydrodynamic formalism, and our main result does not depend on the microscopic theory from which the equation of state and transport coefficients can be obtained. (We will use a specific equation of state to illustrate our results.) The validity of our two-fluid hydrodynamics relies on the following assumptions: the mean free paths from collisions of the microscopic constituents of each fluid among themselves are much smaller than the system size (that's why each of the components can be treated as a fluid), and the mean free path from inter-collisions, i.e., between constituents of one fluid with the other, is larger than those from intra-collisions (that's why the system does not simply become a single fluid), but still smaller than the system size (otherwise the fluids would not interact).
The main motivation for this study is as follows: If the two coexisting fluids move relative to each other, an instability can be expected at a certain critical velocity, termed two-stream instability or counterflow instability [18][19][20]. Here we are only discussing homogeneous systems, i.e., the fluids interpenetrate each other and we do not consider interfaces between the fluids, which also exhibit counterflow instabilities [21]. Instabilities of relativistic ideal fluids in this homogeneous scenario were discussed systematically in Refs. [22,23]. It was found that dynamical instabilities only occur at very large velocities, where there is already an energetic instability, which sets in at a smaller counterflow. Here, dynamical instability refers to an exponentially increasing mode with a certain growth time. Since we shall only work with linearized hydrodynamic equations, this growth time is valid for the initial stage of the instability before nonlinear effects become important. An energetic instability, on the other hand, refers to a negative energy, but there is not necessarily a growth time associated with it. The negative energy indicates that the system wants to relax to a state of lower energy. For this instability to turn into an actual unstable process, i.e., for it to turn dynamical, the coupling to an environment is needed, for instance by dissipative processes. Demonstrating and explaining this effect for the two-stream instability in a two-fluid system is the main goal of this paper. In our hydrodynamic approach, the energetic instability is associated with a change of direction of a sound mode from upstream to downstream due to an increasing counterflow between the fluids. In a more microscopic approach, this flip over is accompanied by the energy of a quasiparticle excitation turning negative [23]. This is completely analogous to Landau's original argument for the instability of a superflow, which thus is an example for an energetic instability. The relation between energetic and dynamical instabilites is also discussed for instance in Refs. [10,24,25]. Another example of an energetic instability is a negative canonical energy of rotating fluids [26], and we shall point out further intriguing parallels of our observations to rotational instabilities in neutron stars below.
We consider a two-fluid system in which one of the fluids is dissipative. Dissipation is taken into account by adding terms of first order in space-time derivatives to the stress-energy tensor. In a single fluid, this gives rise to three dissipation channels, characterized by three transport coefficients: heat conductivity, shear viscosity, and bulk viscosity. Two-fluid systems allow for more coefficients. This is known for instance from a superfluid at finite temperature, which has one dissipative and one non-dissipative component, resulting in three bulk viscosity coefficients [27,28]. For our purpose it is sufficient to consider only the dissipative coefficients that are already present in the absence of the second fluid and thus ignore additional dissipative effects due to the counterflow between the fluids. We also do not discuss two-fluid systems where both fluids would be dissipative in the absence of the other fluid. Moreover, for our main result, we shall restrict ourselves to zero temperature, i.e., only viscosity, not heat conduction, will play a role. In this sense, we are constructing the "minimal" scenario which shows our main result: dissipation triggers a dynamical instability exactly at the critical velocity at which, in the absence of dissipation, an energetic instability would occur.
It is known that in first-order hydrodynamics unphysical instabilities occur even for a single fluid [29], which are cured by going to second order [30] (or possibly, as recently suggested [31,32], by generalizing the common first-order formulations by Eckart and Landau/Lifshitz). We shall therefore start with a careful discussion of a single relativistic dissipative fluid (in an arbitrary rest frame) and identify the unphysical modes before moving on to two coupled fluids. Then, by working in the rest frame of the dissipative fluid and at zero temperature, we ensure that the modes with unphysical instabilities do not mix with the modes we are interested in. As a further confirmation that the two-stream instabilities considered here are not artifacts of first-order relativistic hydrodynamics, we shall also discuss the non-relativistic limit, where no unphysical instabilities occur and where the same conclusion is reached. The observation of the coincidence of an energetic instability and a dynamical instability due to dissipative effects has also been made in the study of a holographic superfluid [33].
Our system may be considered as a toy version of a finite-temperature superfluid. In that case, dissipative and nondissipative components are connected by a microscopic equation of state, and the dissipative component is identical to the entropy current. In the present study we keep the coupling between the fluids completely general, allowing also for entrainment -a nontrivial mixing of the currents and their conjugate momenta. In the future it would be interesting to repeat the calculation using an actual equation of state for a finite temperature superfluid, for instance along the lines of Refs. [34,35]. The generality of our setup allows us to decide whether the non-dissipative fluid behaves as a normal fluid or a superfluid. In a superfluid the hydrodynamic equations are supplemented by a Josephson equation which shuts off certain modes that exist in a normal fluid. Our main result does not depend on whether or not this constraint is taken into account, i.e., it holds for both normal fluids and superfluids.
The results of our study are mostly of general value, but they also have concrete implications. For instance, it has been suggested that two-stream instabilities may provide the trigger for a collective unpinning of vortices from the ion lattice in the inner crust of a neutron star, which is necessary to explain sudden jumps in the rotation frequency, socalled pulsar glitches, see for instance Refs. [36][37][38][39]. By comparing the estimated lag between the rotating superfluid and the crust from glitch data with the regions in parameter space where the two-stream instability sets in it has been argued that the instability only affects very large multipole moments [22]. This estimate, however, was based on the ideal fluid approximation. Our results show that dynamical instabilities in realistic systems set in at lower counterflow velocities, which suggests that the instability may operate at lower multipole moments and thus on larger length scales. Therefore, further studies which connect our general results to concrete systems such as neutron star crusts, but also to two-fluid systems in the laboratory, would be of great interest.
Putting the result in a wider context, we note that our observation is not confined to a linear relative flow between two fluids. In Fig. 1 we show the analogy between our system and the case of a rotating star where certain oscillation modes suffer the so-called Chandrasekhar-Friedman-Schutz (CFS) instability [26,42], which manifests itself in the emission of gravitational waves. In this figure we have sketched the mechanism for f -modes -"fundamental" modes of the star whose instabilities are driven by dissipation or gravitational radiation [43,44]. However, the phenomenologically more interesting instability is predicted for r-modes [45] because it sets in at much smaller rotation frequencies.
(For a recent review about f -mode and r-mode instabilities and their significance for gravitational wave emission of isolated neutron stars see Ref. [46].) The r-modes only exist in a rotating star and, in the absence of any damping mechanism, become unstable at infinitesimally small rotation frequencies. Translated to our two-fluid system, r-modes would find an analogy if there is a mode that only propagates in the presence of a counterflow and which, if dissipation relative flow v sound modes for v > v c sound modes for v = 0 f-modes for rotating star Ω f-modes for c negative energy modes Figure 1: Analogy between the two-stream instability (upper row, discussed in this paper), and the CFS instability in neutron stars, here schematically shown for f -modes (lower row). Left column: The thick (black) arrow indicates the relative velocity v of the second fluid, viewed in the rest frame of the first fluid, or the rotation of the star with frequency Ω, viewed from an inertial frame. Middle column: In the absence of a relative motion of the fluids or a rotation of the star, a given mode can propagate with the same speed in two opposite directions, indicated by two thin (red and blue) arrows of the same length. Right column: At a sufficiently large counterflow between the fluids or a sufficiently large rotation frequency of the star, one of the modes flips over. This point marks the onset of an energetic instability. (In the astrophysics literature, instabilities indicated by negative energies are often referred to as "secular instabilities" [40,41].) If the system is now allowed to interact with its environment by dissipative effects (both cases) or losing angular momentum through the emission of gravitational waves (star), the energetic instability turns dynamical.
is switched on, becomes immediately unstable. We demonstrate the existence of such a mode if entrainment between the two fluids is taken into account. Intriguingly, this observation provides us with a system that is much simpler than the rotating star coupled to gravitational waves and yet shows qualitatively the same instability. This suggests that the r-mode instability is a more general phenomenon than previously thought and not only a consequence of the complicated details of rotating stars and general relativity. It would be interesting to exploit the analogy to the two-fluid system further in the future.

A. Setup
To prepare the two-fluid calculation and establish our notation, we start by discussing a single dissipative fluid in first-order hydrodynamics. We are interested in oscillatory modes, which will be obtained in the usual way by evaluating the hydrodynamic equations to linear order in the fluctuations. With the two-fluid application in mind, we shall discuss these modes in a general rest frame. This generalizes the well-known results of Ref. [29], which are obtained in the rest frame of the fluid. The present formalism for relativistic dissipative fluids has of course been frequently applied in the literature, see for instance Refs. [28,[47][48][49].
Our starting point is the conserved current j µ and the stress-energy tensor T µν , where n is the number density in the rest frame of the fluid, and v µ is the four-velocity, which obeys v µ v µ = 1. We work in Minkowski space with metric convention g µν = (1, −1, −1, −1). We have adopted a formulation in which the current does not receive dissipative corrections. This is referred to as the Eckart frame [50] (in contrast to the Landau/Lifshitz frame formulation [51]). Up to first order in derivatives, the ideal and dissipative contributions are where and P are energy density and pressure in the rest frame of the fluid, and ∆ µν ≡ g µν − v µ v ν . Here, q µ is four-transverse to v µ (i.e., q · v ≡ q µ v µ = 0), and τ µν is four-transverse to v µ , traceless, and symmetric. The explicit form of the dissipative contribution in terms of the transport coefficients is given by with temperature T , shear viscosity η, bulk viscosity ζ, and heat conductivity κ. For the explicit calculation it is useful to write the heat and shear terms in the stress-energy tensor as where we have used v γ ∂ µ v γ = 0 and abbreviated . The hydrodynamic equations are given by the conservation of current, energy, and momentum, It is instructive to compute the entropy production, which also serves as a consistency check for the signs of the dissipative contributions. Treating the chemical potential p and the temperature T as our independent thermodynamic variables 1 and denoting the entropy density by s, we compute where all derivatives with respect to p are taken at fixed T and vice versa. We have used the thermodynamic relations dP = ndp + sdT , d = pdn + T ds, n = ∂P ∂p , s = ∂P ∂T , and + P = pn + sT , which imply dn = ∂n ∂p dp + ∂s ∂p dT and ds = ∂n ∂T dp + ∂s ∂T dT . With ∂ µ j µ = 0 and the definition of the dissipative part of the stress energy tensor (2b), the contraction of Eq. (6b) with the four-velocity yields where the transversality of q µ has been employed. The terms on the right-hand side are rewritten with the help of the explicit forms of the dissipative terms (3), With these expressions and ∂ µ T µν = 0 Eq. (7) becomes where we have defined the entropy current Each of the terms in the entropy production (9) is non-negative provided that η, ζ, κ > 0, as it should be. For instance, in the rest frame of the fluid we can easily write the entropy production terms in the form of squares, where we have kept derivatives of the three-velocity v, but set v = 0 otherwise. The three-velocity is defined via Here and in the following, the components of the three-velocity are denoted by v i , which should not be confused with the spatial components of the four-velocity (which in general have an additional factor γ). One can also write the entropy production (9) in terms of the so-called thermodynamic forces ∂ µ v ν and [∂ µ T − T (v · ∂)v µ ]/T and the corresponding thermodynamic fluxes π µν = τ ∆ µν − τ µν and q µ , which can be compared to the non-relativistic form, see for instance Ref. [9] and Eq. (55b) in our discussion of the non-relativistic limit.
In order to compute the sound modes of the system we introduce harmonic fluctuations in chemical potential, temperature, and three-velocity about their equilibrium values, We insert this ansatz into Eqs. (6) to obtain the non-dissipative contributions [we omit the exponential e i(ωt−k·x) , which multiplies every term on the right-hand sides] Here we have used ∂ µ j µ = 0, have defined the enthalpy density and have abbreviated To compute the dissipative contributions we linearize in the fluctuations δp, δT , δv. In this approximation, and using Eqs. (2b) and (4), we can write where ≡ ∂ µ ∂ µ . This yields the spatial and temporal components (again omitting the exponential) B. Sound modes at T = 0 in the rest frame of the fluid As a warm-up, let us start with the simplest case: zero temperature, T = s = κ = 0, and vanishing fluid velocity, v = 0. In this case, the temporal component of the energy-momentum conservation, ∂ µ T µ0 = 0, is automatically fulfilled. We contract the spatial components ∂ µ T µi by k i and δ ij −k ikj to obtain (three-)longitudinal and (three-)transverse components with respect to the direction of propagation of the mode. Together with the current conservation, this yields where δv ⊥,i = (δ ij −k ikj )δv j , and is the squared sound speed. From Eqs. (19a) and (19b) we find the longitudinal modes (where δv ⊥ = 0) where we have denoted the attenuation constant by (Here, at zero temperature, w = pn.) The attenuation constant has units of time. Note, however, that the damping time scale from the imaginary part in the k 2 contribution is set by (k 2 Γ 0 ) −1 and thus depends on the wavenumber k. The modes (21) are damped since Im ω = k 2 Γ 0 > 0. Imaginary parts of ω with the opposite sign would indicate unstable modes. From Eq. (19c) we find one transverse mode, where δv ⊥ = 0 but δp = k · δv = 0, with This is a purely diffusive mode which does not propagate.

C. Sound modes at nonzero T in general frame
Let us now discuss the general single-fluid case. In particular we allow for a nonzero velocity v. For a single fluid without boundaries such as walls of a capillary and without coupling to a second fluid, this nonzero velocity should of course not lead to any new physics since it is nothing but a Lorentz boost of the whole system. Nevertheless, the calculation will turn out to be useful as a preparation for the two-fluid system.
The equations ∂ µ j µ = ∂ µ T µν = 0 are 5 scalar equations for the 5 scalar variables δp, δT , δv. Since we allow for a general frame, we now have a second three-vector in the problem, namely the fluid velocity v, in addition to the wavevector k. We align the z-axis with k and the y-axis such that the fluid velocity v is in the y-z plane. Instead of the components δv x , δv y , δv z we work with δv x , k · δv, and v · δv, and instead of ∂ µ T µi = 0 we work with where θ is the angle between k and v. Since the velocity oscillations of these modes are transverse to both k and v we term them T ± . Their expansion up to second order in k is shown in Table I. We see that by setting v = T = 0 in the mode T + we recover Eq. (23). In other words, the purely diffusive mode T + now becomes propagating in the boosted fluid. This seems to be a trivial observation, but we shall see that these modes become particularly interesting in the case of two coupled fluids with entrainment, where they can couple to other modes, see Sec. III C. The mode T − has an unphysical instability, which is an artifact of first-order hydrodynamics [29]. Going to second-order hydrodynamics cures this unphysical instability [30] 2 . The mode T − did not appear in the previous subsection, where we worked in the limit v = T = 0 from the beginning. It is useful to keep this in mind for our discussion of the two-fluid system, where we shall work in the zero-temperature limit and in the rest frame of the dissipative fluid, such that the unphysical instability plays no role. The remaining 4 coupled equations give rise to a matrix of coefficients whose determinant is a polynomial of degree 6 in ω, which factorizes into a quadratic and a quartic polynomial. The quadratic polynomial gives exactly the same dispersion as for the modes T ± (24). The four solutions of the quartic polynomial are obviously complicated in general. Two of them only exist at nonzero temperatures, i.e., they are not connected to any of the modes found in the previous subsection. They are (for brevity we only give the dispersion of this mode for v = 0) and The mode (25) shows an (unphysical) instability, just like the mode T − . The two modes (25) and (26) are included in Table I in the last two rows. Table I: Coefficients of the single-fluid modes in a low-momentum expansion of ω(k) up to k 2 , with v-dependence shown only where the expression is short enough, i.e., for L ± and the mode in the 6th row only the v = 0 result is given. The modes in the 4th and 6th row show an unphysical instability in first-order hydrodynamics. The modes of the bottom two rows only exist at nonzero temperature, and play no role in our discussion of the two-fluid system, where we restrict ourselves to zero temperature.
Finally, there are two stable, damped modes that are generalizations to nonzero temperature of the two modes (21).
Here, all derivatives with respect to the chemical potential p are taken at fixed temperature T and vice versa, while the derivative with respect to in Eq. (27a) is taken at fixed entropy per particle; see for instance appendix E of Ref. [52] for a derivation of this change of thermodynamic variables. At zero temperature, the sound speed (27a) reduces to the expression (20). We term these two modes L ± , see first two rows of Table I. In Fig. 2 we plot the 3 stable modes that exist at T = 0, i.e., L + , L − , and T + , as a function of v. It is useful to understand this figure before moving to the two-fluid system. The modes L + and L − correspond to sound waves propagating in opposite directions, parallel and anti-parallel to the fluid velocity v. (Equivalently, we can reinstate the dependence on θ and plot a single mode, say L + , for θ = 0 and θ = π.) For an external observer, the mode L + is sped up by the fluid flow, while the mode L − is slowed down until it flips over at v = c. This inversion of direction is not accompanied by any dynamical instability in the case of a single fluid: the attenuation is positive for all fluid velocities. This is expected since the fluid is stable at rest and isolated and thus no dynamical instabilities should be created by simply moving the fluid. The figure also illustrates, as already mentioned above, that the mode T + only appears to propagate if the fluid moves, and we see that its damping time is determined solely by the shear viscosity η. This is in contrast to the modes L ± , whose damping times are given by a certain linear combination of shear and bulk viscosity.
These results can also be obtained by applying a Lorentz boost Λ µν with 3-velocity v to k µ = (ω, k) with ω = ck + iΓk 2 . Then, denoting the transformed 4-momentum by k µ = (ω , k ) with ω = c k + iΓ k 2 , the relation Λ µν k ν = k µ yields the transformed sound speed c and the transformed attenuation Γ . For boosts in the direction of the wave propagation we find Obviously, the transformed speed c is simply given by the relativistic addition of velocities. One can check that the velocity dependence of the attenuation shown in the right panel of Fig. 2 is indeed given by Γ .

A. Setup
We set up the two-fluid system by first introducing the currents and stress-energy tensor for two ideal, coupled fluids at zero temperature. Here we follow the formalism used for instance in Refs. [5,13,23,53,54]. As a second step, we shall then add the dissipative terms from the previous section for one of the fluids.
As our independent hydrodynamic variables we choose the conjugate momenta p µ n , where the index n = 1, 2 labels the two fluids. We denote the components of the conjugate momenta by p µ n = (µ n , p n ), where the temporal components µ n are the chemical potentials measured in the laboratory frame, and the spatial components are related to the three-velocities v n by p n = µ n v n (neither here nor anywhere in the following do we use n as a summation index). The chemical potentials in the rest frames of the fluids are p n = (p n,µ p µ n ) 1/2 = µ n (1 − v 2 n ) 1/2 . Conjugate momenta and four-velocities are thus related by 3 In a superfluid, temporal and spatial components of the conjugate momentum are not independent since the conjugate momentum can be written as the gradient of a scalar field, p µ = ∂ µ ψ. This ensures the irrotationality of p (in a nonrelativistic superfluid v is irrotational, but this is not true in general). From a microscopic point of view, ψ is the phase of a condensate, for instance a condensed scalar field in a bosonic theory [13] or the Cooper-pair condensate in a fermionic superfluid. The central quantity that encodes the microscopic physics is the so-called generalized pressure Ψ. At zero temperature, it depends on the Lorentz scalars that can be constructed from the conjugate momenta, Ψ = Ψ(p 2 1 , p 2 2 , p 2 12 ), where we have abbreviated p 2 12 ≡ p 1 · p 2 . The conserved currents are obtained by Consequently, we can write where In general, the currents are not proportional to their own conjugate momentum, but, if A = 0, receive a contribution from the conjugate momentum of the other fluid. This effect is called entrainment or Andreev-Bashkin effect [55], and A is the entrainment coefficient. In the present formalism entrainment arises if the generalized pressure depends on p 2 12 . Of course, even in the absence of entrainment the two fluids can be coupled to each other, for instance through terms in the pressure proportional to p 1 p 2 . We shall set up the hydrodynamic equations in full generality, including entrainment terms. Then, for our main result we ignore entrainment, which simplifies the calculation significantly, see Sec. III B. In Sec. III C, where we discuss the analogue of the r-mode instability, entrainment is included and turns out to be crucial.
In the absence of dissipation, the two-fluid stress-energy tensor is The continuity equations for the two conserved currents are and we can write the divergence of the (ideal) stress-energy tensor as with the vorticities To derive Eq. (35) we have used the continuity equations (34) and dΨ = j µ 1 dp 1,µ + j µ 2 dp 2,µ . We see that from the total energy-momentum conservation alone only the sum of the vorticity terms j 1,µ ω µν 1 and j 2,µ ω µν 2 (plus dissipative terms which are not included for now) is required to vanish. It is usually argued, for instance starting from a variational principle [54], that instead there are two separate vorticity equations, j 1,µ ω µν 1 = j 2,µ ω µν 2 = 0. This corresponds to using two separate Euler equations in a non-relativistic two-fluid formalism [4,7]. If one of the fluids, say fluid 2, is a superfluid, i.e., p µ 2 = ∂ µ ψ, then its vorticity vanishes, ω µν 2 = 0. In this case, the energy-momentum conservation becomes equivalent to the vorticity equation j 1,µ ω µν 1 = 0 if fluid 1 is non-dissipative or j 1,µ ω µν 1 + ∂ µ T µν diss = 0 in the presence of dissipation. In Sec. III B this is the situation we shall consider, i.e., we shall assume fluid 1 to be dissipative and fluid 2 to be a superfluid. The main results of that section can also be obtained by simply considering two separate vorticity equations (one of them supplemented by the dissipative terms). The main result of Sec. III C, however, relies on using only the total energy-momentum conservation and not assuming the non-dissipative fluid to be a superfluid. We shall comment on this point in more detail in Sec. III C.
For the explicit evaluation of the current conservation equations (34) we need the derivatives of the coefficients B 1 , B 2 , A, i.e., the second derivatives of the generalized pressure, which we abbreviate by The coefficients b 1 and b 2 are nonzero even in the absence of a coupling between the fluids and are related to the sound speeds of the fluids. The other coefficients characterize the coupling between the fluids with a 1 , a 2 , a 12 being nonzero only in the presence of entrainment, while d describes the non-entrainment coupling.
As in the case of a single fluid, we introduce fluctuations and linearize the conservation equations. As our independent fluctuations we choose the chemical potentials and velocities in the laboratory frame, δv 1 , δv 2 , δµ 1 , δµ 2 . Then, the divergence of the currents becomes [we omit the exponential e i(ωt−k·x) ] and the same with 1 ↔ 2 for ∂ µ j µ 2 . For the energy-momentum conservation we need the temporal and spatial components where (1 ↔ 2) stands for repeating all preceding terms with subscripts 1 and 2 exchanged. In the non-dissipative part one can recover the two contributions from the vorticities (35): all terms containing fluctuations of fluid 1, i.e., δµ 1 and δv 1 , combine to j 1,µ ω µν 1 , while the terms containing fluctuations of fluid 2 combine to the vorticity terms of fluid 2.
Dissipation is added as explained in the introduction: we use the single-fluid dissipative terms for one of the fluids and leave the second fluid unchanged. In particular, we only consider the dissipative coefficients that already appear in single-fluid first-order hydrodynamics, and thus ignore any additional dissipative coefficients that occur due to the relative motion of the two fluids. Since in this section we have restricted ourselves to zero temperature, the only relevant transport coefficients are η and ζ. Choosing without loss of generality fluid 1 to be dissipative, the zero-temperature limit of of Eqs. (18) yields the dissipative contribution,

B. Dynamical instability triggered by dissipation
We perform the calculation of the sound modes in the rest frame of the dissipative fluid, i.e., v 1 = 0 (of course, the fluctuations δv 1 are allowed to be nonzero). In this frame, there are two particular directions: the direction of propagation of the sound wave k and the direction of the (uniform) velocity of the ideal fluid v 2 . As in the case of a single fluid we choose the y-z-plane such that both k and v 2 lie in this plane, and rather than working with the x, y, and z components we work with the x component and the two contractions with k and v 2 . Therefore, we have 8 scalar fluctuations, As mentioned above, let us assume the non-dissipative fluid to be a superfluid, such that p µ 2 = (∂ t ψ, −∇ψ). Hence, ∂ t p 2 = −∇µ 2 , which, within our plane wave ansatz, results in the constraint for the fluctuations ωδp 2 = kδµ 2 or, equivalently, ωµ 2 δv 2 = (k − ωv 2 )δµ 2 . As a consequence, superfluids only allow for longitudinal modes, where p oscillates in the direction of k. We can thus express the fluctuations δv 2 in terms of δµ 2 , which reduces the number of independent variables to 5, namely δµ 1 , δµ 2 , k · δv 1 , v 2 · δv 1 , δv 1x . We find that in this setup the temporal component of the energy-momentum component is automatically fulfilled, ∂ µ T µ0 = 0. This leaves us with the 5 scalar equations , matching the number of variables. We neglect entrainment for now, a 1 = a 2 = a 12 = A = 0. As for a single fluid, we find that the equation ∂ µ T µx = 0 decouples since it is the only equation that contains δv 1x . It yields the mode T + for the dissipative fluid, in which the only nonzero oscillation is δv 1x . Since this mode is completely unaffected by the second fluid, it is not relevant for the following discussion. Notice that the unphysical mode T − does not appear because we work in the rest frame of the dissipative fluid. We are left with 4 equations for 4 variables, which read Here we have eliminated B n and b n in favor of more physical quantities: First, we can relate B n to the densities n n = j n · v n in the respective fluid rest frames. With Eqs. (29) and (31) In the absence of a coupling between the fluids, c 1 and c 2 are the sound speeds of the two fluids, measured in their respective rest frames. The coupling leads to a different (much more complicated) form of the sound speeds, but also affects the thermodynamic quantities B n and b n themselves, i.e., even the quantities c 1 and c 2 as defined in Eq. (44) change as the coupling is varied. Third, we have introduced the dimensionless "mixed susceptibilities" and the same for 1 ↔ 2, which are only nonzero in the presence of a coupling between the fluids. The attenuation constant Γ 0 is as defined in Eq. (22), and additionally we have abbreviated 4 With entrainment we would have and the same with 1 ↔ 2 for n 2 .
The system of equations (42) is further simplified by the observation that v 2 · δv 1 only appears in the last equation (i.e., the 4th column of the 4 × 4 matrix has a nonzero entry only in its last component). This equation can thus be ignored for the calculation of the dispersions. It can later be used to compute v 2 · δv 1 , if needed. We are left with the determinant of a 3 × 3 matrix to determine the dispersion relations ω(k). This determinant is a quartic polynomial in ω, yielding 4 modes. In the absence of a coupling between the fluids, these modes simply correspond to the two modes L ± for each fluid. When the fluids are coupled, the dispersions are very complicated in general. Nevertheless, our main result can be extracted from the general result and stated in a very compact form. To this end, as in Sec. II, we expand the modes in the momentum up to second order, ω ck + iΓk 2 . We expect that one of the modes flips over at a certain relative velocity between the fluids, say v 2 = v 20 such that c(v 20 ) = 0. This is qualitatively the same effect as for a boosted single fluid, see Fig. 2, and constitutes an energetic instability. However, due to the presence of the second fluid we now find that at exactly this point the attenuation constant of this mode changes its sign too, and as a consequence the mode changes from being damped to being dynamically unstable. This observation can be formulated analytically in the vicinity of the critical velocity v 20 . To linear order, the sound speed close to v 20 is with v 20 = c 2 where θ is the angle between the direction of propagation k and the direction of flow of the non-dissipative fluid v 2 . The corresponding attenuation is to linear order In particular, close to the critical velocity Γ does not depend on the angle θ. The dynamical instability sets in first (i.e., at the smallest relative velocity) for modes (anti-)parallel to v 2 . In this case, cos 2 θ = 1 and thus v 20 = c 2 , and we see that c changes its sign from negative (upstream) to positive (downstream) as we move through the critical point by increasing the relative velocity, while Γ changes its sign from positive (damped) to negative (unstable). The time scale for the exponential growth of this mode is given by Eq. (49). It depends on the bulk and shear viscosity coefficients, but also on the thermodynamic properties of the two fluids and their coupling strength. We recall that for this result entrainment plays no role. Of course, nonlinear effects will set in as the amplitude of the mode grows, in the current linear approach we can make no prediction about the fate of the system in the dynamically unstable regime.
To further illustrate this result we plot all four modes and their attenuation as a function of the relative velocity. For this numerical evaluation we need to choose numerical parameters for the quantities appearing in Eq. (42). If we were to choose random numbers we would have no control over whether the two fluids can coexist, and thermodynamic instabilities -which we are not interested in -might manifest themselves in unstable hydrodynamics modes and thus obscure the interpretation of our results. Therefore, we employ the following form of the generalized pressure, which is based on a relativistic field-theoretical model for two complex scalar fields with masses m 1 , m 2 , self-coupling constants λ 1 > 0, λ 2 > 0, a non-entrainment coupling h and an entrainment coupling g [23]. The pressure (50) is valid in the parameter space where condensates of both fields coexist, see Ref. [23] for a systematic discussion of the parameter space. In all following results, we choose our parameters such that at vanishing relative velocity this phase of miscible fluids is stable. Strictly speaking, the pressure (50) refers to a system of two superfluids because the vacuum expectation values (= Bose-Einstein condensates) of the two scalar fields constitute the fluids. However, for our purpose it is sufficient to consider Ψ as an example for a microscopic equation of state even though we assume one of the fluids to be dissipative. We do not need to provide any microscopic form or specific numbers for the viscosity coefficients because the attenuation Γ scales with Γ 0 (22) and otherwise does not depend on η and ζ. We show the result for the 4 upstream/downstream modes in Fig. 3 with parameters m 1 = m 2 ≡ m, p 1 = 12m, p 2 = 10m, λ 1 = 0.2, λ 2 = 0.3. Since we have ignored entrainment in the hydrodynamic calculation, we obviously have to set the microscopic entrainment coupling to zero, g = 0, for consistency. The results are shown for weak coupling, h = −0.03 (left panels) and strong coupling h = −0.12 (right panels). Other than creating a stable system of coexisting fluids at vanishing counterflow, this particular choice of parameters has no specific physical motivation. Somewhat different results are obtained for h > 0, but the main qualitative conclusions are the same and thus we restrict ourselves to h < 0. The main observations are: • The sound speeds themselves show a dynamical instability at large velocities: where the real parts of two branches coincide, their imaginary parts (not shown in the figure) have opposite signs, i.e., one of the modes grows exponentially. Since the results linear in k are not affected by dissipative effects, this is the same observation as already made in Refs. [22,23]. As a consequence, the system is dynamically unstable for all v 2 > v 20 : First due to Re Γ < 0 with a growth time for the instability τ = (k 2 |Re Γ|) −1 , i.e., at the onset the growth time is infinite. Then due to Im c < 0 with growth time τ = (k|Im c|) −1 , while there is no instability in the k 2 terms in this regime, i.e., Re Γ > 0. Interestingly, just before and just after this regime Re Γ diverges, indicating an extreme instability with infinitesimally short growth time. Then, at the highest velocities, the instability is given again by the second-order contribution, Re Γ < 0.
• There is a "role reversal" between damped and undamped modes as a function of the counterflow velocity: For small v 2 there is a pair of upstream modes (c > 0) and a pair of downstream modes (c < 0). Each pair has a strongly damped and weakly damped mode. As v 2 increases, the strongly damped mode becomes weakly damped and vice versa. For small coupling, this role reversal happens in the stable regime, i.e., before the dynamical instability sets in.

C. Analogue of the r-mode instability
The dynamical instability just discussed sets in at a certain nonzero critical velocity. As explained in Sec. I, this is analogous to modes in rotating stars, such as f -modes, which become unstable at a nonzero angular velocity of the star. In contrast, r-modes -whose restoring force is the Coriolis force -only exist in a rotating star and can become unstable at arbitrarily small angular velocities. Exploiting this analogy, it seems natural to ask whether there is an unstable mode in the present two-fluid system with the same properties.
As before, we work in the rest frame of the dissipative fluid. The conservation equations for the two currents and the stress-energy tensor yield the 6 equations ∂ µ j µ 1 = ∂ µ j µ 2 = ∂ µ T µ0 = v 2i ∂ µ T µi = k i ∂ µ T µi = ∂ µ T µx = 0. Now, in contrast to the previous subsection, we do not assume the non-dissipative fluid to be a superfluid, i.e., we allow for both longitudinal and transverse oscillations and thus there are a priori no constraints on the 8 variables (41). Neither do we consider two separate vorticity equations (one of which would be supplemented by the dissipative terms). This is in accordance with all conservation laws of the system, but seems to contradict the usual approach of separating the vorticity equations (or Euler equations in the non-relativistic limit), as already briefly discussed below Eq. (36). The idea behind separate vorticity equations is that the exchange of energy and momentum between the two fluids occurs on larger length and time scales than the exchange within the fluids. By using a single vorticity equation we assume that this difference in scales is not too large, but still sufficiently large to justify the two-fluid picture with independent velocity fields. It is an interesting general question whether such a regime exists, and we leave a more detailed study of this question for the future. Here, our main purpose is to point out a solution of the hydrodynamic equations that shows a dynamical instability for arbitrarily small counterflow. At the very least, this is an interesting toy version of the r-mode instability in neutron stars. Additionally, it would be very interesting to see whether it can be realized in a real-world two-fluid system with counterflow.
In order to calculate the possible modes of the system, we first note that, as above, the equation ∂ µ T µx = 0 decouples because it only contains the fluctuations δv 1x and δv 2x , and these fluctuations appear nowhere else. We are thus left with 5 equations for 6 fluctuation variables. We find propagating modes if we impose the additional constraint δµ 1 = 0 or δµ 2 = 0. In the following we shall restrict ourselves to the case δµ 1 = 0, i.e., we are looking for modes in which the chemical potential of the dissipative fluid does not oscillate. (The main qualitative result can also be obtained with δµ 2 = 0.) Such constraints on the oscillations can in principle be realized experimentally. An example is the so-called fourth sound in a superfluid, where the velocity of the normal component is clamped by an external force and thus cannot oscillate [27,56,57].
We thus consider 5 hydrodynamic equations plus the additional equation δµ 1 = 0. The determinant of the matrix of this system of equations yields a fourth-order polynomial from which four modes are calculated. For the main observation it is crucial to include entrainment. This renders the calculation very tedious, and we have not found any compact analytical results without making use of an explicit equation of state. Let us therefore first plot the results using the generalized pressure (50) and a particular parameter set, and then show the main result in an analytical form for this specific equation of state. In Fig. 4, the upstream/downstream modes are shown for the parameters p 1 = p 2 ≡ p, m 1 = m 2 = 0, λ 1 = 0.2, λ 2 = 0.3, h = 0.05 and three different entrainment couplings g = 0, g = −0.002/p 2 , g = +0.002/p 2 . The parameters are chosen such that at zero counterflow the system of coexisting fluids is stable. Again, we do not need to specify any values for the dissipative coefficients since the attenuation of the modes we consider scales with Γ 0 . In the presence of entrainment, Γ 0 as defined in Eq. (22) depends on the relative velocity v 2 ; for the normalization in Fig. 4 we use Γ 0 (v 2 = 0). We plot only three of the four modes since the fourth mode is a purely diffusive mode which does not couple to the other modes. First of all we observe the phenomenon already discussed in Sec. III B: there is one mode that becomes dynamically unstable at a nonzero critical velocity. This is independent of entrainment. The new observation is related to the mode that, in the absence of entrainment, is given by ω = v 2 cos θ. This is the mode T + in Table I. Since we work in the rest frame of the dissipative fluid, it is the ideal fluid that is in motion and thus we expect this mode to be undamped if there is no coupling to the fluid at rest. This expectation is borne out in the left plots, where the entrainment is set to zero. We see that even in the presence of a non-entrainment coupling h, which is nonzero in all plots of the figure, the mode does not couple to any other modes and thus its attenuation vanishes. Only if we switch on an entrainment coupling g this mode couples to the other modes. This mode mixing was already pointed out in Ref. [23]. Here, in the presence of dissipation, we see that this mode acquires a nonzero attenuation if the entrainment coupling is switched on. Depending on the sign of g, the mode is damped (middle plots, g < 0) or becomes dynamically unstable (right plots, positive g > 0). For small g [compared to p −2 1 , p −2 2 , (p 1 p 2 ) −1 ] and small v 2 and setting h = m 1 = m 2 = 0 for simplicity we find ω ck + iΓk 2 with This confirms analytically that this mode only propagates in the presence of a relative velocity between the fluids and becomes dynamically unstable for g > 0 and arbitrarily small relative velocities.

IV. NON-RELATIVISTIC LIMIT
The purpose of this section is to start from the non-relativistic hydrodynamic equations and re-derive the main result of this paper in this limit. As we have seen in Sec. II, first-order relativistic hydrodynamics gives rise to unphysical instabilities. We have argued that our two-fluid calculation does not suffer from this problem since we have worked in the rest frame of the dissipative fluid and at zero temperature, and thus we expect the dynamical instability we have found to survive in higher-order hydrodynamics. Nevertheless, it is useful to confirm our result in the non-relativistic case, where there are no unphysical instabilities. Here we focus on the non-relativistic version of the results of Sec. III B by considering two separate Navier-Stokes equations for the two fluids and neglecting entrainment, i.e., we do not set up a non-relativistic version of the analogue of the r-mode instability from Sec. III C.

A. Single fluid
As in the relativistic case, let us establish our formalism and notation for a single fluid before moving on to two coupled fluids. In the non-relativistic limit, the hydrodynamic equations for a single fluid are [51] ∂ρ ∂t (For a derivation of these equations from the relativistic formalism see for instance Ref. [58].) Here, ρ is the mass density of the fluid, g = ρv is the momentum density with the three-velocity v, and = 0 + ρv 2 /2 is the energy density with the energy density in the rest frame of the fluid 0 . The energy flux q and the stress tensor Π ij contain dissipative contributions and are given by with To rewrite Eqs. (52) we use the thermodynamic relations 0 + P = µρ + sT and d 0 = µdρ + T ds, where µ is the chemical potential per unit mass. (This definition of µ is common and useful in the non-relativistic literature, although here it is a slight abuse of notation since in the relativistic calculation we denoted the chemical potential itself by µ.) Then we obtain the equivalent, perhaps more familiar, set of equations in the form of continuity equation, entropy production, and Navier-Stokes equation, Next, we introduce the fluctuations in µ, T and v as in Eqs. (13) and use dP = ρdµ + sdT to find 0 where we have neglected the right-hand side of Eq. (55b), which is of second order in the fluctuations. As explained above Eq. (24), we work with the independent velocity fluctuations δv x , v · δv, k · δv. The transverse component of Eq. (56c) decouples from the rest as usual and gives a mode where only δv x is nonzero with dispersion The remaining modes are found by Eqs. (56a), (56b) together with the two equations obtained from contracting Eq. (56c) with k and v, where is the squared speed of sound in the rest frame of the fluid and, in analogy to the relativistic case, we have defined the attenuation constants The determinant of the 4 × 4 matrix yields a fourth-order polynomial in ω, such that in total there are 5 modes, all of which are stable. At T = 0, the second column and the second row of the 4 × 4 matrix vanish, and we find the 3 solutions Together with the mode (57) these are 4 (stable) modes. They are to be compared to the 4 modes L ± and T ± in Table I, one of which suffers an unphysical instability. We see that the only effect of the nonzero velocity is to add v · k to the dispersion. In particular, the attenuation is independent of v. We can obtain the modes (61) from the fully relativistic modes shown in Fig. 2 by taking the small-velocity limit, i.e., neglecting terms quadratic in the velocity, which includes neglecting terms of order vc. [See also Eq. (28), which shows that the attenuation becomes independent of the velocity if terms of order v 2 , vc and higher are neglected.] B. Two fluids For the system of two coupled fluids we restrict ourselves to zero temperature and neglect entrainment. The pressure is a function of both chemical potentials such that dP = ρ 1 dµ 1 +ρ 2 dµ 2 . We consider separate Navier-Stokes equations for the two fluids such that the relevant system of equations is where, as in the relativistic case, we have chosen fluid 1 to be dissipative, and where g n = ρ n v n . (In the presence of entrainment g 1 would receive a contribution from v 2 and vice versa.) Introducing the fluctuations as usual, working in the rest frame of the dissipative fluid, v 1 = 0, and contracting the Navier-Stokes equations with k and v 2 yields the 6 equations where, in analogy to the relativistic calculation, we have introduced "mixed susceptibilities" and the same for 1 ↔ 2. We find one undamped mode ω = v 2 · k and one purely diffusive mode ω = iΓ 1 k 2 . The remaining modes are determined by the 4 × 4 subdeterminant in the upper left corner, which yields a polynomial of fourth order in ω. We expand the dispersion as usual for small momenta, ω = ck + iΓk 2 , and find that one of these four modes flips over (sign change of c) at v 20 = c 2 cos θ , while the corresponding attenuation constant behaves in the vicinity of this point as Comparing with Eq. (49) we have rederived the dynamical instability in the limit of small c 2 . This confirms that the observation that the energetic instability turns dynamical in the presence of dissipation is not an artifact of first-order relativistic hydrodynamics.