11institutetext: Instituto de Física Teórica UAM/CSIC, Calle Nicolás Cabrera
13-15, 28049 Madrid, Spain22institutetext: Departamento de Física Teórica,
Universidad Autónoma de Madrid, Campus de Cantoblanco, 28049 Madrid,
Spain33institutetext: Wilczek Quantum Center, School of Physics and Astronomy,
Shanghai Jiao Tong University, Shanghai 200240, China & Shanghai Research
Center for Quantum Sciences, Shanghai 201315.

# A Holographic Superfluid Symphony

Daniel Areán 3 Matteo Baggioli 1,2 Sebastian Grieninger 1 and Karl Landsteiner
<EMAIL_ADDRESS><EMAIL_ADDRESS><EMAIL_ADDRESS><EMAIL_ADDRESS>
###### Abstract

We study the hydrodynamic excitations of backreacted holographic superfluids
by computing the full set of quasinormal modes (QNMs) at finite momentum and
matching them to the existing hydrodynamic theory of superfluids.
Additionally, we analyze the behavior of the low-energy excitations in real
frequency and complex momentum, going beyond the standard QNM picture.
Finally, we carry out a novel type of study of the model by computing the
support of the hydrodynamic modes across the phase diagram. We achieve this by
determining the support of the corresponding QNMs on the different operators
in the dual theory, both in complex frequency and complex momentum space. From
the support, we are able to reconstruct the hydrodynamic dispersion relations
using the hydrodynamic constitutive relations. Our analysis rules out a role-
reversal phenomenon between first and second sound in this model, contrary to
results obtained in a weakly coupled field theory framework.

††preprint: IFT-UAM/CSIC-21-80

##  1 Introduction

The holographic superfluid model Hartnoll:2008vx ; Hartnoll:2008kx ;
Herzog:2009xv ; Hartnoll:2009sz , often mislabelled as holographic
superconductor,111See Domenech:2010nf for a nice discussion about the
differences in the context of field theory and holography. is one of the most
simple and widely studied setups in the context of applied holography
ammon2015gauge ; zaanen2015holographic ; Hartnoll:2016apf ; Baggioli:2019rrs .
Its potential application in the understanding of high-Tc superconductors
Zaanen2004 has made this system the object of many thorough analyses in the
years since its inception (see Horowitz:2010gk ; Cai:2015cya for a review).
Thus, the lack of a complete study of the quasinormal mode excitations
corresponding to hydrodynamic modes of the holographic superfluid
Hartnoll_2008 is surprising. Until now, the study of the hydrodynamic modes
has been limited to the so-called probe limit, where all the gravitational
degrees of freedom (and in particular energy and momentum) are kept frozen
Amado:2009ts ; Amado:2013xya ; Amado:2013aea ; Esposito:2016ria .
In this work we fill this gap in the literature by carrying out the first
complete study of the quasinormal modes of the holographic s-wave superfluid
with backreaction Hartnoll_2008 .222The p-wave superfluid model Gubser:2008wv
; Ammon:2008fc ; Ammon:2009fe ; Cai:2013pda ; Cai:2013aca with backreaction
may be found in Ammon:2009xh and the analysis of the associated
transport/hydrodynamics in Erdmenger:2011tj ; Erdmenger:2012zu . At the same
time, a holographic model for d-wave superconductivity was introduced in
Chen:2010mk ; Benini:2010pr ; Kim:2013oba and with backreaction explored in
Ge:2012vp .  In particular, we numerically compute the spectrum of the lowest
lying quasinormal modes for the original setup Hartnoll_2008 , and also for
the modified system Andrade:2014xca ; Kim:2015dna ; Baggioli:2015zoa ;
Baggioli:2015dwa where momentum dissipation is introduced using the so-called
holographic axion model Andrade:2013gsa ; Baggioli:2021xuv . We successfully
compare our results to the predictions of superfluid hydrodynamics
Herzog:2011ec extending the results of Herzog:2009md ; Yarom:2009uq .

This study of the QNMs will allow us to carry out a novel analysis of the
support of the hydrodynamic modes of the holographic superfluid. Superfluids
display a characteristic coexistence of two different sound modes
doi:10.1063/1.3248499 which can be easily captured within the two-fluid Tisza-
Landau model tisza1938transport ; landau1941theory . In addition to the
standard density fluctuations, giving rise to the common first sound
excitation, superfluids support also a second sound mode which is driven by
temperature fluctuations. First sound can be visualized as the in-phase motion
of the normal fluid component and the superfluid one; on the contrary, second
sound is displayed as the out-of-phase collective dynamics of the two. In
Schmitt:2014eka , a field theoretical study of the nature of first and second
sound as density or entropy waves or a mixture of the two was performed. Even
though in superfluid 4He, the nature of first and second sound is preserved
for all temperatures, this seems not to be the case for ultra-cold atoms and
weakly-interacting bosons Schmitt:2014eka . Here, we carry out a similar
analysis within the holographic superfluid which stands as a putative toy
model for a strongly coupled superfluid system: we determine the support of
the hydrodynamic modes of the system on excitations of the different operators
present in the dual theory.

In order to analyze the nature of the hydrodynamic excitations of the
holographic superfluid we define a quantitative measure of the support that a
quasinormal mode excitation has on the different dual operators. We then apply
this procedure to the sound and diffusive modes present in our gravity dual.
Our results rule out any role reversal between first and second sound and
place the holographic model in the same phenomenological class of standard
superfluid 4He. This is in contrast to the weakly coupled bosonic field theory
results of Schmitt:2014eka described above. Naturally, an obvious difference
between both setups is given by the inherently strongly coupled nature of the
holographic description. Moreover, Schmitt:2014eka is neglecting dissipating
effects in their description.

In a final part of this work we study the response of the system from a
different perspective. One is typically interested in the poles of (the
Fourier-transformed) Green’s functions in the complex frequency and real
momentum space. These are the QNMs. Nevertheless one can also consider the
poles in the space of real frequency and complex momentum. These modes carry
information about the response of the system to a perturbation of the form
$e^{-i\omega t+ikx}$ where $\omega$ now is a real number and $k$ a complex
one. Although this approach is not usual in holography, it is a standard
practice in the study of propagation of waves, in particular in relation to
absorption and reflection properties of media (e.g. EM waves in a medium
jackson_classical_1999 ). Once momentum is decomposed as
$k=\mathrm{Re}(k)\,+\,i\,\mathrm{Im}(k)$, its real part defines the
propagation wavelength in real space while its imaginary part the so-called
penetration length $\lambda\equiv 1/\mathrm{Im}(k)$. In other words, this
second setup corresponds to considering the decay of a wave in space rather
than in time.

The curves $\omega(k)$ described by the QNMs in the complex plane feature a
complicated structure (e.g. the presence of different sheets, and
singularities). Hence, a priori it is not clear that the behavior of the modes
at real frequency and complex momentum can be obtained straightforwardly from
that of the standard QNMs. Indeed, notice that these two types of modes encode
responses to perturbations with different boundary conditions. Consider for
instance an elastic rod immersed in a viscous liquid. Considering complex
frequency and real momentum corresponds to keeping the extreme of the rod
fixed and creating an oscillating wave on it. The wave in that case will not
propagate in space but just oscillate (and decay) in time. On the contrary,
one could imagine exciting a wave at one of the extremes of the rod and
studying its propagation along the rod. That is equivalent to considering real
frequency and complex momentum.

In the context of holography, the modes at real frequency and complex momentum
have been overlooked and only few works considered this second setup
Amado:2007pv ; Amado:2008ji ; Forcella:2014dwa ; Blake:2014lva ;
Sonner:2017jcf . On a different note, recently, several studies of the modes
in the complex frequency/complex momentum plane have been performed in
relation to the so-called pole-skipping phenomenon Grozdanov:2017ajz and the
convergence radius of the linearized hydrodynamics expansion Grozdanov:2019kge
. In the last part of our work, we study the excitations of a holographic
superfluid at real frequency and complex momentum and discuss their physical
interpretation in detail.

#### Structure of the paper.

In section 2, we first summarize the main results of the hydrodynamic
description of a relativistic superfluid. We then introduce the holographic
model and set the framework of the analysis we carry out in the remainder
sections. In section 3, we compute the quasinormal modes of the backreacted
holographic superfluid model numerically and use the results to match
holographic and hydrodynamic predictions for the system. In section 4, we
perform a new analysis of the support of the various hydrodynamic modes in
terms of the dual field theory operators. Finally, in section 5, we compute
the excitations of the system and their main properties by considering real
frequency and complex momentum. We conclude in section 6 with a summary of our
results and some comments for the future.

##  2 Superfluids in hydrodynamics and holography

In this section, we first review the main features of the hydrodynamic
description of a conformal superfluid. Next, we introduce the holographic dual
of a superfluid and set up the tools for the analysis performed in the
remaining of this work.

###  2.1 Superfluid hydrodynamics, a brief summary

We shall briefly summarize the main features of superfluid hydrodynamics
doi:10.1063/1.1703944 ; Schmitt:2014eka ; 1974anh…..3…..P ; Nicolis:2011cs ;
Son:2002zn ; Bhattacharya:2011tra , with a focus on the properties of the
different hydrodynamic modes. For connections between superfluid
hydrodynamics, holography and effective field theory see Herzog:2011ec ;
Herzog:2009md ; Yarom:2009uq ; Son:2002zn ; Andersen:2002nd ; Nicolis:2011cs ;
Sonner:2010yx ; Escobedo:2010uv ; Berezhiani:2020umi ; Landry:2020ire ;
Landry:2021kko ; Pajer:2018egx ; Delacretaz:2019brr .

Let us start by considering a relativistic charged fluid with conserved
momentum Kovtun:2012rj ; Son:2006em . The longitudinal spectrum contains two
sets of hydrodynamic modes. First, charge and energy conservation imply the
presence of a thermoelectric diffusive mode $\omega=-iD_{q}k^{2}$ whose
diffusion constant is explicitly given later on in eq. (5). Additionally, we
have a pair of propagating first sound modes $\omega=\pm
c_{1}\,k-\,i\,\Gamma_{1}\,k^{2}$ which result from the interplay of
longitudinal momentum and energy. In a conformal system, their propagation
speed is given by

$c_{1}^{2}=\frac{\partial p}{\partial\epsilon}=\frac{1}{2}$ (1)

where $\epsilon$ and $p$ are the energy density and the thermodynamic
pressure, respectively. Additionally, the attenuation constant reads

$\Gamma_{1}=\frac{\eta}{2\,\chi_{\pi\pi}}\,,$ (2)

with $\eta$ being the shear viscosity and $\chi_{\pi\pi}=\epsilon+p$ the
momentum susceptibility.
Below the critical temperature, a superfluid may be viewed in the two-fluid
picture as composed of two components: a normal component as the
(relativistic) fluid just described; and a condensed or superfluid component
that flows without dissipation. Each component has its own velocity and
density. In the normal phase, the interplay of charge, energy and longitudinal
momentum gives rise to a propagating mode: first sound, and a diffusive mode.
In the superfluid phase $T<T_{c}$, first sound remains in the spectrum, the
diffusive mode gets a purely imaginary gap, and the Goldstone fluctuations,
coupled to those of charge and energy give rise to a new propagating mode:
second sound $\omega=\pm c_{2}\,k-\,i\,\Gamma_{2}\,k^{2}$, known as second
sound doi:10.1063/1.3248499 . Its speed and attenuation constant are given by

$\displaystyle
c_{2}^{2}\,=\,\frac{\rho_{s}\,\xi^{2}}{\left(\mu\,\rho_{n}+sT\right)\,\frac{\partial\xi}{\partial
T}}\,,$ (3)
$\displaystyle\Gamma_{2}\,=\,\frac{\mu\,\rho_{s}}{2\,\chi_{\pi\pi}\,\left(\mu\,\rho_{n}+sT\right)}\,\eta\,+\,\frac{\mu\,\chi_{\pi\pi}}{2\,T^{3}\,\rho^{2}\,(\partial\xi/\partial
T)_{|_{\mu}}}\,\kappa\,+\,\frac{\rho_{s}\,\chi_{\pi\pi}}{2\,\mu\left(\mu\rho_{n}+sT\right)}\,\zeta\,.$
(4)

Here we have defined: the superfluid density $\rho_{s}$, the normal density
$\rho_{n}=\rho-\rho_{s}$, the reduced entropy $\xi=s/\rho$, the diffusivity
$\kappa$ and the superfluid diffusivity $\zeta$.
Notice the only difference with the results presented in Herzog:2011ec , and
valid for $d=4$, is the numerical constant in front of the term proportional
to the shear viscosity $\eta$.

In the normal phase, $\rho_{s}=0$, the first sound remains in the spectrum and
the second sound becomes the thermoelectric diffusive mode Hartnoll:2014lpa ;
Kovtun:2014nsa

$\omega\,=\,-iD_{q}k^{2}\,=-i\,k^{2}\,\frac{\kappa\,\mu\,\chi_{\pi\pi}}{T^{3}\,\rho^{2}\,(\partial\xi/\partial
T)_{|_{\mu}}}\,.$ (5)

Once momentum conservation is broken, the first sound in the normal phase gets
substituted by the energy diffusion mode $\omega=-iD_{e}\,k^{2}$ whose
diffusion constant obeys the standard Einstein relation. At the same time, in
the superfluid phase, the second sound mode is replaced by the so-called
fourth sound first predicted by PhysRev.73.608 ; PhysRev.113.962 where we keep
the normal component stationary. The dispersion relation of fourth sound reads
$\omega=\pm c_{4}\,k-\,i\,\,\Gamma_{4}\,k^{2}$, where the velocity and
attenuation constant are given by

$\displaystyle
c_{4}^{2}\,=\,\frac{\rho_{s}}{\mu\,\left(\frac{\partial_{\rho}}{\partial_{\mu}}\right)}\,,\qquad\Gamma_{4}\,=\,\frac{\kappa}{2\,T\,\chi_{qq}}\,+\,\frac{\rho_{s}}{2\,\mu}\,\zeta\,.$
(6)

For completeness, we list here all the Kubo formulas used to extract the
various transport and thermodynamic parameters entering in the dispersion
relations above

$\displaystyle\eta=\lim_{\omega\rightarrow
0}\frac{1}{\omega}\mathrm{Im}\,\mathcal{G}_{T_{xy}T_{xy}}(\omega,0)\,,$ (7)
$\displaystyle\frac{\kappa}{T}=\lim_{\omega\rightarrow 0}\lim_{k\rightarrow
0}\frac{1}{k}\mathrm{Im}\,\mathcal{G}_{J_{x}J_{t}}(\omega,k)\,,$ (8)
$\displaystyle\zeta=\lim_{\omega\rightarrow
0}\omega\,\mathrm{Im}\,\mathcal{G}_{\phi\phi}(\omega,0)\,,$ (9)

together with the expression for the electric conductivity

$\sigma(\omega)=\frac{i}{\omega}\mathcal{G}_{J_{x}J_{x}}(\omega,0)=\left(\frac{i}{\omega}+\delta(\omega)\right)\left[\frac{\rho_{n}^{2}}{\mu\,\rho_{n}+s\,T}+\frac{\rho_{s}}{\mu}\right]+\sigma_{0}\,.$
(10)

For slow momentum dissipation this two-point function is modified to

$\sigma(\omega)\,=\,\frac{i}{\omega}\,\mathcal{G}_{J_{x}J_{x}}(\omega,0)\,=\,\frac{\rho_{n}^{2}}{\mu\,\rho_{n}+s\,T}\,\frac{1}{\Gamma\,-\,i\,\omega}\,+\frac{i}{\omega}\,\frac{\rho_{s}}{\mu}+\sigma_{0}\,.$
(11)

In this discussion, $\phi$ refers to the superfluid Goldstone mode which is
dual to the phase of the complex bulk scalar (see next section for the
details). Also, $\mathcal{G}_{AB}$ is the retarded Greens function of the
operators $A$ and $B$. These Kubo formulas will be used in section 3 when we
analyze the transport properties of holographic superfluids.

###  2.2 The holographic model

We consider the standard bulk action for a s-wave holographic superfluid
Hartnoll:2008vx ; Hartnoll_2008

$S\,=\,\frac{1}{2\kappa_{4}^{2}}\int\mathrm{d}^{d+1}x\,\sqrt{-g}\left[R-2\Lambda-\frac{1}{4}F_{mn}F^{mn}-|D\psi|^{2}-M^{2}|\psi|^{2}\right]\,,$
(12)

where $\psi$ is a complex scalar, $F\equiv\mathrm{d}A$, and the covariant
derivative is defined as $D_{m}\psi\equiv(\partial_{m}-i\,q\,A_{m})\,\psi$. We
consider the 4-dimensional bulk case $d=3$, and fix $\Lambda=-3$ and
$2\kappa_{4}^{2}=1$. The AdS radius $L$ has also been set to one for
convenience.
The resulting equations of motion are

$\displaystyle
R_{mn}-\frac{1}{2}\left(R-2\Lambda-\frac{1}{4}F^{2}-|D\psi|^{2}-M^{2}|\psi|^{2}\right)g_{mn}$
$\displaystyle=\frac{1}{2}F_{mk}F_{n}^{\;\;k}+\frac{1}{2}\left[D_{m}\psi(D_{n}\psi)^{\star}+D_{n}\psi(D_{m}\psi)^{\star}\right]\,,$
(13a)
$\displaystyle\frac{1}{\sqrt{-g}}\,\partial_{m}\left(\sqrt{-g}\,F^{mn}\right)-iq(\psi^{\star}D^{n}\psi-\psi
D^{n}\psi^{\star})=0\,,$ (13b)
$\displaystyle\frac{1}{\sqrt{-g}}\,D_{m}\left(\sqrt{-g}\,D^{m}\psi\right)-M^{2}\psi=0\,.$
(13c)

We consider the following ansatz for the background

$\displaystyle\mathrm{d}s^{2}=\frac{1}{z^{2}}\left[-u(z)\,\mathrm{d}t^{2}-2\,\mathrm{d}t\mathrm{d}z+g(z)\,(\mathrm{d}x^{2}+\mathrm{d}y^{2})\right],$
(14) $\displaystyle
A=A_{t}(z)\,\mathrm{d}t\,,\qquad\psi=\psi_{1}(z)-i\,\psi_{2}(z)\,,$ (15)

in Eddington-Finkelstein coordinates $\\{t,z,x,y\\}$ with the AdS boundary
located at $z=0$ and the black brane horizon at $z_{h}=1$.
Assuming this ansatz, the temperature and chemical potential of the dual field
theory are given as

$T\,=\,-\frac{u^{\prime}(1)}{4\pi}\,,\qquad\mu\,=\,A_{t}(0)-A_{t}(1)\,.$ (16)

For numerical convenience, we impose $A_{t}(1)=0$. Furthermore, we choose the
value of the background condensate to be real by setting $\psi_{2}^{(2)}$ to
zero (see (42) for the definition of the expectation values)333We add a factor
of 2 to the expectation values of the scalars $\psi$ to account for the non-
canonical normalization in the action (12).

$\langle\mathcal{O}^{\psi}_{\text{cond}}\rangle=2\,\psi^{(2)}_{1}.$ (17)

The analysis in sections 3 and 4 will rely on the study of the fluctuations
about the background discussed above. In particular, we switch on the full set
of longitudinal perturbations
$g_{mn}=\hat{g}_{mn}+\epsilon\,h_{mn}/z^{2},\,A_{m}=\hat{A}_{m}+\epsilon\,a_{m},\,\psi_{i}=\hat{\psi}_{i}+\epsilon\,\delta\psi_{i}$

$\\{h_{tt},\,h_{tx},\,h_{xx},\,h_{yy},\,a_{t},\,a_{x},\,\delta\psi_{1},\,\delta\psi_{2}\\}\,,$
(18)

where we choose the momentum to be aligned with the $x$-direction without loss
of generality. From now on, we will add to all symbols referring to background
quantities a hat $\hat{A}$. Moreover, we assume the radial gauge444Note that
since we work in Eddington-Finkelstein coordinates, the bulk solution for the
scalar field is complex even though we impose the $\psi^{(2)}_{2}$ to be zero.

$h_{zi}=0=a_{z}.$ (19)

We decompose all our fluctuations as

$\mathfrak{f}(t,z,x)\,=\,e^{-i\omega t+ikx}\,\mathfrak{f}(z)\,.$ (20)

The quasinormal modes and all the transport coefficients can be obtained
following standard methods. To compute the thermodynamic derivatives
accurately, for example the one appearing in eq. (3), we re-write the
derivative in terms of dimensionless quantities; we then compute the
background on a Chebychev grid in that variable (for example $T/\mu$) and use
the spectral differentiation matrices to compute the thermodynamic derivative
(see appendix B.4 in Ammon:2020rvg for more details).

In order to introduce momentum dissipation and study the properties of fourth
sound, we supplement the action (12) by a linear axion term Andrade:2013gsa

$S_{2}\,=\,-\frac{1}{2\kappa_{4}^{2}}\,\int\mathrm{d}^{d+1}x\,\sqrt{-g}\,\partial_{m}\phi^{I}\partial^{m}\phi^{I}$
(21)

where the index $I$ runs only along the spatial directions $I=x,y$, and the
background profile for the scalars is

$\hat{\phi}^{I}\,=\,\alpha\,x^{I}\,.$ (22)

In the dual field theory, this corresponds to switching on an operator which
breaks translational invariance explicitly. For more details about this model
we refer to the recent review Baggioli:2021xuv . In presence of momentum
dissipation, we add the fluctuation of the longitudinal axion field
$\delta\phi^{x}$ to the fluctuations in eq. (18).

Finally, throughout all the manuscript we set $q=1$ and $M^{2}=-2$.

##  3 Transport & hydrodynamic modes

In this section, we explicitly verify that the hydrodynamic framework
described in the previous section (see Herzog:2011ec for the original
reference, albeit in one dimension more) is indeed the correct hydrodynamic
description of the holographic superfluid model Hartnoll:2008vx .
Surprisingly, even though partial checks may be found in the literature
Herzog:2011ec ; Yarom:2009uq ; Herzog:2009md , a complete matching of the two
pictures considering the full backreaction Hartnoll:2008kx in the
gravitational model has not appeared yet.

Figure 1: The normal and superfluid densities $\rho_{n}$ and $\rho_{s}$ as
function of the temperature. Left: Translational invariant theory. Right:
Model with broken translations, $\alpha/\mu=1/10$.

We start by plotting the normal and superfluid densities as function of
$T/T_{c}$ in the superfluid phase. The normal and superfluid density,
respectively, may be extracted numerically from eq. (10) (translationally
invariant case) and eq. (11) (broken translations). While in the latter case,
the superfluid density is simply given by the $i/\omega$ pole in the
conductivity, we have to use that $\rho_{n}=\rho-\rho_{s}$ and solve eq. (10)
for $\rho_{s}$ in the translationally invariant case. The results are shown in
fig. 1 for both the translational invariant model and the model where
translations are broken explicitly using axion fields Baggioli:2021xuv ;
Baggioli:2016rdj . As it should, the superfluid density vanishes at the
critical temperature for both models. In the translational invariant case, we
also observe that the normal density vanishes at zero temperature. This last
feature holds in superfluid Helium and BCS superconductors, while a residual
normal density at zero temperature has been reported in some high-Tc
superconductors PMID:27535534 . Interestingly, Gouteraux:2019kuy has shown
that in quantum critical superfluids the nature of the IR fixed point of the
theory is crucial for the fate of the normal density: it vanishes for fixed
points with an emergent Lorentz symmetry as in our translation invariant case
Horowitz:2009ij , while some quantum critical points of the hyperscaling-
Lifshitz type can support a nonvanishing normal density Gouteraux:2020asq .
The model breaking translations (right panel of fig. 1) shows a slower decay
of the normal density for $T\to 0$ and a fit of the numerical data shows a
scaling of $\rho_{n}/\rho\sim 0.05265+3.14843(T/\mu)^{1.2583}$ for
$T/\mu\in[5\cdot 10^{-7},5\cdot 10^{-6}]$ indicating a non-zero normal density
in the zero temperature limit. It would be interesting to characterize the IR
geometry of our model and analyze it in the light of the results of
Gouteraux:2019kuy . We leave a more accurate analysis of this feature for the
future.

We then move to the study of the sound modes. On the gravity side these
correspond to quasinormal modes of the black hole geometry. As we review in
section 4.1, with further technical details in appendix D, we obtain them by
numerically solving the linearized equations of motion for the fluctuations
(18). This will allow us to compute the speed and attenuation of the different
sound modes and compare them against the hydrodynamic predictions of section
2.1. The propagation speed and attenuation constant of the first sound mode
are shown in fig. 2. The propagation speed is fixed by conformal symmetry to
the expected value $c_{1}^{2}=1/2$ and it is continuous across the phase
transition at $T=T_{c}$. Its dimensionless attenuation constant $\Gamma_{1}T$
is shown in the right panel of fig. 2 and it is in perfect agreement with the
hydrodynamic formula in eq. (2) in which all the parameters are independently
computed using thermodynamic relations and transport Kubo formulas. One should
notice that the dimensionless attenuation constant vanishes at zero
temperature, indicating that at $T=0$ first sound is an infinite living
excitation. This is tantamount to saying that the shear viscosity of the
system vanishes at $T=0$. Finally, we observe that the attenuation constant of
first sound is continuous at $T=T_{c}$ even though its first derivative is
not.

Figure 2: The propagation speed and the attenuation constant of first sound as
function of $T/T_{c}$ above and below the critical point. The red lines are
the hydrodynamic formulas (1), (2) and the colored circles the numerical data
extracted from the QNMs. Blue color indicates the normal phase. The plot for
the speed is not extended into the normal phase since the speed is trivially
constant in the full range of $T$ due to the conformal symmetry.

We perform an analogous analysis for second sound in fig. 3. We remind the
reader that the second sound mode exists only in the superfluid phase below
the critical temperature $T<T_{c}$. Interestingly, we observe a very non-
monotonic behavior for the second sound speed which vanishes both at the
critical temperature $T=T_{c}$ and at zero temperature $T=0$, with a maximum
around $T=0.6\,T_{c}$. Apparently, this behavior is not universal and it
rather depends on the charge and the conformal dimension of the dual scalar
operator Herzog:2009md . Again the hydrodynamic formula (3) fits the numerical
data very well.

The attenuation constant of second sound is more interesting. First, it
vanishes at zero temperature. Second, across the phase transition it is
connected to the imaginary part of the gapped scalar modes present in the
normal phase hydroSC and responsible for the superfluid instability at
$T=T_{c}$. Nevertheless, the data display a clear and sharp jump around
$T=T_{c}$ where our hydrodynamic description is no longer applicable as is
evident in the right panel of fig. 3. Moreover, the attenuation constant of
second sound is not continuously connected to the diffusion constant which we
present in the inset in the same figure 3. Our numerical observation was
confirmed by the authors of Donos:2021pkk in a simpler setup (without chemical
potential) using analytic arguments.

Figure 3: The propagation speed and the attenuation constant of second sound
as function of $T/T_{c}$ above and below the critical point. The red lines are
the hydrodynamic formulas (3), (4) and the colored circles the numerical data
extracted from the QNMs. Blue color indicates the normal phase. The inset in
the right panel shows the diffusion constant of the purely diffusive mode in
the normal phase together with its hydrodynamic expression.

Finally, when momentum is dissipated and the normal component is prevented
from flowing, a new excitation called fourth sound appears in the superfluid
PhysRev.113.962 ; PhysRev.73.608 . We obtained the speed of propagation and
the attenuation constant numerically in fig. 4 and compared them with the
hydrodynamic expressions in eq. (6). As expected the speed of fourth sound
interpolates between the speed of second sound at the critical temperature
$T_{c}$ to the speed of first sound for $T\rightarrow 0$ PhysRev.113.962 ;
PhysRev.137.A1383 . Additionally, notice how the attenuation constant of
fourth sound is much larger than those of first and second sounds. Finally, we
observe that the speed and attenuation constant of fourth sound have the same
qualitative behavior as those of the sound mode obtained in the probe limit in
Amado:2013xya . We do expect the two to coincide in the limit of very fast
momentum dissipation, $\alpha/T\gg 1$.

To conclude this first part, let us iterate that the hydrodynamic description
presented in the previous section perfectly matches the numerical results
obtained from the holographic superfluid model both in the normal and
superfluid phases with and without momentum dissipation.

Figure 4: The propagation speed and the attenuation constant of fourth sound
in function of $T/T_{c}$ below the critical temperature for $\alpha/\mu=0.1$.
The red lines are the hydrodynamic formulas (6) and the colored circles the
numerical data extracted from the QNMs. The dashed lines in the left panel are
the speed of first and second sound for $\alpha/\mu=0$.

##  4 Tomography of superfluid sound modes

After verifying the hydrodynamic theory in our holographic model explicitly,
we now move on to investigate the so-called support of the hydrodynamic modes.

###  4.1 Theory

A quasinormal mode is a solution to the linearized equations of motion in a
black hole background with in-going boundary conditions on the horizon and
normalizable ones at the boundary of AdS Berti:2009kk ; Nollert_1999 . In
general, such a solution will consist of a collection of fields which we
denote by $\Phi^{I}$. For a single QNM all these fields oscillate with the
same frequency and decay exponentially with the same decay rate. At late times
and close to equilibrium, we may write a generic QNM solution as a
superposition of modes

$\phi^{I}(z,t,\vec{x})=\int\frac{\mathrm{d}^{d}k}{(2\pi)^{d}}\sum_{n}A^{I}_{n}(z,\vec{k})e^{-i\Omega_{n}(\vec{k})t}e^{-\Gamma_{n}(\vec{k})t}e^{i\vec{k}\vec{x}}\,.$
(23)

The bulk wave functions of a QNM $A^{I}_{n}(\vec{k},z)$ do not contain the
leading (non–normalizable) terms of the asymptotic boundary condition. We can
therefore multiply the field $\phi^{I}$ with the appropriate power of $z$ such
that $A_{n}^{I}(0,\vec{k})$ is the coefficient of the non-normalizable mode.
According to the holographic dictionary the boundary value of $\Phi_{n}^{I}$
encodes the expectation value of an operator $\langle O^{I}(t,\vec{x})\rangle$
(see appendix A). Since the solution (23) evaluated at the boundary represents
a propagating and attenuated sound wave in the dual field theory, it is
instructive to think about the quasi-normal mode as corresponding to a
hydrodynamic sound wave. One should notice that sound waves consist of more
than one perturbation, for example the perturbations of the energy density
$\epsilon$ and the pressure $p$ (as is clear from the thermodynamic formula
for the speed of sound $c_{s}^{2}=\partial p/\partial\epsilon$). Hence, we
expect the QNM solution dual to a sound mode to consist of more than one field
$\Phi^{I}$ and taking the form of eq. (23).

To simplify the analysis further, we assume that rotational symmetry is not
broken. Then the quasinormal frequencies are function of $|\vec{k}|^{2}$ only.
From now on we write $k$ for $|\vec{k}|$. Since the collection $\Phi^{I}$
includes vector fields and metric components the dependence of the amplitudes
may be more complicated in general.

Demanding that the fields $\Phi^{I}$ be real results in the following
condition for the quasi-normal frequencies and amplitudes

$(\phi^{I})^{*}=\phi^{I}=\int\frac{\mathrm{d}^{d}k}{(2\pi)^{d}}\sum_{m}A_{m}^{I*}(z,-\vec{k})e^{i\Omega_{m}(k)t}e^{-\Gamma_{m}(k)t}e^{i\vec{k}\vec{x}}\,.$
(24)

Combining this constraint with (23) leads to the condition that the
quasinormal fequencies come in pairs $(n,m)$ with

$\displaystyle\Omega_{m}=-\Omega_{n}\,,\qquad\Gamma_{m}=\Gamma_{n}\,,\qquad
A_{m}^{I*}(z,-\vec{k})=A_{n}(z,\vec{k})\,,$

or are purely damped with vanishing $\Omega_{n}$ and
$A_{n}^{I*}(z,-\vec{k})=A_{n}^{I}(z,\vec{k})$.

In view of this structure, the minimal way of constructing a real QNM solution
for $\Phi^{I}$ is to excite the mode $n$ at wave number $\vec{k}$ and the
mirror mode $m$ with $\Omega_{m}=-\Omega_{n}$ at wave number $-\vec{k}$.

$\Phi^{I}(z,t,\vec{x})=e^{-\Gamma_{n}t}\left[A^{I}(z,\vec{k})\,e^{-i\Omega
t+i\vec{k}\vec{x}}+A^{I}(z,\vec{k})^{*}e^{i\Omega
t-i\vec{k}\vec{x}}\right]\,.$ (25)

If we evaluate this at the boundary, we get the space-time evolution of the
set of operators $O^{I}$ as

$\langle O^{I}(t,\vec{x})\rangle=e^{-\Gamma_{n}t}\langle
O^{I}(0,\vec{k})\rangle\cos(\Omega_{n}t-\vec{k}\vec{x}+\alpha^{I}_{n})\,,$
(26)

where we take $\langle O^{I}(0,\vec{k})\rangle$ to be positive since the phase
$\alpha^{I}_{n}$ accounts for the sign of $\langle O^{I}(t,\vec{x})\rangle$.
We need to take into account, however, that (26) represents a solution to the
linearized equations of motion. This means that we can multiply it with an
arbitrary complex number and obtain another valid solution. Therefore only
relative amplitudes and phases contain physical information. Moreover, since a
generic solution will involve operators $O^{I}$ with different dimensions, we
should only compare dimensionless expectation values that can be easily
obtained by normalizing with the corresponding power of the temperature,
namely

$\mathcal{O}^{I}(t,\vec{x})=T^{-\Delta}\langle O^{I}(t,\vec{x})\rangle\,.$
(27)

We also denote the dimensionless absolute values of the amplitudes as
$R_{n}^{I}(\vec{k})=T^{-\Delta}|\langle O^{I}_{n}(0,\vec{k})\rangle|$ (note
that all $R^{I}_{n}\in\mathbb{R}^{+}_{0}$). We can then write our normalized
operators as

$\mathcal{O}^{I}(t,\vec{x})=e^{-\Gamma_{n}t}R_{n}^{I}(\vec{k})\cos(\Omega_{n}t-\vec{k}\vec{x}+\alpha^{I}_{n})\,.$
(28)

A measure of the relative support on the different operators $O^{I}$ of a
given QNM is thus given by the ratio of their amplitudes $R_{n}^{I}$. Namely,
$r_{n}^{JI}=R^{J}_{n}/R^{I}_{n}$ measures the relative support of a mode on
the operator $O^{I}$ with respect to the operator $O^{J}$. In addition, the
relative phases $\phi^{IJ}=\alpha^{I}-\alpha^{J}$ encode the phase delay
between the excitations of operators $O^{I}$ and $O^{J}$. 555We could also
consider taking the logarithm of the ratios $r^{IJ}_{n}$ and defining
$\tilde{r}^{IJ}_{n}=\log(R^{J}_{n})-\log(R^{I}_{n})$. This has the advantage
that the $\tilde{r}^{JI}_{n}$ are antisymmetric in $JI$ as are the relative
phases $\phi^{IJ}$.

A useful way of thinking about these ratios and relative phases of the
expectation values, would be to consider a (thought) experiment where a device
couples directly to the operator $O^{I}$ and only to that one. This device
would act as a source of $O^{I}$ alone among all other operators. But the
operator mixing inherent to the system would result in a state where operators
other than $O^{I}$ are also excited. Then, the ratios $r_{n}^{JI}$ and phases
$\phi^{JI}$ will precisely characterize the response of the operator $O^{J}$
when a source is switched on for the operator $O^{I}$.

A (thorough) study of relative amplitudes and phases of QNMs is lacking in the
literature, although closely related analyses of the residues of holographic
Green’s functions were carried out in Amado:2007yr ; Amado:2008ji ;
Davison:2011uk . We will do this here for a particularly interesting example,
namely the first and second sound modes of a holographic superconductor.
Notice that a similar study of relative phases and amplitudes for sound modes
of a weakly coupled superfluid has been carried out in Alford:2013koa , and we
will eventually compare our results to it.

Finally, we shall comment on the difference between the ratios of amplitudes
that we will determine here, and the residues of a certain mode $\omega(k)$.
The residue of a mode quantifies the contribution of that mode to a given
response: i.e. the weight of its corresponding pole in the spectral
representation of a given Green’s function. Instead, the amplitudes we discuss
here represent the weight of a given operator in the collective excitation
corresponding to a specific mode. In a sense, these observables provide
complementary information about the system.

###  4.2 Results

For concreteness, we study the longitudinal spectrum of excitations, to which
first and second sounds belong. This amounts to solving for the set of
fluctuations in eq. (18), as we describe at the end of section 2.2. Moreover,
we will consider solutions corresponding to QNMs: these are dual to
configurations in which all the sources of the dual operators vanish, while
(some of) the expectation values do not. In appendix A, we work out the map
between fluctuations and dual operators. As a result, in the dual theory we
will be dealing with the following set of operators corresponding to the
expectation values of the fluctuations

$\displaystyle\leavevmode\hbox to4.71pt{\vbox
to4.71pt{\pgfpicture\makeatletter\hbox{\hskip 2.35277pt\lower-2.35277pt\hbox
to0.0pt{\pgfsys@beginscope\pgfsys@invoke{
}\definecolor{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@rgb@stroke{0}{0}{0}\pgfsys@invoke{
}\pgfsys@color@rgb@fill{0}{0}{0}\pgfsys@invoke{
}\pgfsys@setlinewidth{0.4pt}\pgfsys@invoke{ }\nullfont\hbox
to0.0pt{\pgfsys@beginscope\pgfsys@invoke{
}{}{{}}{}{{{}}{}{}{}{}{}{}{}{}}\pgfsys@beginscope\pgfsys@invoke{
}\color[rgb]{1,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{1,0,0}\pgfsys@color@rgb@stroke{1}{0}{0}\pgfsys@invoke{
}\pgfsys@color@rgb@fill{1}{0}{0}\pgfsys@invoke{
}\definecolor[named]{pgffillcolor}{rgb}{1,0,0}\definecolor[named]{pgffillcolor}{rgb}{1,0,0}\pgfsys@color@rgb@fill{1}{0}{0}\pgfsys@invoke{
}{}\pgfsys@moveto{0.0pt}{0.0pt}\pgfsys@moveto{2.15277pt}{0.0pt}\pgfsys@curveto{2.15277pt}{1.18895pt}{1.18895pt}{2.15277pt}{0.0pt}{2.15277pt}\pgfsys@curveto{-1.18895pt}{2.15277pt}{-2.15277pt}{1.18895pt}{-2.15277pt}{0.0pt}\pgfsys@curveto{-2.15277pt}{-1.18895pt}{-1.18895pt}{-2.15277pt}{0.0pt}{-2.15277pt}\pgfsys@curveto{1.18895pt}{-2.15277pt}{2.15277pt}{-1.18895pt}{2.15277pt}{0.0pt}\pgfsys@closepath\pgfsys@moveto{0.0pt}{0.0pt}\pgfsys@fillstroke\pgfsys@invoke{
} \pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope
\pgfsys@invoke{\lxSVG@closescope
}\pgfsys@endscope{}{}{}\hss}\pgfsys@discardpath\pgfsys@invoke{\lxSVG@closescope
}\pgfsys@endscope\hss}}\lxSVG@closescope\endpgfpicture}}\,\,\langle\delta
T^{tt}\rangle\,\,:\,\,\,\text{energy density},$ $\displaystyle\leavevmode\hbox
to4.71pt{\vbox to4.71pt{\pgfpicture\makeatletter\hbox{\hskip
2.35277pt\lower-2.35277pt\hbox to0.0pt{\pgfsys@beginscope\pgfsys@invoke{
}\definecolor{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@rgb@stroke{0}{0}{0}\pgfsys@invoke{
}\pgfsys@color@rgb@fill{0}{0}{0}\pgfsys@invoke{
}\pgfsys@setlinewidth{0.4pt}\pgfsys@invoke{ }\nullfont\hbox
to0.0pt{\pgfsys@beginscope\pgfsys@invoke{
}{}{{}}{}{{{}}{}{}{}{}{}{}{}{}}\pgfsys@beginscope\pgfsys@invoke{
}\color[rgb]{0,0,1}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,1}\pgfsys@color@rgb@stroke{0}{0}{1}\pgfsys@invoke{
}\pgfsys@color@rgb@fill{0}{0}{1}\pgfsys@invoke{
}\definecolor[named]{pgffillcolor}{rgb}{0,0,1}\definecolor[named]{pgffillcolor}{rgb}{0,0,1}\pgfsys@color@rgb@fill{0}{0}{1}\pgfsys@invoke{
}{}\pgfsys@moveto{0.0pt}{0.0pt}\pgfsys@moveto{2.15277pt}{0.0pt}\pgfsys@curveto{2.15277pt}{1.18895pt}{1.18895pt}{2.15277pt}{0.0pt}{2.15277pt}\pgfsys@curveto{-1.18895pt}{2.15277pt}{-2.15277pt}{1.18895pt}{-2.15277pt}{0.0pt}\pgfsys@curveto{-2.15277pt}{-1.18895pt}{-1.18895pt}{-2.15277pt}{0.0pt}{-2.15277pt}\pgfsys@curveto{1.18895pt}{-2.15277pt}{2.15277pt}{-1.18895pt}{2.15277pt}{0.0pt}\pgfsys@closepath\pgfsys@moveto{0.0pt}{0.0pt}\pgfsys@fillstroke\pgfsys@invoke{
} \pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope
\pgfsys@invoke{\lxSVG@closescope
}\pgfsys@endscope{}{}{}\hss}\pgfsys@discardpath\pgfsys@invoke{\lxSVG@closescope
}\pgfsys@endscope\hss}}\lxSVG@closescope\endpgfpicture}}\,\,\langle\delta
T^{tx}\rangle\,:\,\,\,\text{longitudinal momentum},$
$\displaystyle\leavevmode\hbox to4.71pt{\vbox
to4.71pt{\pgfpicture\makeatletter\hbox{\hskip 2.35277pt\lower-2.35277pt\hbox
to0.0pt{\pgfsys@beginscope\pgfsys@invoke{
}\definecolor{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@rgb@stroke{0}{0}{0}\pgfsys@invoke{
}\pgfsys@color@rgb@fill{0}{0}{0}\pgfsys@invoke{
}\pgfsys@setlinewidth{0.4pt}\pgfsys@invoke{ }\nullfont\hbox
to0.0pt{\pgfsys@beginscope\pgfsys@invoke{
}{}{{}}{}{{{}}{}{}{}{}{}{}{}{}}\pgfsys@beginscope\pgfsys@invoke{
}\color[rgb]{0,1,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,1,0}\pgfsys@color@rgb@stroke{0}{1}{0}\pgfsys@invoke{
}\pgfsys@color@rgb@fill{0}{1}{0}\pgfsys@invoke{
}\definecolor[named]{pgffillcolor}{rgb}{0,1,0}\definecolor[named]{pgffillcolor}{rgb}{0,1,0}\pgfsys@color@rgb@fill{0}{1}{0}\pgfsys@invoke{
}{}\pgfsys@moveto{0.0pt}{0.0pt}\pgfsys@moveto{2.15277pt}{0.0pt}\pgfsys@curveto{2.15277pt}{1.18895pt}{1.18895pt}{2.15277pt}{0.0pt}{2.15277pt}\pgfsys@curveto{-1.18895pt}{2.15277pt}{-2.15277pt}{1.18895pt}{-2.15277pt}{0.0pt}\pgfsys@curveto{-2.15277pt}{-1.18895pt}{-1.18895pt}{-2.15277pt}{0.0pt}{-2.15277pt}\pgfsys@curveto{1.18895pt}{-2.15277pt}{2.15277pt}{-1.18895pt}{2.15277pt}{0.0pt}\pgfsys@closepath\pgfsys@moveto{0.0pt}{0.0pt}\pgfsys@fillstroke\pgfsys@invoke{
} \pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope
\pgfsys@invoke{\lxSVG@closescope
}\pgfsys@endscope{}{}{}\hss}\pgfsys@discardpath\pgfsys@invoke{\lxSVG@closescope
}\pgfsys@endscope\hss}}\lxSVG@closescope\endpgfpicture}}\,\,\langle\delta
T^{xx}\rangle:\,\,\,\text{pressure in the $x$ -- direction},$
$\displaystyle\leavevmode\hbox to4.71pt{\vbox
to4.71pt{\pgfpicture\makeatletter\hbox{\hskip 2.35277pt\lower-2.35277pt\hbox
to0.0pt{\pgfsys@beginscope\pgfsys@invoke{
}\definecolor{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@rgb@stroke{0}{0}{0}\pgfsys@invoke{
}\pgfsys@color@rgb@fill{0}{0}{0}\pgfsys@invoke{
}\pgfsys@setlinewidth{0.4pt}\pgfsys@invoke{ }\nullfont\hbox
to0.0pt{\pgfsys@beginscope\pgfsys@invoke{
}{}{{}}{}{{{}}{}{}{}{}{}{}{}{}}\pgfsys@beginscope\pgfsys@invoke{
}\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@invoke{
}\pgfsys@color@gray@fill{0}\pgfsys@invoke{
}\definecolor[named]{pgffillcolor}{rgb}{0,0,0}\definecolor[named]{pgffillcolor}{rgb}{0,0,0}\pgfsys@color@gray@fill{0}\pgfsys@invoke{
}{}\pgfsys@moveto{0.0pt}{0.0pt}\pgfsys@moveto{2.15277pt}{0.0pt}\pgfsys@curveto{2.15277pt}{1.18895pt}{1.18895pt}{2.15277pt}{0.0pt}{2.15277pt}\pgfsys@curveto{-1.18895pt}{2.15277pt}{-2.15277pt}{1.18895pt}{-2.15277pt}{0.0pt}\pgfsys@curveto{-2.15277pt}{-1.18895pt}{-1.18895pt}{-2.15277pt}{0.0pt}{-2.15277pt}\pgfsys@curveto{1.18895pt}{-2.15277pt}{2.15277pt}{-1.18895pt}{2.15277pt}{0.0pt}\pgfsys@closepath\pgfsys@moveto{0.0pt}{0.0pt}\pgfsys@fillstroke\pgfsys@invoke{
} \pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope
\pgfsys@invoke{\lxSVG@closescope
}\pgfsys@endscope{}{}{}\hss}\pgfsys@discardpath\pgfsys@invoke{\lxSVG@closescope
}\pgfsys@endscope\hss}}\lxSVG@closescope\endpgfpicture}}\,\,\langle\delta
T^{yy}\rangle:\,\,\,\text{pressure in the $y$ -- direction},$
$\displaystyle\leavevmode\hbox to4.71pt{\vbox
to4.71pt{\pgfpicture\makeatletter\hbox{\hskip 2.35277pt\lower-2.35277pt\hbox
to0.0pt{\pgfsys@beginscope\pgfsys@invoke{
}\definecolor{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@rgb@stroke{0}{0}{0}\pgfsys@invoke{
}\pgfsys@color@rgb@fill{0}{0}{0}\pgfsys@invoke{
}\pgfsys@setlinewidth{0.4pt}\pgfsys@invoke{ }\nullfont\hbox
to0.0pt{\pgfsys@beginscope\pgfsys@invoke{
}{}{{}}{}{{{}}{}{}{}{}{}{}{}{}}\pgfsys@beginscope\pgfsys@invoke{
}\color[rgb]{0.6,0.4,0.2}\definecolor[named]{pgfstrokecolor}{rgb}{0.6,0.4,0.2}\pgfsys@color@rgb@stroke{0.6}{0.4}{0.2}\pgfsys@invoke{
}\pgfsys@color@rgb@fill{0.6}{0.4}{0.2}\pgfsys@invoke{
}\definecolor[named]{pgffillcolor}{rgb}{0.6,0.4,0.2}\definecolor[named]{pgffillcolor}{rgb}{0.6,0.4,0.2}\pgfsys@color@rgb@fill{0.6}{0.4}{0.2}\pgfsys@invoke{
}{}\pgfsys@moveto{0.0pt}{0.0pt}\pgfsys@moveto{2.15277pt}{0.0pt}\pgfsys@curveto{2.15277pt}{1.18895pt}{1.18895pt}{2.15277pt}{0.0pt}{2.15277pt}\pgfsys@curveto{-1.18895pt}{2.15277pt}{-2.15277pt}{1.18895pt}{-2.15277pt}{0.0pt}\pgfsys@curveto{-2.15277pt}{-1.18895pt}{-1.18895pt}{-2.15277pt}{0.0pt}{-2.15277pt}\pgfsys@curveto{1.18895pt}{-2.15277pt}{2.15277pt}{-1.18895pt}{2.15277pt}{0.0pt}\pgfsys@closepath\pgfsys@moveto{0.0pt}{0.0pt}\pgfsys@fillstroke\pgfsys@invoke{
} \pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope
\pgfsys@invoke{\lxSVG@closescope
}\pgfsys@endscope{}{}{}\hss}\pgfsys@discardpath\pgfsys@invoke{\lxSVG@closescope
}\pgfsys@endscope\hss}}\lxSVG@closescope\endpgfpicture}}\,\,\langle\delta
J^{t}\rangle\,\,\,:\,\,\,\text{charge density},$
$\displaystyle\leavevmode\hbox to4.71pt{\vbox
to4.71pt{\pgfpicture\makeatletter\hbox{\hskip 2.35277pt\lower-2.35277pt\hbox
to0.0pt{\pgfsys@beginscope\pgfsys@invoke{
}\definecolor{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@rgb@stroke{0}{0}{0}\pgfsys@invoke{
}\pgfsys@color@rgb@fill{0}{0}{0}\pgfsys@invoke{
}\pgfsys@setlinewidth{0.4pt}\pgfsys@invoke{ }\nullfont\hbox
to0.0pt{\pgfsys@beginscope\pgfsys@invoke{
}{}{{}}{}{{{}}{}{}{}{}{}{}{}{}}\pgfsys@beginscope\pgfsys@invoke{
}\color[rgb]{1,0,1}\definecolor[named]{pgfstrokecolor}{rgb}{1,0,1}\pgfsys@color@rgb@stroke{1}{0}{1}\pgfsys@invoke{
}\pgfsys@color@rgb@fill{1}{0}{1}\pgfsys@invoke{
}\definecolor[named]{pgffillcolor}{rgb}{1,0,1}\definecolor[named]{pgffillcolor}{rgb}{1,0,1}\pgfsys@color@rgb@fill{1}{0}{1}\pgfsys@invoke{
}{}\pgfsys@moveto{0.0pt}{0.0pt}\pgfsys@moveto{2.15277pt}{0.0pt}\pgfsys@curveto{2.15277pt}{1.18895pt}{1.18895pt}{2.15277pt}{0.0pt}{2.15277pt}\pgfsys@curveto{-1.18895pt}{2.15277pt}{-2.15277pt}{1.18895pt}{-2.15277pt}{0.0pt}\pgfsys@curveto{-2.15277pt}{-1.18895pt}{-1.18895pt}{-2.15277pt}{0.0pt}{-2.15277pt}\pgfsys@curveto{1.18895pt}{-2.15277pt}{2.15277pt}{-1.18895pt}{2.15277pt}{0.0pt}\pgfsys@closepath\pgfsys@moveto{0.0pt}{0.0pt}\pgfsys@fillstroke\pgfsys@invoke{
} \pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope
\pgfsys@invoke{\lxSVG@closescope
}\pgfsys@endscope{}{}{}\hss}\pgfsys@discardpath\pgfsys@invoke{\lxSVG@closescope
}\pgfsys@endscope\hss}}\lxSVG@closescope\endpgfpicture}}\,\,\langle\delta
J^{x}\rangle\,\,:\,\,\,\text{longitudinal current},$
$\displaystyle\leavevmode\hbox to4.71pt{\vbox
to4.71pt{\pgfpicture\makeatletter\hbox{\hskip 2.35277pt\lower-2.35277pt\hbox
to0.0pt{\pgfsys@beginscope\pgfsys@invoke{
}\definecolor{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@rgb@stroke{0}{0}{0}\pgfsys@invoke{
}\pgfsys@color@rgb@fill{0}{0}{0}\pgfsys@invoke{
}\pgfsys@setlinewidth{0.4pt}\pgfsys@invoke{ }\nullfont\hbox
to0.0pt{\pgfsys@beginscope\pgfsys@invoke{
}{}{{}}{}{{{}}{}{}{}{}{}{}{}{}}\pgfsys@beginscope\pgfsys@invoke{
}\color[rgb]{1,0.5,0}\definecolor[named]{pgfstrokecolor}{rgb}{1,0.5,0}\pgfsys@color@rgb@stroke{1}{0.5}{0}\pgfsys@invoke{
}\pgfsys@color@rgb@fill{1}{0.5}{0}\pgfsys@invoke{
}\definecolor[named]{pgffillcolor}{rgb}{1,0.5,0}\definecolor[named]{pgffillcolor}{rgb}{1,0.5,0}\pgfsys@color@rgb@fill{1}{0.5}{0}\pgfsys@invoke{
}{}\pgfsys@moveto{0.0pt}{0.0pt}\pgfsys@moveto{2.15277pt}{0.0pt}\pgfsys@curveto{2.15277pt}{1.18895pt}{1.18895pt}{2.15277pt}{0.0pt}{2.15277pt}\pgfsys@curveto{-1.18895pt}{2.15277pt}{-2.15277pt}{1.18895pt}{-2.15277pt}{0.0pt}\pgfsys@curveto{-2.15277pt}{-1.18895pt}{-1.18895pt}{-2.15277pt}{0.0pt}{-2.15277pt}\pgfsys@curveto{1.18895pt}{-2.15277pt}{2.15277pt}{-1.18895pt}{2.15277pt}{0.0pt}\pgfsys@closepath\pgfsys@moveto{0.0pt}{0.0pt}\pgfsys@fillstroke\pgfsys@invoke{
} \pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope
\pgfsys@invoke{\lxSVG@closescope
}\pgfsys@endscope{}{}{}\hss}\pgfsys@discardpath\pgfsys@invoke{\lxSVG@closescope
}\pgfsys@endscope\hss}}\lxSVG@closescope\endpgfpicture}}\,\,\langle\delta\mathrm{H}\rangle\,\,\,\,:\,\,\,\text{Higgs
mode (condensate fluctuations)},$ $\displaystyle\leavevmode\hbox
to4.71pt{\vbox to4.71pt{\pgfpicture\makeatletter\hbox{\hskip
2.35277pt\lower-2.35277pt\hbox to0.0pt{\pgfsys@beginscope\pgfsys@invoke{
}\definecolor{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@rgb@stroke{0}{0}{0}\pgfsys@invoke{
}\pgfsys@color@rgb@fill{0}{0}{0}\pgfsys@invoke{
}\pgfsys@setlinewidth{0.4pt}\pgfsys@invoke{ }\nullfont\hbox
to0.0pt{\pgfsys@beginscope\pgfsys@invoke{
}{}{{}}{}{{{}}{}{}{}{}{}{}{}{}}\pgfsys@beginscope\pgfsys@invoke{
}\color[rgb]{0.666667,0.666667,1}\definecolor[named]{pgfstrokecolor}{rgb}{0.666667,0.666667,1}\pgfsys@color@rgb@stroke{0.666667}{0.666667}{1}\pgfsys@invoke{
}\pgfsys@color@rgb@fill{0.666667}{0.666667}{1}\pgfsys@invoke{
}\definecolor[named]{pgffillcolor}{rgb}{0.666667,0.666667,1}\definecolor[named]{pgffillcolor}{rgb}{0.666667,0.666667,1}\pgfsys@color@rgb@fill{0.666667}{0.666667}{1}\pgfsys@invoke{
}{}\pgfsys@moveto{0.0pt}{0.0pt}\pgfsys@moveto{2.15277pt}{0.0pt}\pgfsys@curveto{2.15277pt}{1.18895pt}{1.18895pt}{2.15277pt}{0.0pt}{2.15277pt}\pgfsys@curveto{-1.18895pt}{2.15277pt}{-2.15277pt}{1.18895pt}{-2.15277pt}{0.0pt}\pgfsys@curveto{-2.15277pt}{-1.18895pt}{-1.18895pt}{-2.15277pt}{0.0pt}{-2.15277pt}\pgfsys@curveto{1.18895pt}{-2.15277pt}{2.15277pt}{-1.18895pt}{2.15277pt}{0.0pt}\pgfsys@closepath\pgfsys@moveto{0.0pt}{0.0pt}\pgfsys@fillstroke\pgfsys@invoke{
} \pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope
\pgfsys@invoke{\lxSVG@closescope
}\pgfsys@endscope{}{}{}\hss}\pgfsys@discardpath\pgfsys@invoke{\lxSVG@closescope
}\pgfsys@endscope\hss}}\lxSVG@closescope\endpgfpicture}}\,\,\langle\,\partial_{t}\varphi\rangle\,\,:\,\,\,\text{time
derivative of the phase fluctuation (Goldstone)}\,,$
$\displaystyle\leavevmode\hbox to4.71pt{\vbox
to4.71pt{\pgfpicture\makeatletter\hbox{\hskip 2.35277pt\lower-2.35277pt\hbox
to0.0pt{\pgfsys@beginscope\pgfsys@invoke{
}\definecolor{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@rgb@stroke{0}{0}{0}\pgfsys@invoke{
}\pgfsys@color@rgb@fill{0}{0}{0}\pgfsys@invoke{
}\pgfsys@setlinewidth{0.4pt}\pgfsys@invoke{ }\nullfont\hbox
to0.0pt{\pgfsys@beginscope\pgfsys@invoke{
}{}{{}}{}{{{}}{}{}{}{}{}{}{}{}}\pgfsys@beginscope\pgfsys@invoke{
}\color[rgb]{1,1,0}\definecolor[named]{pgfstrokecolor}{rgb}{1,1,0}\pgfsys@color@cmyk@stroke{0}{0}{1}{0}\pgfsys@invoke{
}\pgfsys@color@cmyk@fill{0}{0}{1}{0}\pgfsys@invoke{
}\definecolor[named]{pgffillcolor}{rgb}{1,1,0}\definecolor[named]{pgffillcolor}{rgb}{1,1,0}\pgfsys@color@cmyk@fill{0}{0}{1}{0}\pgfsys@invoke{
}{}\pgfsys@moveto{0.0pt}{0.0pt}\pgfsys@moveto{2.15277pt}{0.0pt}\pgfsys@curveto{2.15277pt}{1.18895pt}{1.18895pt}{2.15277pt}{0.0pt}{2.15277pt}\pgfsys@curveto{-1.18895pt}{2.15277pt}{-2.15277pt}{1.18895pt}{-2.15277pt}{0.0pt}\pgfsys@curveto{-2.15277pt}{-1.18895pt}{-1.18895pt}{-2.15277pt}{0.0pt}{-2.15277pt}\pgfsys@curveto{1.18895pt}{-2.15277pt}{2.15277pt}{-1.18895pt}{2.15277pt}{0.0pt}\pgfsys@closepath\pgfsys@moveto{0.0pt}{0.0pt}\pgfsys@fillstroke\pgfsys@invoke{
} \pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope
\pgfsys@invoke{\lxSVG@closescope
}\pgfsys@endscope{}{}{}\hss}\pgfsys@discardpath\pgfsys@invoke{\lxSVG@closescope
}\pgfsys@endscope\hss}}\lxSVG@closescope\endpgfpicture}}\,\,\langle\,\eta_{1}\rangle\,\,\,\,\,:\,\,\,\text{real
part of the charged scalar operator}\,,$ $\displaystyle\leavevmode\hbox
to4.71pt{\vbox to4.71pt{\pgfpicture\makeatletter\hbox{\hskip
2.35277pt\lower-2.35277pt\hbox to0.0pt{\pgfsys@beginscope\pgfsys@invoke{
}\definecolor{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@rgb@stroke{0}{0}{0}\pgfsys@invoke{
}\pgfsys@color@rgb@fill{0}{0}{0}\pgfsys@invoke{
}\pgfsys@setlinewidth{0.4pt}\pgfsys@invoke{ }\nullfont\hbox
to0.0pt{\pgfsys@beginscope\pgfsys@invoke{
}{}{{}}{}{{{}}{}{}{}{}{}{}{}{}}\pgfsys@beginscope\pgfsys@invoke{
}\color[rgb]{0.5,0,0.5}\definecolor[named]{pgfstrokecolor}{rgb}{0.5,0,0.5}\pgfsys@color@rgb@stroke{0.5}{0}{0.5}\pgfsys@invoke{
}\pgfsys@color@rgb@fill{0.5}{0}{0.5}\pgfsys@invoke{
}\definecolor[named]{pgffillcolor}{rgb}{0.5,0,0.5}\definecolor[named]{pgffillcolor}{rgb}{0.5,0,0.5}\pgfsys@color@rgb@fill{0.5}{0}{0.5}\pgfsys@invoke{
}{}\pgfsys@moveto{0.0pt}{0.0pt}\pgfsys@moveto{2.15277pt}{0.0pt}\pgfsys@curveto{2.15277pt}{1.18895pt}{1.18895pt}{2.15277pt}{0.0pt}{2.15277pt}\pgfsys@curveto{-1.18895pt}{2.15277pt}{-2.15277pt}{1.18895pt}{-2.15277pt}{0.0pt}\pgfsys@curveto{-2.15277pt}{-1.18895pt}{-1.18895pt}{-2.15277pt}{0.0pt}{-2.15277pt}\pgfsys@curveto{1.18895pt}{-2.15277pt}{2.15277pt}{-1.18895pt}{2.15277pt}{0.0pt}\pgfsys@closepath\pgfsys@moveto{0.0pt}{0.0pt}\pgfsys@fillstroke\pgfsys@invoke{
} \pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope
\pgfsys@invoke{\lxSVG@closescope
}\pgfsys@endscope{}{}{}\hss}\pgfsys@discardpath\pgfsys@invoke{\lxSVG@closescope
}\pgfsys@endscope\hss}}\lxSVG@closescope\endpgfpicture}}\,\,\langle\,\eta_{2}\rangle\,\,\,\,\,:\,\,\,\text{imaginary
part of the charged scalar operator}\,,$ $\displaystyle\leavevmode\hbox
to4.71pt{\vbox to4.71pt{\pgfpicture\makeatletter\hbox{\hskip
2.35277pt\lower-2.35277pt\hbox to0.0pt{\pgfsys@beginscope\pgfsys@invoke{
}\definecolor{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@rgb@stroke{0}{0}{0}\pgfsys@invoke{
}\pgfsys@color@rgb@fill{0}{0}{0}\pgfsys@invoke{
}\pgfsys@setlinewidth{0.4pt}\pgfsys@invoke{ }\nullfont\hbox
to0.0pt{\pgfsys@beginscope\pgfsys@invoke{
}{}{{}}{}{{{}}{}{}{}{}{}{}{}{}}\pgfsys@beginscope\pgfsys@invoke{
}\color[rgb]{0.666667,0.666667,0}\definecolor[named]{pgfstrokecolor}{rgb}{0.666667,0.666667,0}\pgfsys@color@rgb@stroke{0.666667}{0.666667}{0}\pgfsys@invoke{
}\pgfsys@color@rgb@fill{0.666667}{0.666667}{0}\pgfsys@invoke{
}\definecolor[named]{pgffillcolor}{rgb}{0.666667,0.666667,0}\definecolor[named]{pgffillcolor}{rgb}{0.666667,0.666667,0}\pgfsys@color@rgb@fill{0.666667}{0.666667}{0}\pgfsys@invoke{
}{}\pgfsys@moveto{0.0pt}{0.0pt}\pgfsys@moveto{2.15277pt}{0.0pt}\pgfsys@curveto{2.15277pt}{1.18895pt}{1.18895pt}{2.15277pt}{0.0pt}{2.15277pt}\pgfsys@curveto{-1.18895pt}{2.15277pt}{-2.15277pt}{1.18895pt}{-2.15277pt}{0.0pt}\pgfsys@curveto{-2.15277pt}{-1.18895pt}{-1.18895pt}{-2.15277pt}{0.0pt}{-2.15277pt}\pgfsys@curveto{1.18895pt}{-2.15277pt}{2.15277pt}{-1.18895pt}{2.15277pt}{0.0pt}\pgfsys@closepath\pgfsys@moveto{0.0pt}{0.0pt}\pgfsys@fillstroke\pgfsys@invoke{
} \pgfsys@invoke{\lxSVG@closescope }\pgfsys@endscope
\pgfsys@invoke{\lxSVG@closescope
}\pgfsys@endscope{}{}{}\hss}\pgfsys@discardpath\pgfsys@invoke{\lxSVG@closescope
}\pgfsys@endscope\hss}}\lxSVG@closescope\endpgfpicture}}\,\,\langle\,\delta\phi^{x}\rangle\,:\,\,\,\text{longitudinal
axion operator}\,.$ (29)

In eq. (42) we express these expectation values in terms of the asymptotics of
the fluctuations (18). Notice that in the superfluid phase we are rewriting
the operator dual to the fluctuation of the complex scalar in terms of its
modulus and phase, which we denote by $\delta H$ and $\varphi$ respectively.
Moreover, for the model breaking translations, we switch on an extra
fluctuation, $\delta\phi^{x}$. The meaning of all of these quantities is well
known apart from $\delta\phi^{x}$, corresponding to the longitudinal axion
operator. The colored circles in (29) display the color scheme used in all the
figures in this section.

In the remainder of this section we will present results for the amplitudes
and phases of the expectation values (29) for the different QNMs of the
system. These are all computed in Fourier space (see eq. (20)), at finite
frequency and momentum. Finally, we restrict our analysis of the amplitudes
and phases to the following hydrodynamic modes

  1. 1.

first and second sounds in the superfluid phase (with preserved translations);

  2. 2.

fourth sound and energy diffusion in the superfluid phase with broken
translations;

  3. 3.

first sound and a purely diffusive mode in the normal phase (with preserved
translations);

  4. 4.

thermo-electric diffusion modes (energy and charge diffusion in the decoupling
limit) in the normal phase with broken translations.

Additionally, in the normal phase we shall also consider the gapped scalar
mode that drives the superfluid phase transition at $T=T_{c}$.

Finally, let us explain our normalization of the amplitudes. We may formulate
the linearized equations for the fluctuations in terms of a generalized
eigenvalue problem (see appendix D for more details)

$(\bm{A}(k)\,\omega-\bm{B}(k)\,)\,\bm{x}(k)=0,$ (30)

where $\bm{A}$ and $\bm{B}$ are differential operators and $\bm{x}$ is a
vector consisting of the fluctuations. Solving this generalized eigenvalue
problem numerically, we find a set of eigenfunctions $\bm{x}_{n}$ (consisting
of the fluctuations) for each eigenvalue $\omega_{n}$ (the quasinormal
frequency). We will refer to the set of eigenfunctions as eigenvector from now
on. First, we make all expectation values of the eigenvectors consisting of
the fluctuations (29) dimensionless by dividing by the appropriate power of
the temperature (depending on the conformal dimension of the corresponding
operator). Then, we normalize the sum of all absolute values of the
dimensionless expectations values of the eigenvectors to one. This means that
an amplitude of $0.5$ implies for example that $50\%$ of the support of such a
mode is controlled by this type of fluctuation. With this normalization, the
amplitudes are invariant under the scaling symmetry inherent to the linear
analysis which acts on all fluctuations as $\langle\delta
f^{i}\rangle\rightarrow\lambda\langle\delta f^{i}\rangle$.

Figure 5: The amplitudes of the various field theory operators (29) in the
normal phase ($T>T_{c}$) at fixed $k/T=0.01$ for the diffusive mode (left),
first sound (center), and the gapped scalar modes (right). The color scheme is
that illustrated in (29).

In order to make the reader familiar with the language, we start our analysis
in the normal phase which corresponds to considering a standard relativistic
charged fluid in the field theory dual. There, the two hydrodynamic modes are
first sound and the purely diffusive mode. Nevertheless, for completeness we
will consider also the gapped scalar modes which are the responsible for the
superfluid instability at $T\leq T_{c}$ hydroSC . Notice that in the normal
phase the background value of the scalar is zero and therefore the
fluctuations of the phase are not a well defined object. Thus, in the normal
phase, we will consider the expectation value of the fluctuations of the
imaginary part of the scalar operator dual to the complex bulk field $\psi$.

The results for the amplitudes of the normal fluid are shown in fig. 5. As
expected, we find that the gapped scalar mode is solely carried by the
fluctuations of the scalar charge operator. The purely diffusive mode (left
panel of fig. 5) is carried by the fluctuations of the charge density
$\langle\delta J^{t}\rangle$, as one can derive analytically from
hydrodynamics (see appendix (62)). The hydrodynamic treatment shows that in
first order hydrodynamics this mode consists solely of $\langle\delta
J^{t}\rangle$ (and of the by a Ward identity related $\langle\delta
J^{t}\rangle$). Moreover, the second order correction includes $\delta T^{tt}$
with a relative amplitude of $\langle\delta T^{tt}\rangle/(\langle\delta
J^{t}\rangle\,T)=k^{2}/T\,\rho/(\mu\rho+s\,T)$. We shall point out that our
result that the diffusive mode is solely carried by $\langle\delta
J^{t}\rangle$ is not in contradiction with the findings of Hartnoll:2014lpa ;
Davison:2015taa where it was observed that the diffusive mode is carried by
the ‘incoherent charge’ $\delta Q^{\text{diff}}\equiv\delta
J^{t}-\rho/(\mu\rho+sT)\,\delta T^{tt}$ 666We thank Blaise Goutéraux and
Richard Davison for discussions about this point.. As we detail in appendix C,
our computation describes the situation where only the diffusive mode is
excited in the system. The authors of Hartnoll:2014lpa ; Davison:2015taa
consider a generic perturbation in which all modes are excited instead. We
show in that appendix that one can relate both situations via the hydrodynamic
framework and match them to the holographic computation.

Finally, first sound is dominated by the fluctuations of energy density and is
almost entirely supported on the fluctuations of energy density, pressure and
longitudinal momentum. The ratio between pressure and energy fluctuations is
approximately $0.5$ as we expect from the Ward identities that constrain the
expectation values given in eq. (46). Note that the pressure fluctuations
contribute approximately equally via $\langle\delta T^{xx}\rangle$ and
$\langle\delta T^{yy}\rangle$ to first sound as is obvious from the
tracelessness of the energy-momentum tensor and the Ward identities (46). We
furthermore observe that the amplitudes of the fluctuation of longitudinal
momentum are approximately 0.7 time the energy fluctuations (or 1.4 times the
pressure fluctuations) which may be seen from eq. (46) or equivalently from
eq. (58) using the conformal speed. From the inset of fig. 3 we know that in
the range of temperatures shown there, the diffusion constant takes values in
the range $D_{q}\,T\approx\\{0.04,0.1\\}$. Using the Ward identity (47), we
find that the fluctuations of the longitudinal current in the purely diffusive
mode should be approximately $\\{4\cdot 10^{-4},10^{-3}\\}$ of the charge
fluctuation as evident in the left panel in fig. 5 and in eq. (53). The simple
relation eq. (53) is Fick’s law for diffusion. This discussion shows that we
can use the amplitudes and the Ward identities to compute the diffusion
constant in the hydrodynamic dispersion relation and vice versa. In the case
of first sound a similar analysis holds albeit it is a bit more complicated.
We refer the interested reader to appendix B. In the normal phase, the scalar
sector decouples and the fluctuations of the scalar operator do not enter into
either of the sound modes in the normal phase.

Figure 6: The amplitudes of the various operators (29) for the thermo-electric
diffusive modes in the normal phase with momentum dissipation. We fixed
$k/T=0.01$ and $T/\mu=0.025$. The color scheme is that illustrated in (29).

As already mentioned, when the axion fields are switched on and momentum is
not a conserved quantity anymore, first sound is destroyed at large wavelength
and the only two hydrodynamic modes in the normal phase are the thermo-
electric diffusive modes Kovtunbound . We show the amplitudes for these two
modes in fig. 6. At $\alpha=0$, the left panel corresponds to the diffusive
mode shown in fig. 5. Indeed, for $\alpha/T\ll 1$, such a mode is totally
dominated by charge fluctuations. Interestingly, around $\alpha/T\approx 8$,
we observe a smooth crossover to a regime where the original diffusive mode is
now dominated by energy fluctuations. Moreover, the ratio between pressure and
energy fluctuations is approximately $0.5$ for all $\alpha/T$, as we already
observed for the first sound mode which is a consequence of the tracelessness
of the energy-momentum tensor (since $\langle\delta
T^{tt}\rangle=\langle\delta T^{xx}\rangle+\langle\delta T^{yy}\rangle$).
Notice that in that large $\alpha/T\gg 1$ regime the contribution from charge
fluctuations drops to zero indicating a drastic change in the nature of the
diffusive mode.

As shown on the right panel of fig. 6, the situation is different for the
other diffusive mode present in the normal phase of the model with momentum
dissipation. For small values of momentum dissipation, $\alpha/T\lesssim 2$,
the momentum chosen for the extraction of the amplitudes, $k/T=0.01$, is too
large and the mode is not diffusive anymore but already sound-like, i.e.
propagating. There, we find the same features as for first sound in fig. 5.
Around $\alpha/T\approx 2$, momentum is not a well-defined quantity anymore,
and its contribution to the mode drops rapidly to zero. At very large value of
momentum dissipation, the mode in the right panel of fig. 6 is totally
decoupled from everything else and it is dominated by energy and pressure.
Note that for $\alpha/T\gtrsim 200$, we enter the regime where the scalar
modes are unstable since the critical temperature (which is a function of the
momentum dissipation strength) is now larger the temperature we chose.
Similarly to the diffusive mode, we observe that the ratio between energy and
pressure fluctuations is approximately 0.5 for all $\alpha/T$.

Figure 7: Top: The amplitudes of the various operators for the two
thermoelectric diffusion modes at $T/\mu=5$. Bottom: The corresponding
diffusion constants. The blue points are numerical data obtained from the
amplitude data in the upper panel using the Fick’s law (53) as explained in
appendix B. The red lines are the analytical formulas (31) valid in the
incoherent limit, which is reached around $\alpha/T\sim 10$. At that point,
the support of the two modes decouple into energy density and charge density
fluctuations. See also the figures in appendix C.2.

Before moving on to the broken phase, let us investigate the normal phase with
broken translations further; in particular, we are interested in the two
thermoelectric diffusion modes. At intermediate values of $T/\mu$ and
$\alpha/T$, energy diffusion and charge diffusion are coupled and therefore
their support is mixed between the various operators. This is no longer true
in the so-called incoherent limit Hartnoll:2014lpa , at which the momentum
dissipation rate becomes the dominant scale in the system, i.e.
$\alpha/T,\;\alpha/\mu\gg 1$. In such a regime, we do expect the two diffusive
modes to decouple again Davison:2015bea and the corresponding diffusion
constants acquire the simple values

$D_{e}\,=\,\frac{\kappa}{c_{v}}\,,\qquad
D_{q}\,=\,\frac{\sigma}{\chi_{\rho\rho}}\,,$ (31)

where $\kappa$ is the thermal conductivity, $c_{v}$ the specific heat,
$\sigma$ the electric conductivity and $\chi_{\rho\rho}$ the charge
susceptibility. For the linear axion model considered in this work, the
formulas for the diffusion constants above are known analytically and can be
found in the literature (see for example Baggioli:2017ojd ). In fig. 7, we
repeat the analysis of fig. 6 but this time going into the aforementioned
incoherent limit. As expected, in this case, the contributions to the support
of the two diffusive modes from charge and energy fluctuations crosses around
$\alpha/T\sim 10$, which corresponds exactly to the location at which the
diffusion constants approach the values mentioned in eq. (31), i.e. the edge
of the incoherent limit. Going towards $\alpha/T\rightarrow\infty$, the
support of the two modes becomes solely dominated by either charge or charge
and energy fluctuations, as expected from their decoupling. This result is
explained in more detail in appendix C.2 where we show the relative
coefficients in fig. 14 explicitly. Note that we extracted the diffusion
constants from the amplitude data at one fixed $k/T$ using eq. (53). Not only
that, but the amplitudes of these two diffusive modes can be derived
analytically at any value of the parameters from hydrodynamics (see appendix
(62)), matching perfectly the numerical results presented here.

Figure 8: The amplitudes of the various operators (29) for first (left panel)
and second (right panel) sound normalized to $T$ at fixed $k/T=0.01$. The
inset is a non-logarithmic presentation of the same curves to show the ratio
of approximately $0.5$ between pressure and energy fluctuations. The color
scheme is that illustrated in (29).

We now move on to the analysis of the superfluid phase of our system. We study
the amplitudes of the operators (29) for the two hydrodynamic modes in the
longitudinal sector of this phase: first and second sound. The results are
plotted in fig. 8. Notice that in this condensed phase we represent the
fluctuations of the scalar operator in terms of its amplitude (the Higgs mode
$\delta H$) and the time derivative of its phase which carries information
about the norm of the Goldstone boson since in our background
$|A_{\mu}-i\,\partial_{\mu}\varphi|\sim\mu-i\,\partial_{t}\varphi+\mathcal{O}(\varphi^{2})$.
As one can see on the left panel of fig. 8, the behavior of the amplitudes of
the first sound mode is rather featureless and dominated by energy, pressure,
and longitudinal momentum fluctuations across the phase transition and down to
the lowest temperatures. The fluctuations of the Goldstone mode are highly
suppressed, and those of charge and current decrease as T is lowered. We show
the results for second sound on the right panel of 8. This mode is dominated
by the fluctuations of the superfluid condensate (the Higgs mode). In view of
these results one can rule out a role reversal between first and second sound
in which the nature of the fluctuations supporting those modes is interchanged
as the temperature is lowered across the superfluid phase as in Alford:2013koa
. Finally, notice that the fluctuations of the scalar operator are not shown
above the critical temperature in the left panel of fig. 8. They decouple from
the other operators and therefore do not contribute to the support of the
first sound mode above $T_{c}$.

Figure 9: The amplitudes of the various operators (29) for the diffusive mode
(left) and $4^{th}$ sound (right) in the superfluid phase with momentum
dissipation. We fixed $\alpha/T=6$ and $k/T=0.01$. The color scheme is that
illustrated in (29).

We can now perform the same analysis in the superfluid phase in presence of
momentum dissipation. In this case, first sound is replaced by a diffusive
mode and second sound is replaced by fourth. In fig. 9 we show the results for
the amplitudes of the fluctuations for the diffusive mode (left) and fourth
sound (right). Except at temperatures near the critical temperature, where the
Higgs mode takes over, the diffusive mode is dominated by energy and pressure
fluctuations. Fourth sound is always dominated by energy fluctuations.
Interestingly, at low temperature the ratios between energy, longitudinal
momentum, and pressure fluctuations asymptote to the same ratios observed for
first sound since the speed of fourth sound tends to the speed of first sound
as $T\to 0$ there (while the attenuation constant tends to zero). Finally, we
shall point out that we also observe a purely imaginary, gapped non-
hydrodynamic mode which is solely supported by the scalar fluctuations. This
higher quasi-normal mode might be the Higgs mode discussed in Behrle_2018 ;
Bhaseen:2012gg ; Li:2013fhw ; Pekker_2015 .

####  4.2.1 About the relative phases

As explained in detail in section 4.1, the value of the amplitudes of the
various fluctuations (29) of the different (hydrodynamic) modes is not the
only information we can extract from our analysis. We can also compute the
phase differences between the various fluctuations present in a specific
(hydrodynamic) mode. We used the phase of the charge density fluctuation
$\Theta_{\delta J^{t}}$ as reference phase and computed the phases of the
other fluctuations with respect to $\Theta_{\delta J^{t}}$. Note that the
phases are defined modulo $2\pi$. The phase differences of the various
fluctuations are listed in table 1. The phases differences remain almost
constant throughout the phase diagram. However, it is important to notice that
the “ideal” phase differences as listed in 1 are only reached in the $k\to 0$
limit. As we show in eq. (54) and eq. (59), the relative phases between
certain fluctuations are intimately related to damping, i.e. for vanishing
phase differences we reach ideal hydrodynamics. In other words, the relative
phases give us insights about the damping in the system and we may compute the
diffusion and attenuation constants, respectively, from them.

without momentum dissipation |  $\delta T^{tt}$ |  $\delta T^{tx}$ |  $\delta T^{xx}$ |  $\delta T^{yy}$ |  $\delta J^{x}$ |  $\delta\mathrm{H}$ |  $\partial_{t}\varphi$ |  $\delta\phi^{x}$
---|---|---|---|---|---|---|---|---
diffusive mode ($T>T_{c}$) | $0$ | $0$ | $\pi$ | $\pi$ | $\pi/2$ |  |  |
$1^{st}$ sound with $\mathrm{Re}(\omega)>0$ ($T>T_{c}$) | $\pi$ | $\pi$ | $0$ | $0$ | $\pi$ |  |  |
$1^{st}$ sound with $\mathrm{Re}(\omega)>0$ ($T<T_{c}$) | $\pi$ | $\pi$ | $\pi$ | $\pi$ | $0$ | $0$ | 0 |
$2^{nd}$ sound with $\mathrm{Re}(\omega)>0$ ($T<T_{c}$) | $\pi/2$ | $\pi/2$ | $\pi/2$ | $\pi/2$ | $0$ | $0$ | 0 |
with momentum dissipation |  |  |  |  |  |  |  |
$4^{nd}$ sound with $\mathrm{Re}(\omega)>0$ ($T<T_{c}$) | $\pi$ | $\pi$ | $\pi$ | $\pi$ | $0$ | $0$ | 0 | $\pi/2$
diffusion ($T<T_{c}$) | $\pi$ | $\pi/2$ | $\pi$ | $\pi$ | $-\pi/2$ | $\pi$ | $0$ | $\pi/2$

Table 1: The phase differences modulo 2$\pi$ with respect to $\delta J^{t}$
for the various fluctuations for the normal fluid ($T>T_{c}$) and superfluid
$(T<T_{c})$ phases for temperatures sufficiently far from the phase transition
and in the ideal limit $k\to 0$. The gray cells indicate that such a
fluctuation is not present in the corresponding mode. In the case of the sound
modes, the phases of the “partner” mode may be obtained by taking into account
that the phases add up to 0 (for parity even operators) and $\pi$ (for parity
odd operators), respectively. Concretely, the computation of the phases for
the momentum dissipation case has been performed for $\alpha/T=6$.

####  4.2.2 Momentum dependence of the amplitudes

So far, all the fluctuations amplitudes have been computed for a single value
of momentum, $k/T=0.01$. In this section, we briefly discuss the momentum
dependence of the amplitudes presented in section 4.2. The results are shown
in fig. 10.
For first sound, both the amplitudes and the relative phases are constant in a
large range of momenta (we considered them until $k/T\approx 12$). In the case
of second sound, we observe that the biggest amplitude, the Higgs
fluctuations, remains dominant until $k/T\approx 10$ and the relative phases
approximately constant. Away from the small momentum limit, we notice that the
contribution from the fluctuations of the superfluid condensate becomes
smaller and the mode is eventually dominated by the energy fluctuations which
are the biggest contribution at large $k$. At the same time, at $k/T\sim 1$,
we observe a crossover between the fluctuations of the charge density (which
is the first subleading contribution at small $k$, even thought it accounts
for less than $10\%$ of the total weight) and that of energy and longitudinal
momentum which become more relevant at large wave-vector.
All in all, we can state that, in the hydrodynamic limit $k/T\ll 1$, our
results depend only minimally on the specific value of the momentum $k$.

Figure 10: The momentum dependence of the normalized amplitudes (29) at
$T/\mu=0.016\,(T/T_{c}=0.804)$. Top panels: First sound. Central panels:
Second sound; The inset shows a zoom on small momenta (and we omit plotting
the amplitude of the Higgs fluctuation in orange). Bottom panels: Second sound
beyond the linear part of the dispersion relation ($k/T\gg 1$).

##  5 Alternative approach via complex momentum modes

###  5.1 Theory

Quasinormal modes can be understood as the response of the system to a
perturbation that happened in the past. Thus, quasinormal modes describe how
the system relaxes back to equilibrium at late times.

A different way of studying the response is by perturbing the system with a
periodic source of frequency $\omega$ and looking for the corresponding wave
with complex wave number $q=k+i\kappa$. If we were to excite only one
particular mode then the response would be

$\Phi^{I}(z,t,x)=\tilde{A}^{I}_{n}(\omega,z)e^{-i\omega
t+ik_{n}x}e^{-\kappa_{n}|x|}\,.$ (32)

The complex momentum is now a function of the real frequency
$k_{n}=k_{n}(\omega)$ and $\kappa_{n}=\kappa_{n}(\omega)$. This was
investigated in Amado:2008ji where it was shown that one can isolate the modes
that decay exponentially by choosing the integration in the complex momentum
plane appropriately. Only a few other works considered this choice within the
holographic framework Amado:2007pv ; Forcella:2014dwa ; Blake:2014lva .

We can now proceed with the same attitude as in the previous section and
investigate the amplitudes and the relative phases. Whether these are
different from the ones of the quasinormal modes is an open question and it is
our main motivation for this last analysis. A naive expectation would be that
the phases and amplitude ratios are the same as long as the mode in question
has not gone through some crossing point with another mode in the complex
plane. Generally, such a crossing point has complex frequency and complex
momentum so its neither visible from the quasinormal mode (real momentum) nor
from the complex momentum mode (real frequency).

###  5.2 Results

Before discussing the amplitudes and the phases of the various operators, let
us briefly comment on the nature of the low-energy modes for complex momentum
and real frequency. Let us start with the equation for a generic attenuated
sound wave

$\omega^{2}\,=\,v^{2}\,k^{2}\,-\,i\,\omega\,2\,\Gamma\,k^{2}\,+\,\dots\,,$
(33)

where as always $v,\Gamma$ are the propagation speed and the attenuation
constant, respectively. Let us now consider eq. (33) in terms of a real-valued
frequency and a complex valued momentum. At low frequency, it is easy to
verify that

$\displaystyle\mathrm{Re}(k)\,=\,\pm\,\frac{\omega}{v}\,+\,\mathcal{O}\left(\omega^{3}\right)\,,$
$\displaystyle\mathrm{Im}(k)\,=\,\pm\,\frac{\Gamma}{v^{3}}\,\omega^{2}\,+\,\mathcal{O}\left(\omega^{4}\right)\,.$
(34)

The same analysis can be repeated with a damped diffusive mode (such as the
diffusive mode in the superfluid phase)

$\omega\,=\,-\,i\,\gamma\,-\,i\,D\,k^{2}$ (35)

yielding

$\displaystyle\mathrm{Re}(k)\,=\,\pm\,\frac{\omega}{2\,\sqrt{\gamma\,D}}\,+\,\mathcal{O}\left(\omega^{3}\right)\,,$
$\displaystyle\mathrm{Im}(k)\,=\,\pm\,\sqrt{\frac{\gamma}{D}}\,+\,\mathcal{O}\left(\omega^{2}\right)\,.$
(36)

The equations (34) and (36) represent the low-momentum behavior of the low-
energy modes in this new inverted picture. As expected, our numerical results
shown in fig. 11 are in perfect agreement with this picture.

Figure 11: The three lowest modes in the superfluid phase. The temperature is
fixed at $T/\mu=0.019$ ($T/T_{c}=0.92$). Second sound (green), first sound
(red), the former damped diffusive mode (black). The dashed lines in the top
panel show the expectation from eq. (34).

We can now analyze the amplitudes of the fluctuations in (29) at constant
frequency. We show them in fig. 12 for a fixed small value of (real) frequency
$\omega/T=0.01$. The results for first sound are exactly identical to those
obtained at complex frequency and real momentum (inset in the left panel of
fig. 8). In the second sound case shown on the right panel of fig. 12 (to be
compared with the right panel of fig. 8), the dominant contributions
corresponding to the condensate and charge density fluctuations are identical
while the subleading contributions are slightly different between the real and
imaginary frequency modes. Still the qualitative behavior in the two
frameworks is very similar.

Figure 12: The amplitudes of the fluctuations (29) for first (left) and second
sound (right) at fixed $\omega/T=0.01$ in the superfluid phase ($T<T_{c}$)
without axion fields.

##  6 Conclusions

In this work we have revisited the holographic S-wave superfluid model
Hartnoll:2008kx focusing on the hydrodynamic description of the low-energy
modes and taking into account the backreaction on the gravity side. Despite
the long history of this model, a complete analysis of this sort was missing.
Through this approach we have explicitly verified the match between
relativistic superfluid hydrodynamics and the holographic model, both in terms
of the dispersion relations and, most importantly, in terms of all the first
order transport coefficients. Our study has revealed that all the hydrodynamic
transport coefficients are continuous at the phase transition. Moreover,
differently from what Amado:2013xya found in the probe limit, even the second
sound attenuation constant is smoothly connected to the gapped scalar modes in
the normal phase. We also analyzed the low energy modes for real frequency and
complex momentum, finding agreement with the hydrodynamic description.
Finally, we extended the analysis to a setup where momentum was not conserved
and succesfully matched the peculiar fourth sound mode with its hydrodynamic
counterpart. Interestingly, we find a non-zero normal density in the setup
with broken translations, similarly to the results of Gouteraux:2019kuy ;
Gouteraux:2020asq for hyperscaling-Lifschitz type geometries.

After demonstrating the validity of the hydrodynamic description of the
holographic superfluid, we carried out a study of the support of the various
hydrodynamic modes in terms of the dual field theory operators. In a nutshell,
we have quantified “how much” a single boundary field theory operator
contributes to a specific hydrodynamic mode. The analysis has been performed
both in the superfluid and normal phases and also in presence of momentum
dissipation. In contrast to the field theory results of Alford:2013koa , we do
not observe a role-reversal phenomenon between first and second sound in terms
of the support of those modes. From the amplitudes and phases, we can
reconstruct the dispersion relations of the hydrodynamic modes. This allows us
to compute the speed of sound and the attenuation/diffusion constants from our
data of the amplitudes and relative phases at a fixed value of the momentum.
We notice that the absolute values of the amplitudes give us insights about
the speed of the sound modes, and the relative phases between the fluctuations
are related to damping. Furthermore, the amplitudes and relative phases give
us an (computational) easy recipe to compute the hydrodynamic dispersion
relations since it is sufficient to determine them at only one point of the
dispersion relation (within the radius of convergence). Beyond hydrodynamics,
we observe a purely imaginary gapped mode in the superfluid phase which is
solely supported by the scalar fluctuations; this could be the Higgs mode
discussed in Behrle_2018 ; Li:2013fhw ; Pekker_2015 .

In the future, it would be interesting to enlarge the analysis in this work to
scenarios like those in the following list

  * •

Consider more complex symmetry breaking patterns such as the pseudo-
spontaneous breaking of the $U(1)$ symmetry (as initiated in Donos:2021pkk ).
Even though the explicit breaking of a $U(1)$ vector symmetry is definitely
not a physical option, this could serve as a toy-model to improve our
understanding of the pseudo-spontaneous breaking of translations in the
homogeneous models Baggioli:2014roa ; Andrade:2017cnc ; Alberte:2017cch ;
Amoretti:2016bxs ; Amoretti:2017frz ; Amoretti:2019kuf ; Ammon:2019apj ;
Ammon:2020xyv and the possible relation between phase relaxation and pseudo-
spontaneous breaking Baggioli:2019abx ; Ammon:2019wci ; Amoretti:2018tzw ;
Amoretti:2019cef ; Baggioli:2020haa ; Baggioli:2020edn ; Baggioli:2020nay .
Finally, if these two breaking mechanisms are analogous or even equivalent,
one could try to tackle the phenomenology of superfluid current relaxation
induced by vortices Davison:2016hno within this effective construction.

  * •

Since we investigate the match of hydrodynamics and our holographic model, it
would be interesting to determine the radius of convergence of linearized
superfluid hydrodynamics throughout the phase diagram following the recent
studies Grozdanov:2019kge ; Grozdanov2019 ; Withers:2018srf ; Abbasi:2020ykq ;
Abbasi:2020xli ; Jansen:2020hfd ; Baggioli:2020loj ; Arean:2020eus ;
Wu:2021mkk ; Baggioli:2021ujk ; Grozdanov:2021gzh ; Grozdanov:2020koi ;
Jeong:2021zsv in other holographic models.

  * •

Compute the second order contributions to equilibrium transport, as initiated
in Grieninger:2021rxd ; Kovtun:2018dvd , in superfluids and systems with
second order phase transitions in general.

  * •

Work out the hydrodynamics of holographic superfluids in presence of a
background superfluid velocity Arean:2010zw ; Arean:2010wu ; Arean:2011gz .
The excitations in the probe limit have been already presented in
Amado:2013aea and the dynamical Landau instability in Lan:2020kwn . It would
be interesting to extend them to the fully backreacted scenario and compare
them to the hydrodynamic description.

  * •

Holographic superfluid models in presence of disorder have been constructed in
the past years Arean:2013mta ; Arean:2015sqa . It would be interesting to
study the transport of those models and compare it with our results in the
simpler axion model.

We leave these questions for the near future.

## Acknowledgements

We thank Martin Ammon, Andrea Amoretti, Seán Gray, Hyun-Sik Jeong, Wei-Jia Li,
Yan Liu, and Ya-Wen Sun for useful discussions on related topics. We thank
Blaise Goutéraux and Richard Davison for several useful comments on the first
version of this manuscript and for useful discussions about the interpretation
of our results. D.A. and S.G. are supported by the ‘Atracción de Talento’
program (2017-T1/TIC-5258, Comunidad de Madrid) and through the grants
SEV-2016-0597 and PGC2018-095976-B-C21. M.B. acknowledges the support of the
Shanghai Municipal Science and Technology Major Project (Grant
No.2019SHZDZX01).

##  Appendix A The pocket dictionary

In this appendix, following the standard holographic renormalization
procedure, we compute the one- and two-point functions relevant for the
analysis in the main text. We start by computing the variation of the
renormalized on-shell action. To the action (12) we add the standard
counterterms

$S_{\text{ct}}={1\over
2\kappa_{4}^{2}}\int\mathrm{d}^{3}x\sqrt{-\gamma}\,\left(2K-4-|\psi|^{2}-R^{(\gamma)}\right)\,.$
(37)

A straightforward computation leads to the following expression for the
variation of the renormalized on-shell action (see also Herzog:2009md )

$\displaystyle\\!\delta S^{\rm o-s}_{\rm ren}$
$\displaystyle\\!=\\!-\\!\int\\!\mathrm{d}^{3}x\sqrt{-\gamma}\bigg{\\{}\\!n_{m}F^{mn}\delta
A_{n}+\\!2n_{m}\left[\delta\psi
g^{mn}\left(\partial_{n}\\!+\\!iq\,A_{n}\right)\bar{\psi}\\!+\\!\delta\bar{\psi}g^{mn}\left(\partial_{n}\\!-\\!iq\,A_{n}\right)\psi\right]$
$\displaystyle-\delta\gamma^{mn}\\!\left(\\!K_{mn}\\!-(K+2+\\!{|\psi|^{2}\over
2}+\\!{R^{(\gamma)}\over
2})\,\gamma_{mn}\\!-\\!R^{(\gamma)}_{mn}\\!\right)+2\psi\delta\bar{\psi}+2\bar{\psi}\delta\psi\bigg{\\}}\bigg{|}_{z=\epsilon}\\!\\!\\!\\!\\!\\!\\!,$
(38)

where $\gamma_{mn}$ denotes the induced metric at a $z=\epsilon$ hypersurface,
while $K_{mn}$ is the extrinsic curvature at that hypersurface, $K$ is its
trace, and $R^{(\gamma)}_{mn}$ the Ricci tensor.777We write the induced metric
and extrinsic curvature in the four-indices notation, see e.g. McNeesnotes .

In order to read the one-point functions from (38) we need the boundary
expansions of the bulk fields. These take the generic form

$\displaystyle
g_{\mu\nu}=g_{\mu\nu}^{(0)}(1+o(z))+g_{\mu\nu}^{(3)}\,z^{3}(1+o(z))\,,$ (39)
$\displaystyle
A_{\mu}=A_{\mu}^{(0)}+A_{\mu}^{(1)}z+o(z^{2})\,,\qquad\psi=\psi^{(1)}z+\psi^{(2)}z^{2}+o(z^{3})\,.$
(40)

For the background solution (indicated with a hat symbol), we have

$\hat{\psi}^{(1)}=0\,,\qquad\hat{\psi}^{(2)}=\langle\mathcal{O}^{\psi}_{\text{cond}}\rangle/2\,,\qquad\hat{A}_{t}^{(0)}=\mu\,,\qquad\hat{A}_{t}^{(1)}=-\rho,\qquad\hat{g}_{tt}^{(3)}=-u_{3}\,.$
(41)

Hence, we are taking the leading contribution of the scalar towards the
boundary to be real. As usual this will correspond to picking a real value for
the symmetry-breaking condensate.

Equations (38-41) result in the following expressions for the expectation
values of the dual field theory operators up to first order in fluctuations

$\displaystyle\langle\Psi\rangle={1\over\sqrt{-\gamma_{(b)}}}{\delta S^{\rm
o-s}_{\rm ren}\over\delta\bar{\psi}_{({\rm
b})}}=\langle\mathcal{O}^{\psi}_{\text{cond}}\rangle+2\psi^{(2)}(x,t)\,,$
(42a)
$\displaystyle\langle\bar{\Psi}\rangle={1\over\sqrt{-\gamma_{(b)}}}{\delta
S^{\rm o-s}_{\rm ren}\over\delta\psi_{({\rm
b})}}=\langle\mathcal{O}^{\psi}_{\text{cond}}\rangle^{*}+2\psi^{(2)*}(x,t)\,,$
(42b) $\displaystyle\langle J^{t}\rangle={1\over\sqrt{-\gamma_{(b)}}}{\delta
S^{\rm o-s}_{\rm ren}\over\delta A^{t}_{({\rm b})}}=\rho-a_{t}^{(1)}(x,t)\,,$
(42c) $\displaystyle\langle J^{i}\rangle={1\over\sqrt{-\gamma_{(b)}}}{\delta
S^{\rm o-s}_{\rm ren}\over\delta A^{i}_{({\rm
b})}}=a_{i}^{(1)}(x,t)\,,\quad(i=x,y)\,,$ (42d) $\displaystyle\langle
T_{tt}\rangle={2\over\sqrt{-\gamma_{(b)}}}{\delta S^{\rm o-s}_{\rm
ren}\over\delta\gamma^{tt}_{({\rm b})}}=2u_{3}-2h_{tt}^{(3)}(x,t)\,,$ (42e)
$\displaystyle\langle T_{xx}\rangle={2\over\sqrt{-\gamma_{(b)}}}{\delta S^{\rm
o-s}_{\rm ren}\over\delta\gamma^{xx}_{({\rm
b})}}=u_{3}-h_{tt}^{(3)}(x,t)+3h_{yy}^{(3)}(x,t)\,,$ (42f)
$\displaystyle\langle T_{yy}\rangle={2\over\sqrt{-\gamma_{(b)}}}{\delta S^{\rm
o-s}_{\rm ren}\over\delta\gamma^{yy}_{({\rm
b})}}=u_{3}-h_{tt}^{(3)}(x,t)+3h_{xx}^{(3)}(x,t)\,,$ (42g)
$\displaystyle\langle T_{tx}\rangle=-3h_{tx}^{(3)}(x,t)\,,\qquad\langle
T_{ty}\rangle=-3h_{ty}^{(3)}(x,t)\,,\qquad\langle
T_{xy}\rangle=-3h_{xy}^{(3)}(x,t)\,,$ (42h)

where the energy density is given by $\varepsilon=2\,u_{3}$. The subindex (b)
stands for the finite value of the different fields towards the boundary,
hence once stripped of the powers of $z$ (e.g.
$\delta\gamma_{mn(b)}=\delta\gamma_{mn}/z^{2}$).

It is worth pointing out that the boundary expansions (40) fulfill the
necessary relations for the relevant Ward Identities to be satisfied. First,
they obey the constraint

$h_{xx}^{(3)}(x,t)+h_{yy}^{(3)}(x,t)=0\,,$ (43)

guaranteeing that $\langle T^{\mu}_{\mu}\rangle=0$, as required by conformal
invariance. As for the Ward Identity $\langle\nabla^{\mu}T_{\mu\nu}\rangle=0$,
that results from translational invariance, it boils down to the following
constraints:

$\displaystyle 2\omega h_{tt}^{(3)}(\omega,k)+3kh_{tx}^{(3)}(\omega,k)=0\,,$
(44a) $\displaystyle-3\omega
h_{tx}^{(3)}(\omega,k)+k(-h_{tt}^{(3)}(\omega,k)+3h_{yy}^{(3)}(\omega,k))=0\,,$
(44b)

where we have assumed an harmonic spacetime dependence of the fluctuations as
in (20). Finally, the Ward Identity $\langle\nabla_{\mu}J^{\mu}\rangle=0$,
stemming from the global $U(1)$ symmetry, is satisfied once

$k\,a_{x}^{(1)}+\omega\,a_{t}^{(1)}=0\,,$ (45)

holds. By solving the equations (13) towards the boundary we have checked that
indeed the constraints (44, 45) are satisfied.

From the Ward identities and the tracelessness of the energy-momentum tensor,
we find the following relations between the expectation values of the
fluctuations

$\displaystyle\langle\delta
T^{tt}\rangle=\frac{1}{1-\omega^{2}/k^{2}}\,\langle\delta
T^{yy}\rangle=-\frac{k}{\omega}\,\langle\delta
T^{tx}\rangle=\frac{k^{2}}{\omega^{2}}\,\langle\delta T^{xx}\rangle,$ (46)
$\displaystyle\langle\delta J^{t}\rangle=\frac{k}{\omega}\,\langle\delta
J^{x}\rangle,$ (47)

where the $\delta$ in the expectation values indicates that we are referring
to the fluctuations.
Regarding the axion fields introducing momentum dissipation, the boundary
expansion for the axion fluctuations reads we have (see Andrade:2013gsa for
the derivation)

$\phi^{I}\,=\,\phi_{(0)}^{I}\,+\,\ldots\,+\,\phi_{(3)}^{I}\,z^{3}$ (48)

Hence, in addition to the expectation values (42), we have the expectation
value of the axion fluctuation Andrade:2013gsa

$\langle\,\phi^{I}\rangle\,=\,6\,\phi_{(3)}^{I}\,.$ (49)

for the operator dual to $\phi^{I}$. Since the axions break translational
invariance, the Ward identities discussed above do no longer hold at the level
of the fluctuations. For a discussion see Kim:2016hzi .

Finally, by computing the one-point functions at sources on from (38), we can
determine the form of the correlators relevant for the Kubo-formulas employed
in the main text:

$\displaystyle\langle\Psi\Psi\rangle={\psi^{(2)}\over\psi^{(1)}}+i(q\mu-\omega)\,,\qquad\langle
T_{xy}T_{xy}\rangle=u_{3}-3{h_{xy}^{(3)}\over
h_{xy}^{(0)}}+i\omega\left({k^{2}\over 2}-\omega^{2}\right),$ (50a)
$\displaystyle\langle J^{x}J^{x}\rangle={a_{x}^{(1)}\over
a_{x}^{(0)}}+i\omega\,,\qquad\langle J^{x}J^{t}\rangle=-{a_{x}^{(1)}\over
a_{t}^{(0)}}+ik\,.$ (50b)

##  Appendix B Dispersion Relations from Amplitudes

From the amplitudes computed at a certain value of the momentum, we might
extract the diffusion constants and the hydrodynamic dispersion relations. The
dispersion relation of the purely diffusive mode in the normal phase is for
example given by eq. (5)

$\omega=-i\,D_{q}\,k^{2}.$ (51)

Using the Ward identity eq. (47), we find

$\langle\delta J^{x}\rangle=-i\,D_{q}\,k\,\langle\delta J^{t}\rangle$ (52)

or in other words, the diffusion constant is given by

$D_{q}=\frac{i}{k}\frac{\langle\delta J^{x}\rangle}{\langle\delta
J^{t}\rangle}$ (53)

which is nothing else than the very well-known Fick’s law. The imaginary unit
in eq. (53) simply implies that the amplitudes have a phase difference of
$\pi/2$

$D_{q}=\frac{1}{k}\frac{|\langle\delta J^{x}\rangle|}{|\langle\delta
J^{t}\rangle|}\,e^{i(\Theta_{\delta J^{x}}-\Theta_{\delta
J^{t}}+\pi/2)}\overset{D_{q}\in\mathbb{R}}{=}\frac{1}{|k|}\frac{|\langle\delta
J^{x}\rangle|}{|\langle\delta J^{t}\rangle|}\sin(\Theta_{\delta
J^{t}}-\Theta_{\delta J^{x}}),$ (54)

where we required the diffusion constant to be real in the last step.
Requiring the diffusion constant to be real and positive and thus the cosine
to vanish fixes the relative amplitudes to $\pi/2$ mod $2\pi$. From this
equation it is obvious, that damping is intimately related to a phase
difference between the fluctuations. We will see the same for sound waves.
Note that the momentum dependence of the expectation values is more general,
i.e. they also contain all higher orders in $k$.

For a sound mode with speed $v_{s}$, attenuation constant $\Gamma_{s}$ and
real momentum $k$, we notice with eq. (46) that

$-v_{s}+i\,\Gamma_{s}/2\,k=\frac{\langle\delta T^{tx}\rangle}{\langle\delta
T^{tt}\rangle}.$ (55)

In the case, where we know the speed of sound (first sound), we can simply
take the absolute value of this equation and solve for the attenuation
constant

$\Gamma_{s}=\frac{2}{k}\,\sqrt{\left|\frac{\langle\delta
T^{tx}\rangle}{\langle\delta T^{tt}\rangle}\right|^{2}-v_{s}^{2}}.$ (56)

The easiest way to compute the speed is to use the information about the
relative phases and decompose the ratio of the expectation values as

$\frac{\langle\delta T^{tx}\rangle}{\langle\delta
T^{tt}\rangle}=\frac{|\langle\delta T^{tx}\rangle|}{|\langle\delta
T^{tt}\rangle|}\,e^{i\,(\Theta_{\delta T^{tx}}-\Theta_{\delta T^{tt}})}$ (57)

and we thus find for the speed of sound according to eq. (55) and Euler’s
formula

$v_{s}=\frac{|\langle\delta T^{tx}\rangle|}{|\langle\delta
T^{tt}\rangle|}\,\cos\left(\Theta_{\delta T^{tt}}-\Theta_{\delta
T^{tx}}+\pi\right).$ (58)

Note that in the case of first sound we observe a ratio of approximately 0.7
between $|\langle\delta T^{tx}\rangle|$ and $|\langle\delta T^{tt}\rangle|$
(see for example figure 5 or 8) which yields a conformal speed of
$1/\sqrt{2}$. The attenuation constant is then simply given by

$\Gamma_{s}=\frac{2}{k}\,\frac{|\langle\delta T^{tx}\rangle|}{|\langle\delta
T^{tt}\rangle|}\,\sin\left(\Theta_{\delta T^{tt}}-\Theta_{\delta
T^{tx}}+\pi\right).$ (59)

Instead of the Ward identity between $\langle\delta T^{tt}\rangle$ and
$\langle\delta T^{tx}\rangle$, we can of course derive similar relations with
$\langle\delta T^{xx}\rangle$ and $\langle\delta T^{yy}\rangle$ using the
tracelessness of the energy-momentum tensor. This is in particular important
in the case of 4th sound, where the Ward identity associated with energy
momentum conservation no longer holds. Independently, we may also use eq. (47)
to express the speed and attenuation constant in terms of the expectation
values.

##  Appendix C Eigenvectors in Hydrodynamics

In this appendix, we discuss two simple examples within the hydrodynamic
framework. They serve the purpose of putting our results for the amplitudes in
section 4 into context and connecting them to the existing literature.

###  C.1 The purely diffusive mode in normal fluids

Let us consider the longitudinal sector of a conformal charge fluid with
conserved momentum. Following the same notations as Davison:2015taa , we may
write the Fourier transformed hydrodynamic equations of motion in the form

$\begin{pmatrix}-i\omega&ik&0\\\ ik\,\beta_{1}&\frac{k^{2}\pi
s}{4\,(\mu\,\rho+sT)}-i\omega&0\\\
k^{2}\,\alpha_{1}\,\alpha_{Q}&\frac{ik\,\rho}{\mu\rho+sT}&-i\omega-\frac{k^{2}\,\alpha_{1}\alpha_{Q}\,(\mu\rho+sT)}{\rho}\end{pmatrix}\cdot\begin{pmatrix}\delta
T^{tt}\\\ \delta T^{tx}\\\ \delta J^{t}\end{pmatrix}=\begin{pmatrix}0\\\ 0\\\
0\end{pmatrix},$ (60)

where
$\alpha_{1}=\left(\frac{\partial\mu}{\partial\epsilon}\right)_{\rho}-\frac{\mu}{T}\left(\frac{\partial
T}{\partial\varepsilon}\right)_{\rho},\,\beta_{1}=\left(\frac{\partial
P}{\partial\epsilon}\right)_{\rho},\,\alpha_{Q}=\sigma$. Note that $\delta
T^{tt}$ and $\delta T^{tx}$ are not independent variables but related by the
Ward identities (46). Instead of an independent hydrodynamic variable, $\delta
T^{tx}$ should be viewed as an auxiliary variable to reduce the generalized
eigenvalue problem to first order in $\omega$ as explained below eq. (67).

Solving the hydrodynamic equations as a generalized eigenvalue problem with
respect to the frequency $\omega$, we find a pair of sound modes and a purely
diffusive mode in agreement with the results in the main text, i.e.

$\omega_{1,2}=\pm\sqrt{\beta_{1}}k-i\,\frac{\pi\,s}{8\,(sT+\mu\rho)}k^{2}+o(k^{3})\,,\quad\omega_{3}=i\frac{\alpha_{1}\,(sT+\mu\rho)\,\alpha_{Q}}{\rho}\,k^{2}+o(k^{3})\,,$
(61)

together with their corresponding eigenvectors

$v_{1,2}=\begin{pmatrix}\frac{\mu\rho+sT}{\rho}\\\
\mp\sqrt{\beta_{1}}(\frac{\mu\rho+sT}{\rho})+i\frac{\pi s}{8\rho}k\\\
1\end{pmatrix}+o(k^{2})\,,\quad v_{3}=\begin{pmatrix}0\\\ 0\\\
1\end{pmatrix}+o(k^{2}).$ (62)

Importantly, this means that the eigenvector corresponding to the purely
diffusive mode, $v_{3}$, consists only of $\langle\delta J^{t}\rangle$ to
first order in $k$, exactly as we observe in the left plot of fig. 5.
Furthermore, we see that the ratio of the contributions to first sound from
$\langle\delta T^{tt}\rangle$ and $\langle\delta J^{t}\rangle$ is given by
$\frac{\mu\rho+sT}{\rho}$, exactly as we observe in the middle plot of fig. 5.

Let us now explain why our results are not in contradiction with those in
Davison:2015taa and how they relate to them. In Davison:2015taa it was stated
that the purely diffusive mode is carried by the incoherent current $\delta
Q^{\text{diff}}\equiv\delta J^{t}-\frac{\rho}{sT+\mu\rho}\delta T^{tt}$. To
get to this result, we have to perform a basis transformation and consider

$\begin{pmatrix}\tilde{v}_{1}&\tilde{v}_{2}&\tilde{v}_{3}\end{pmatrix}=\begin{pmatrix}v_{1}&v_{2}&v_{3}\end{pmatrix}^{-1}\begin{pmatrix}\delta
T^{tt}\\\ \delta T^{tx}\\\ \delta J^{t}\end{pmatrix}\,.$ (63)

Notice that, in this basis, the diffusive mode corresponds to
$\tilde{v}_{3}=\delta J^{t}-\frac{\rho}{sT+\mu\rho}\delta T^{tt}$. In the
language of Davison:2015taa , $\tilde{v}_{3}$ is the incoherent current. Note
that the relative coefficient in $\tilde{v}_{3}$ is exactly the inverse of the
relative coefficient between $\delta T^{tt}$ and $\delta J^{t}$ in the sound
mode. For simplicity, we define the coefficients $c_{I}$ appearing in the
decomposition of the various eigenvectors in terms of the basis of operator
fluctuations using:

$\tilde{v}_{n}\,=\,\sum_{I}\,c_{I}\,\delta\mathcal{O}^{I}\,,$ (64)

where $n$ labels the eigenvectors and $I$ the basis of operators. The
coefficients $c_{I}$ are exactly the quantities discussed in Davison:2015taa
and displayed in figure 13, 14 and 15.

Figure 13: The contribution $c_{\delta T^{tt}}$ (eq. (64)) of the $\delta
T^{tt}$ fluctuations to the incoherent current $\tilde{v}_{3}$. The blue
points are numerically extracted by inverting the matrix of the amplitudes of
the different modes presented in the main text as in eq. (63). The red dashed
line is the analytic result $\frac{\rho\,T}{sT+\mu\rho}$.

The remaining vectors in (63) are

$\displaystyle\tilde{v}_{1,2}=\frac{\rho\left(\mp\frac{i\pi
ks}{\sqrt{\beta_{1}}}+8\mu\rho+8sT\right)}{16(\mu\rho+sT)^{2}}\delta
T^{tt}\mp\frac{\rho\left(\pi^{2}k^{2}s^{2}+128\beta_{1}(\mu\rho+sT)^{2}\right)}{256\beta_{1}^{3/2}(\mu\rho+sT)^{3}}\delta
T^{tx}+\ldots$ (65)

and they correspond to the first sound mode. Using the Ward Identity relating
$\delta T^{tt}$ and $\delta T^{tx}$ and assuming that the modes corresponding
to $\tilde{v}_{1,2}$ (sound modes) are not excited, i.e.
$\tilde{v}_{1},\tilde{v}_{2}=0$, one immediately obtains the condition $\delta
T^{tt}=0$. This implies that only charge can fluctuate and leads to
$\tilde{v}_{3}=\delta J^{t}$, which matches exactly our result presented in
the main text shown in fig. 5. In a nutshell, the main physical difference
between the approach of Davison:2015taa and ours lies in considering the
response of the system to different perturbations. While in section 4 we
studied the support of the different hydrodynamic modes when only one of those
modes is switched on, ref. Davison:2015taa considered the question for a
generic perturbation in which all modes are excited. Obviously, as shown in
eq. (63), the information of the two approaches can be extracted from the same
original hydrodynamic matrix (e.g. eq. (60)). As a proof of that, we can
invert the matrix of the amplitudes of the eigenvectors, project it on the
basis of the fluctuations and obtain exactly the $\frac{\rho}{sT+\mu\rho}$
factor derived in Davison:2015taa . The results of this procedure are shown in
fig. 13.

###  C.2 The diffusive modes in a normal fluid with momentum dissipation

We can perform the same analysis for the hydrodynamic modes in the
longitudinal sector when momentum is dissipated at a rate $\sim\alpha/T$. This
section may be viewed as a supplemental discussion of the data presented in
fig. 7. With momentum dissipation, the only conserved quantities are energy
and electric charge. This leads to a $2\times 2$ hydrodynamic system in terms
of energy and charge fluctuations $\delta T^{tt},\delta J^{t}$ which, in
matrix formalism, takes (in the conventions of Hartnoll:2014lpa ) the
following form

$\begin{pmatrix}\frac{k^{2}(\chi(\bar{\alpha}\mu+\bar{\kappa})-\zeta(\mu\sigma+\bar{\alpha}T))}{\chi
c_{\mu}-T\zeta^{2}}-i\omega&\frac{k^{2}\left(\bar{\alpha}\left(c_{\mu}T-\chi\mu^{2}\right)-\bar{\kappa}(\chi\mu+T\zeta)+\mu\sigma(c_{\mu}+\mu\zeta)\right)}{\chi
c_{\mu}-T\zeta^{2}}\\\ \frac{k^{2}(\bar{\alpha}\chi-\sigma\zeta)}{\chi
c_{\mu}-T\zeta^{2}}&\frac{k^{2}(\sigma(c_{\mu}+\mu\zeta)-\bar{\alpha}(\chi\mu+T\zeta))}{\chi
c_{\mu}-T\zeta^{2}}-i\omega\end{pmatrix}\cdot\begin{pmatrix}\delta T^{tt}\\\
\delta J^{t}\end{pmatrix}=\begin{pmatrix}0\\\ 0\end{pmatrix}.$ (66)

Figure 14: Top: The relative amplitudes of $\langle\delta T^{tt}\rangle$ and
$\langle\delta J^{t}\rangle$ for the two thermoelectric diffusion modes at
$T/\mu=5$ (in blue, corresponding to the data in fig. 7) and the hydrodynamic
prediction for the ratio (dashed red line). Bottom: The normalized,
dimensionless coefficients $c_{\delta J^{t}}\,T^{2}$ and $c_{\delta
T^{tt}}\,T^{3}$ in the basis of the eigenvectors from eq. (64). In the
incoherent limit $\alpha\gg\mu$, $\delta J_{+}\to\delta J^{t}$ and $\delta
J_{-}\sim\delta T^{tt}+\gamma_{-}\,\alpha\delta J^{t}$. This implies that
$\delta J^{-}$ is purely $\delta J^{t}$ and the other mode is still a mixture
of $\delta T^{tt}$ and $\delta J^{t}$. The relative coefficient of the bottom
right plot is presented in fig. 15.

We have defined: $\sigma$ the electric conductivity,
$\chi\equiv\partial\rho/\partial\mu$ the charge susceptibility, $\bar{\kappa}$
the thermal conductivity, $\bar{\alpha}$ the thermoelectric conductivity,
$\zeta=\partial s/\partial\mu$ and the specific heat $c_{\mu}=T\partial
s/\partial T$. All these coefficients are known analytically for the linear
axion model, see for example Baggioli:2017ojd .
At this point, we can perform the same analysis as before and compute the
eigenvectors $J_{\pm}$ for the matrix eq. (66). The matrix is diagonal in the
basis of the eigenvectors $\delta J_{\pm}=a_{\pm}\,(\delta
T^{tt}+\gamma_{\pm}\,\alpha\,\delta J^{t}),$ where $a_{\pm}$ is just an overal
normalization and
$\gamma_{\pm}=-\frac{3\epsilon}{4\alpha\,\rho}\left(1\pm\sqrt{1+\frac{16\alpha^{2}\,\rho^{2}}{9\epsilon^{2}}}\right)$
Davison:2015bea , with $\epsilon$ being the energy density and $\alpha$ the
dimensionful strength of momentum dissipation. We present the amplitudes and
coefficients $c_{I}$ in fig. 14. In fig. 15, we explicitly verify that the
relative coefficient matches the analytic expression of Davison:2015bea .

Figure 15: The relative coefficient of $\delta J^{t}$ to $\delta T^{tt}$ in
the basis of the eigenvectors $J_{\pm}$ in the incoherent limit (cf. eq.
(64)). The blue points are numerically extracted by computing the coefficients
in the basis of the eigenvectors of the data presented in the top right panel
of fig. 14. The red dashed line is the analytic result for $\delta
J_{-}\sim\delta T^{tt}+\gamma_{-}\,\alpha\,\delta J^{t}$.

##  Appendix D Numerical methods

In this section, we briefly describe the numerical methods we employed to
solve the differential equations in this work following Grieninger:2020wsb ;
Ammon:2016fru ; Grieninger:2017jxz (for a detailed introduction see
Boyd1989ChebyshevAF ). Both, the background equilibrium solutions as well as
the solutions to the linearized perturbations about the equilibrium state are
obtained by means of pseudo-spectral methods. The idea of spectral methods is
to approximate the numerical solution in terms of basis functions on a
discretized grid. Throughout this work, we choose Chebychev polynomials as
basis functions and a Chebychev-Lobatto grid to discretize the radial
direction. Spectral methods solve the equations of motions globally which is a
big advantage compared to shooting methods or finite-differences where we have
to vary the initial conditions on one side of the domain until we find the
desired boundary values at the other end of the interval. More importantly,
spectral method are highly accurate and have a fast convergence rate and are
thus perfectly suited for problems in numerical holography.

Computing the static background and two-point functions in terms of pseudo-
spectral methods is fairly standard and we thus only comment on how to obtain
the quasi-normal modes. After obtaining the static background solution, we
want to investigate time- and space dependent (linearized) fluctuations about
that background to compute (a) the quasi-normal modes and (b) the retarded
Green’s functions for the transport coefficients. The coupled ordinary
differential equations for the linearized fluctuations are generally of the
form

$(\bm{C}(k)\,\omega^{2}+\bm{A}(k)\,\omega-\bm{B}(k))\,\bm{x}(\omega,k)=0,\ \ \
(\bm{F}(\omega)\,k^{2}+\bm{D}(\omega)\,k-\bm{E}(\omega))\,\bm{x}(\omega,k)=0$
(67)

where $\bm{A}$ to $\bm{F}$ are differential operators and the vector $\bm{x}$
consists of all the fluctuations
$\bm{x}=\\{h_{tt},\,h_{tx},\,h_{xx},\,h_{yy},\,a_{t},\,a_{x},\,\delta\psi_{1},\,\delta\psi_{2},\,\delta\phi_{x}\\}$.
We may recast both problems in terms of a generalized eigenvalue problem with
respect to the quasi-normal frequency $\omega$ or the momentum $k$ by
introducing auxiliary functions of the form
$\tilde{\mathfrak{f}}=\omega\,\mathfrak{f}$ and
$\tilde{\mathfrak{f}}=k\,\mathfrak{f}$ so that we find

$(\tilde{\bm{A}}\,\omega-\tilde{\bm{B}})\,\tilde{\bm{x}}_{1}=0,\ \ \
(\tilde{\bm{D}}\,k-\tilde{\bm{E}})\,\tilde{\bm{x}}_{2}=0$ (68)

where

$\tilde{\bm{x}}_{1}=\\{h_{tt},\,h_{tx},\,h_{xx}+h_{yy},\,h_{xx}-h_{yy},\,a_{t},a_{x},\,\delta\psi_{1},\,\delta\psi_{2},\,\delta\phi_{x},\,\tilde{h}_{xx}+\tilde{h}_{yy}\\}$
(69)

and

$\tilde{\bm{x}}_{2}=\\{h_{tt},\,h_{tx},\,h_{xx},\,h_{yy},\,a_{t},a_{x},\,\delta\psi_{1},\,\delta\psi_{2},\,\delta\phi_{x},\,\tilde{h}_{tt},\,\tilde{h}_{yy},\,\tilde{a}_{t},\,\delta\tilde{\psi}_{1},\,\delta\tilde{\psi}_{2}\\}\,.$
(70)

The tilde in $\tilde{\bm{x}}_{1}$ and $\tilde{\bm{x}}_{2}$ denote the
auxiliary functions.

By solving the generalized eigenvalue problem, we find for each eigenvalue
$\omega_{n}$ (or $k_{n}$) an eigenvector $\tilde{\bm{x}}$. Finally, to resolve
the strong gradients at the horizon at low temperatures accurately we use up
to 350 gridpoints in order to satisfy the equations of motion and constraints
with accuracy better than $10^{-12}$.
The Chebychev coefficients of the slowest convergent solution drop below
$10^{-13}$ in this case.

## References

  * (1) S. A. Hartnoll, C. P. Herzog and G. T. Horowitz, Building a Holographic Superconductor, Phys. Rev. Lett. 101 (2008) 031601 [0803.3295].
  * (2) S. A. Hartnoll, C. P. Herzog and G. T. Horowitz, Holographic Superconductors, JHEP 12 (2008) 015 [0810.1563].
  * (3) C. P. Herzog, Lectures on Holographic Superfluidity and Superconductivity, J.Phys. A42 (2009) 343001 [0904.1975].
  * (4) S. A. Hartnoll, Lectures on holographic methods for condensed matter physics, Class. Quant. Grav. 26 (2009) 224002 [0903.3246].
  * (5) O. Domenech, M. Montull, A. Pomarol, A. Salvio and P. J. Silva, Emergent Gauge Fields in Holographic Superconductors, JHEP 08 (2010) 033 [1005.1776].
  * (6) M. Ammon and J. Erdmenger, Gauge/Gravity Duality.  Cambridge University Press, 2015.
  * (7) J. Zaanen, Y. Liu, Y. Sun and K. Schalm, Holographic Duality in Condensed Matter Physics.  Cambridge University Press, 2015.
  * (8) S. A. Hartnoll, A. Lucas and S. Sachdev, Holographic quantum matter, 1612.07324.
  * (9) M. Baggioli, Applied Holography: A Practical Mini-Course.  SpringerBriefs in Physics. Springer, 2019.
  * (10) J. Zaanen, Why the temperature is high, Nature 430 (Jul, 2004) 512–513.
  * (11) G. T. Horowitz, Introduction to Holographic Superconductors, Lect. Notes Phys. 828 (2011) 313–347 [1002.1722].
  * (12) R.-G. Cai, L. Li, L.-F. Li and R.-Q. Yang, Introduction to Holographic Superconductor Models, Sci. China Phys. Mech. Astron. 58 (2015), no. 6 060401 [1502.00437].
  * (13) S. A. Hartnoll, C. P. Herzog and G. T. Horowitz, Holographic superconductors, Journal of High Energy Physics 2008 (dec, 2008) 015–015.
  * (14) I. Amado, M. Kaminski and K. Landsteiner, Hydrodynamics of Holographic Superconductors, JHEP 05 (2009) 021 [0903.2209].
  * (15) I. Amado, D. Arean, A. Jimenez-Alba, K. Landsteiner, L. Melgar and I. S. Landea, Holographic Type II Goldstone bosons, JHEP 07 (2013) 108 [1302.5641].
  * (16) I. Amado, D. Areán, A. Jiménez-Alba, K. Landsteiner, L. Melgar and I. Salazar Landea, Holographic Superfluids and the Landau Criterion, JHEP 02 (2014) 063 [1307.8100].
  * (17) A. Esposito, S. Garcia-Saenz and R. Penco, First sound in holographic superfluids at zero temperature, JHEP 12 (2016) 136 [1606.03104].
  * (18) S. S. Gubser and S. S. Pufu, The Gravity dual of a p-wave superconductor, JHEP 11 (2008) 033 [0805.2960].
  * (19) M. Ammon, J. Erdmenger, M. Kaminski and P. Kerner, Superconductivity from gauge/gravity duality with flavor, Phys. Lett. B 680 (2009) 516–520 [0810.2316].
  * (20) M. Ammon, J. Erdmenger, M. Kaminski and P. Kerner, Flavor Superconductivity from Gauge/Gravity Duality, JHEP 10 (2009) 067 [0903.1864].
  * (21) R.-G. Cai, S. He, L. Li and L.-F. Li, A Holographic Study on Vector Condensate Induced by a Magnetic Field, JHEP 12 (2013) 036 [1309.2098].
  * (22) R.-G. Cai, L. Li and L.-F. Li, A Holographic P-wave Superconductor Model, JHEP 01 (2014) 032 [1309.4877].
  * (23) M. Ammon, J. Erdmenger, V. Grass, P. Kerner and A. O’Bannon, On Holographic p-wave Superfluids with Back-reaction, Phys. Lett. B 686 (2010) 192–198 [0912.3515].
  * (24) J. Erdmenger, P. Kerner and H. Zeller, Transport in Anisotropic Superfluids: A Holographic Description, JHEP 01 (2012) 059 [1110.0007].
  * (25) J. Erdmenger, D. Fernandez and H. Zeller, New Transport Properties of Anisotropic Holographic Superfluids, JHEP 04 (2013) 049 [1212.4838].
  * (26) J.-W. Chen, Y.-J. Kao, D. Maity, W.-Y. Wen and C.-P. Yeh, Towards A Holographic Model of D-Wave Superconductors, Phys. Rev. D 81 (2010) 106008 [1003.2991].
  * (27) F. Benini, C. P. Herzog, R. Rahman and A. Yarom, Gauge gravity duality for d-wave superconductors: prospects and challenges, JHEP 11 (2010) 137 [1007.1981].
  * (28) K.-Y. Kim and M. Taylor, Holographic d-wave superconductors, JHEP 08 (2013) 112 [1304.6729].
  * (29) X.-H. Ge, S. F. Tu and B. Wang, d-Wave holographic superconductors with backreaction in external magnetic fields, JHEP 09 (2012) 088 [1209.4272].
  * (30) T. Andrade and S. A. Gentle, Relaxed superconductors, JHEP 06 (2015) 140 [1412.6521].
  * (31) K.-Y. Kim, K. K. Kim and M. Park, A Simple Holographic Superconductor with Momentum Relaxation, JHEP 04 (2015) 152 [1501.00446].
  * (32) M. Baggioli and M. Goykhman, Phases of holographic superconductors with broken translational symmetry, JHEP 07 (2015) 035 [1504.05561].
  * (33) M. Baggioli and M. Goykhman, Under The Dome: Doped holographic superconductors with broken translational symmetry, JHEP 01 (2016) 011 [1510.06363].
  * (34) T. Andrade and B. Withers, A simple holographic model of momentum relaxation, JHEP 05 (2014) 101 [1311.5157].
  * (35) M. Baggioli, K.-Y. Kim, L. Li and W.-J. Li, Holographic Axion Model: a simple gravitational tool for quantum matter, 2101.01892.
  * (36) C. P. Herzog, N. Lisker, P. Surowka and A. Yarom, Transport in holographic superfluids, JHEP 08 (2011) 052 [1101.3330].
  * (37) C. P. Herzog and A. Yarom, Sound modes in holographic superfluids, Phys. Rev. D 80 (2009) 106002 [0906.4810].
  * (38) A. Yarom, Fourth sound of holographic superfluids, JHEP 07 (2009) 070 [0903.1353].
  * (39) R. J. Donnelly, The two-fluid theory and second sound in liquid helium, Physics Today 62 (2009), no. 10 34–39 [https://doi.org/10.1063/1.3248499].
  * (40) L. Tisza, Transport phenomena in helium ii, Nature 141 (1938), no. 3577 913–913.
  * (41) L. Landau, Theory of the superfluidity of helium ii, Physical Review 60 (1941), no. 4 356.
  * (42) A. Schmitt, Introduction to Superfluidity: Field-theoretical approach and applications, 1404.1284.
  * (43) J. D. Jackson, Classical electrodynamics.  Wiley, New York, NY, 3rd ed. ed., 1999.
  * (44) I. Amado, C. Hoyos-Badajoz, K. Landsteiner and S. Montero, Absorption lengths in the holographic plasma, JHEP 09 (2007) 057 [0706.2750].
  * (45) I. Amado, C. Hoyos-Badajoz, K. Landsteiner and S. Montero, Hydrodynamics and beyond in the strongly coupled N=4 plasma, JHEP 07 (2008) 133 [0805.2570].
  * (46) D. Forcella, A. Mezzalira and D. Musso, Electromagnetic response of strongly coupled plasmas, JHEP 11 (2014) 153 [1404.4048].
  * (47) M. Blake, A. Donos and D. Tong, Holographic Charge Oscillations, JHEP 04 (2015) 019 [1412.2003].
  * (48) J. Sonner and B. Withers, Universal spatial structure of nonequilibrium steady states, Phys. Rev. Lett. 119 (2017), no. 16 161603 [1705.01950].
  * (49) S. Grozdanov, K. Schalm and V. Scopelliti, Black hole scrambling from hydrodynamics, Phys. Rev. Lett. 120 (2018), no. 23 231601 [1710.00921].
  * (50) S. Grozdanov, P. K. Kovtun, A. O. Starinets and P. Tadić, Convergence of the Gradient Expansion in Hydrodynamics, Phys. Rev. Lett. 122 (2019), no. 25 251601 [1904.01018].
  * (51) E. P. Gross, Hydrodynamics of a superfluid condensate, Journal of Mathematical Physics 4 (1963), no. 2 195–207 [https://doi.org/10.1063/1.1703944].
  * (52) S. J. Putterman, Superfluid hydrodynamics, vol. 3.  Jan., 1974.
  * (53) A. Nicolis, Low-energy effective field theory for finite-temperature relativistic superfluids, 1108.2513.
  * (54) D. T. Son, Low-energy quantum effective action for relativistic superfluids, hep-ph/0204199.
  * (55) J. Bhattacharya, S. Bhattacharyya, S. Minwalla and A. Yarom, A Theory of first order dissipative superfluid dynamics, JHEP 05 (2014) 147 [1105.3733].
  * (56) J. O. Andersen, Effective field theory for Goldstone bosons in nonrelativistic superfluids, cond-mat/0209243.
  * (57) J. Sonner and B. Withers, A gravity derivation of the Tisza-Landau Model in AdS/CFT, Phys. Rev. D 82 (2010) 026001 [1004.2707].
  * (58) M. A. Escobedo and C. Manuel, Effective field theory and dispersion law of the phonons of a non-relativistic superfluid, Phys. Rev. A 82 (2010) 023614 [1004.2567].
  * (59) L. Berezhiani, On Effective Theory of Superfluid Phonons, Phys. Lett. B 805 (2020) 135451 [2001.08696].
  * (60) M. J. Landry, Second sound and non-equilibrium effective field theory, 2008.11725.
  * (61) M. J. Landry, Higher-form symmetries and ’t Hooft anomalies in non-equilibrium systems, 2101.02210.
  * (62) E. Pajer and D. Stefanyszyn, Symmetric Superfluids, JHEP 06 (2019) 008 [1812.05133].
  * (63) L. V. Delacrétaz, D. M. Hofman and G. Mathys, Superfluids as Higher-form Anomalies, SciPost Phys. 8 (2020) 047 [1908.06977].
  * (64) P. Kovtun, Lectures on hydrodynamic fluctuations in relativistic theories, J. Phys. A 45 (2012) 473001 [1205.5040].
  * (65) D. T. Son and A. O. Starinets, Hydrodynamics of r-charged black holes, JHEP 0603 (2006) 052 [hep-th/0601157].
  * (66) S. A. Hartnoll, Theory of universal incoherent metallic transport, Nature Phys. 11 (2015) 54 [1405.3651].
  * (67) P. Kovtun, Fluctuation bounds on charge and heat diffusion, J. Phys. A 48 (2015), no. 26 265002 [1407.0690].
  * (68) J. R. Pellam, Wave transmission and reflection phenomena in liquid helium ii, Phys. Rev. 73 (Mar, 1948) 608–617.
  * (69) K. R. Atkins, Third and fourth sound in liquid helium ii, Phys. Rev. 113 (Feb, 1959) 962–965.
  * (70) M. Ammon, S. Grieninger, J. Hernandez, M. Kaminski, R. Koirala, J. Leiber and J. Wu, Chiral hydrodynamics in strong external magnetic fields, JHEP 04 (2021) 078 [2012.09183].
  * (71) M. Baggioli.  PhD thesis, Barcelona U., 2016.  1610.02681.
  * (72) I. Božović, X. He, J. Wu and A. Bollinger, Dependence of the critical temperature in overdoped copper oxides on superfluid density, Nature 536 (August, 2016) 309—311.
  * (73) B. Goutéraux and E. Mefford, Normal charge densities in quantum critical superfluids, Phys. Rev. Lett. 124 (2020), no. 16 161604 [1912.08849].
  * (74) G. T. Horowitz and M. M. Roberts, Zero Temperature Limit of Holographic Superconductors, JHEP 11 (2009) 015 [0908.3677].
  * (75) B. Goutéraux and E. Mefford, Non-vanishing zero-temperature normal density in holographic superfluids, JHEP 11 (2020) 091 [2008.02289].
  * (76) I. Amado, M. Kaminski and K. Landsteiner, Hydrodynamics of Holographic Superconductors, JHEP 05 (2009) 021 [0903.2209].
  * (77) A. Donos, P. Kailidis and C. Pantelidou, Dissipation in Holographic Superfluids, 2107.03680.
  * (78) K. A. Shapiro and I. Rudnick, Experimental determination of the fourth sound velocity in helium ii, Phys. Rev. 137 (Mar, 1965) A1383–A1391.
  * (79) E. Berti, V. Cardoso and A. O. Starinets, Quasinormal modes of black holes and black branes, Class. Quant. Grav. 26 (2009) 163001 [0905.2975].
  * (80) H.-P. Nollert, Quasinormal modes: the characteristic sound of black holes and neutron stars, Classical and Quantum Gravity 16 (nov, 1999) R159–R216.
  * (81) I. Amado, C. Hoyos-Badajoz, K. Landsteiner and S. Montero, Residues of correlators in the strongly coupled N=4 plasma, Phys. Rev. D77 (2008) 065004 [0710.4458].
  * (82) R. A. Davison and N. K. Kaplis, Bosonic excitations of the $AdS_{4}$ Reissner-Nordstrom black hole, JHEP 1112 (2011) 037 [1111.0660].
  * (83) M. G. Alford, S. K. Mallavarapu, A. Schmitt and S. Stetina, Role reversal in first and second sound in a relativistic superfluid, Phys. Rev. D89 (2014), no. 8 085005 [1310.5953].
  * (84) R. A. Davison, B. Goutéraux and S. A. Hartnoll, Incoherent transport in clean quantum critical metals, JHEP 10 (2015) 112 [1507.07137].
  * (85) P. Kovtun, Fluctuation bounds on charge and heat diffusion, J. Phys. A48 (2015), no. 26 265002 [1407.0690].
  * (86) R. A. Davison and B. Goutéraux, Dissecting holographic conductivities, JHEP 09 (2015) 090 [1505.05092].
  * (87) M. Baggioli and W.-J. Li, Diffusivities bounds and chaos in holographic Horndeski theories, JHEP 07 (2017) 055 [1705.01766].
  * (88) A. Behrle, T. Harrison, J. Kombe, K. Gao, M. Link, J.-S. Bernier, C. Kollath and M. Köhl, Higgs mode in a strongly interacting fermionic superfluid, Nature Physics 14 (May, 2018) 781–785.
  * (89) M. J. Bhaseen, J. P. Gauntlett, B. D. Simons, J. Sonner and T. Wiseman, Holographic Superfluids and the Dynamics of Symmetry Breaking, Phys. Rev. Lett. 110 (2013), no. 1 015301 [1207.4194].
  * (90) W.-J. Li, Y. Tian and H.-b. Zhang, Periodically Driven Holographic Superconductor, JHEP 07 (2013) 030 [1305.1600].
  * (91) D. Pekker and C. Varma, Amplitude/higgs modes in condensed matter physics, Annual Review of Condensed Matter Physics 6 (Mar, 2015) 269–297.
  * (92) M. Baggioli and O. Pujolas, Electron-Phonon Interactions, Metal-Insulator Transitions, and Holographic Massive Gravity, Phys. Rev. Lett. 114 (2015), no. 25 251602 [1411.1003].
  * (93) T. Andrade, M. Baggioli, A. Krikun and N. Poovuttikul, Pinning of longitudinal phonons in holographic spontaneous helices, JHEP 02 (2018) 085 [1708.08306].
  * (94) L. Alberte, M. Ammon, M. Baggioli, A. Jimenez and O. Pujolas, Black hole elasticity and gapped transverse phonons in holography, JHEP 01 (2018) 129 [1708.08477].
  * (95) A. Amoretti, D. Areán, R. Argurio, D. Musso and L. A. Pando Zayas, A holographic perspective on phonons and pseudo-phonons, JHEP 05 (2017) 051 [1611.09344].
  * (96) A. Amoretti, D. Arean, B. Gouteraux and D. Musso, Effective holographic theory of charge density waves, 1711.06610.
  * (97) A. Amoretti, D. Arean, B. Gouteraux and D. Musso, Gapless and gapped holographic phonons, JHEP 01 (2020) 058 [1910.11330].
  * (98) M. Ammon, M. Baggioli, S. Gray and S. Grieninger, Longitudinal Sound and Diffusion in Holographic Massive Gravity, JHEP 10 (2019) 064 [1905.09164].
  * (99) M. Ammon, M. Baggioli, S. Gray, S. Grieninger and A. Jain, On the Hydrodynamic Description of Holographic Viscoelastic Models, Phys. Lett. B 808 (2020) 135691 [2001.05737].
  * (100) M. Baggioli and S. Grieninger, Zoology of solid & fluid holography — Goldstone modes and phase relaxation, JHEP 10 (2019) 235 [1905.09488].
  * (101) M. Ammon, M. Baggioli and A. Jimenez-Alba, A Unified Description of Translational Symmetry Breaking in Holography, JHEP 09 (2019) 124 [1904.05785].
  * (102) A. Amoretti, D. Areán, B. Goutéraux and D. Musso, Universal relaxation in a holographic metallic density wave phase, Phys. Rev. Lett. 123 (2019), no. 21 211602 [1812.08118].
  * (103) A. Amoretti, D. Areán, B. Goutéraux and D. Musso, Diffusion and universal relaxation of holographic phonons, 1904.11445.
  * (104) M. Baggioli and M. Landry, Effective Field Theory for Quasicrystals and Phasons Dynamics, 2008.05339.
  * (105) M. Baggioli, S. Grieninger and L. Li, Magnetophonons & type-B Goldstones from Hydrodynamics to Holography, JHEP 09 (2020) 037 [2005.01725].
  * (106) M. Baggioli, Homogeneous Holographic Viscoelastic Models & Quasicrystals, 2001.06228.
  * (107) R. A. Davison, L. V. Delacrétaz, B. Goutéraux and S. A. Hartnoll, Hydrodynamic theory of quantum fluctuating superconductivity, Phys. Rev. B94 (2016), no. 5 054502 [1602.08171]. [Erratum: Phys. Rev.B96,no.5,059902(2017)].
  * (108) S. Grozdanov, P. K. Kovtun, A. O. Starinets and P. Tadić, The complex life of hydrodynamic modes, Journal of High Energy Physics 2019 (Nov, 2019) 97.
  * (109) B. Withers, Short-lived modes from hydrodynamic dispersion relations, JHEP 06 (2018) 059 [1803.08058].
  * (110) N. Abbasi and S. Tahery, Complexified quasinormal modes and the pole-skipping in a holographic system at finite chemical potential, 2007.10024.
  * (111) N. Abbasi and M. Kaminski, Constraints on quasinormal modes and bounds for critical points from pole-skipping, JHEP 03 (2021) 265 [2012.15820].
  * (112) A. Jansen and C. Pantelidou, Quasinormal modes in charged fluids at complex momentum, 2007.14418.
  * (113) M. Baggioli, How small hydrodynamics can go, 2010.05916.
  * (114) D. Arean, R. A. Davison, B. Goutéraux and K. Suzuki, Hydrodynamic diffusion and its breakdown near AdS2 fixed points, 2011.12301.
  * (115) N. Wu, M. Baggioli and W.-J. Li, On the universality of AdS2 diffusion bounds and the breakdown of linearized hydrodynamics, 2102.05810.
  * (116) M. Baggioli, U. Gran and M. Tornsö, Collective modes of polarizable holographic media in magnetic fields, 2102.09969.
  * (117) S. Grozdanov, A. Starinets and P. Tadić, Hydrodynamic dispersion relations at finite coupling, 2104.11035.
  * (118) S. Grozdanov, Bounds on transport from univalence and pole-skipping, 2008.00888.
  * (119) H.-S. Jeong, K.-Y. Kim and Y.-W. Sun, The breakdown of magneto-hydrodynamics near AdS2 fixed point and energy diffusion bound, 2105.03882.
  * (120) S. Grieninger and A. Shukla, Second order equilibrium transport in strongly coupled $\mathcal{N}=4$ supersymmetric $SU(N_{c})$ Yang-Mills plasma via holography, JHEP 08 (2021) 108 [2105.08673].
  * (121) P. Kovtun and A. Shukla, Kubo formulas for thermodynamic transport coefficients, JHEP 10 (2018) 007 [1806.05774].
  * (122) D. Arean, P. Basu and C. Krishnan, The Many Phases of Holographic Superfluids, JHEP 10 (2010) 006 [1006.5165].
  * (123) D. Arean, M. Bertolini, C. Krishnan and T. Prochazka, Type IIB Holographic Superfluid Flows, JHEP 03 (2011) 008 [1010.5777].
  * (124) D. Arean, M. Bertolini, C. Krishnan and T. Prochazka, Quantum Critical Superfluid Flows and Anisotropic Domain Walls, JHEP 09 (2011) 131 [1106.1053].
  * (125) S. Lan, H. Liu, Y. Tian and H. Zhang, Landau Instability and soliton formations, 2010.06232.
  * (126) D. Arean, A. Farahi, L. A. Pando Zayas, I. S. Landea and A. Scardicchio, Holographic superconductor with disorder, Phys. Rev. D 89 (2014), no. 10 106003 [1308.1920].
  * (127) D. Arean, L. A. Pando Zayas, I. S. Landea and A. Scardicchio, Holographic disorder driven superconductor-metal transition, Phys. Rev. D 94 (2016), no. 10 106003 [1507.02280].
  * (128) R. McNees, Conventions, definitions, identities, and formulas, 2020.  http://jacobi.luc.edu/Useful.html, [Accessed 12-June-2020].
  * (129) K. K. Kim, M. Park and K.-Y. Kim, Ward identity and Homes’ law in a holographic superconductor with momentum relaxation, JHEP 10 (2016) 041 [1604.06205].
  * (130) S. Grieninger, Non-equilibrium dynamics in Holography.  PhD thesis, Jena U., 2020.  2012.10109.
  * (131) M. Ammon, S. Grieninger, A. Jimenez-Alba, R. P. Macedo and L. Melgar, Holographic quenches and anomalous transport, JHEP 09 (2016) 131 [1607.06817].
  * (132) S. Grieninger, Holographic quenches and anomalous transport, Master’s thesis, Jena U., TPI, 2016, [1711.08422].
  * (133) J. P. Boyd, Chebyshev and Fourier Spectral Methods.  Dover Publications Inc., 2003.