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, Spain

# Pseudospectra of Holographic Quasinormal Modes

Daniel Areán 1,2 David García Fariña 1 Karl Landsteiner<EMAIL_ADDRESS><EMAIL_ADDRESS><EMAIL_ADDRESS>
###### Abstract

Quasinormal modes and frequencies are the eigenvectors and eigenvalues of a
non-Hermitian differential operator. They hold crucial significance in the
physics of black holes. The analysis of quasinormal modes of black holes in
asymptotically Anti-de Sitter geometries plays also a key role in the study of
strongly coupled quantum many-body systems via gauge/gravity duality. In
contrast to normal Sturm-Liouville operators, the spectrum of non-Hermitian
(and non-normal) operators generally is unstable under small perturbations.
This research focuses on the stability analysis of the spectrum of quasinormal
frequencies pertaining to asymptotically planar AdS black holes, employing
pseudospectrum analysis. Specifically, we concentrate on the pseudospectra of
scalar and transverse gauge fields, shedding light on their relevance within
the framework of gauge/gravity duality.

††preprint: IFT-UAM/CSIC-23-091

##  1 Introduction

Quasinormal modes are the solutions of linear differential equations on an
open domain. The openness is encoded in outgoing boundary conditions. Such
boundary conditions arise naturally in many branches of physics from the study
of black holes to optics. Quasinormal modes can be thought of as the
eigenmodes of linear differential operators. Due to the open boundary
conditions the differential operator is non-Hermitian and the eigenvalues are
complex.

In gravitational physics, quasinormal modes are of utmost importance in the
theory of black holes. Famously, a black hole horizon acts as a perfectly
absorbing membrane. Therefore, solving a wave equation in a black hole
background leads naturally to outgoing boundary conditions. In asymptotically
flat space, radiation is also outgoing at infinity. The quasinormal modes that
arise in asymptotically flat spacetimes are supposed to describe the late time
ringdown of perturbed black holes. They appear for example in black hole
mergers and are of highest interest for gravitational wave physics
Nollert:1999ji ; Kokkotas:1999bd ; Berti:2009kk ; Franchini:2023eda ;
Jaramillo:2022kuv .

The importance of quasinormal modes in asymptotically Anti de-Sitter (AdS)
spaces derives from gauge/gravity duality Maldacena:1997re ; Aharony:1999ti ;
BigJ ; zaanen2015holographic ; Hartnoll:2018xxg . Gauge/gravity duality is
conjectured to describe certain strongly coupled quantum many-body systems.
While up to date no concrete physical system has been found that indeed is
modelled in all detail by gauge/gravity duality, it has lead to important
insights into areas such as hydrodynamics and transport theory in the
relativistic regime Policastro:2001yc ; Baier:2007ix ; Bhattacharyya:2007vjd .
Some key results are the extremely low specific shear viscosity of holographic
models of the quark gluon plasma Kovtun:2004de  or phase transitions towards
superconducting states Gubser:2008px ; Hartnoll:2008vx  as well as strongly
coupled quantum critical phases Herzog:2007ij . Quasinormal modes are a key
ingredient in this development. They can be understood as the poles of the
retarded Green’s functions of the strongly coupled quantum many-body system
Horowitz:1999jd ; Birmingham:2001pj ; Kovtun:2005ev . The shear viscosity can
be read off from a quasinormal mode corresponding to momentum diffusion
Policastro:2002se1  and second order phase transitions arise as a quasinormal
frequency moves into the upper half plane thus indicating an instability of
the ground state Herzog:2009ci ; Amado:2009ts .

In particular the modern way of understanding relativistic hydrodynamics is
very much influenced by the properties of quasinormal modes in asymptotically
AdS black holes. Relativistic hydrodynamics can be understood as an effective
field theory organised in a derivative expansion Kovtun:2012rj . As an
effective field theory it has a limited range of validity. This limit can be
set as the wavelength at which the first non-hydrodynamic mode decays at a
slower rate than the hydrodynamic mode Amado:2008ji . More recently it is
conjectured that the limit is set by the radius in the complex plane at which
the lowest hydrodynamic quasinormal mode collides with the first non-
hydrodynamic mode Withers:2018srf ; Grozdanov:2019 ; GrozdanovKovtun:2019 .

A question that so far has not been investigated in the context of
gauge/gravity duality is the spectral stability of the quasinormal frequencies
under small perturbations of the background geometry. We will refer to this
spectral instability simply as instability. Since in this work we are not
considering eigenvalues lying on the upper half of the complex plane (which
would denote instability of the background), no confusion should arise from
this wording.

In asymptotically flat space studies on the stability of quasinormal
frequencies have a long history Nollert:1996rf ; Nollert:1998ys ;
Aguirregabiria:1996  and recently the issue has been taken up in
Jaramillo:2020tuu ; Destounis:2021lum ; Sarkar:2023rhp ; Cheung:2021bol ;
Berti:2022xfj ; Konoplya:2022pbc ; Courty:2023rxk ; Daghigh:2020jyk ;
Qian:2020cnz ; alsheikh:tel-04116011 ; Destounis:2023ruj ; Torres:2023nqg
(see also Boyanov:2022ark  for horizonless compact objects). The key new
development of these recent studies is the application of the method of
pseudospectra Trefethen:2005 ; Sjostrand:2019 ; davies:2007 .

The eigenvalues of Hermitian operators are stable in a technical sense. In
fact, this stability is what makes perturbation theory successful in quantum
physics. On the other hand, the eigenvalues of non-Hermitian and non-normal
operators generically are unstable. This means they can be displaced large
distances by even small perturbations of the original operator. Here stability
and smallness of a perturbation have precise meanings and we will give exact
definitions later on.

In the case of quasinormal modes the spectral instability appears as a
consequence of the non-conservative nature of the system. This non-
conservative nature is associated with the damping of the fluctuations as they
fall into the black hole and radiate away to infinity in asymptotically flat
spacetime. AdS at large distances from the black hole horizon acts as a
confining box but it is still true that fluctuations are damped as they fall
into the black hole.

Thus the damping of fluctuations in black hole geometries is independent of
the asymptotics of the spacetime. It is then reasonable to expect that such
quasinormal modes and frequencies are also inherently unstable as occurs in
asymptotically flat space. The motivation for investigating this question is
twofold. Firstly, it would serve as increasing evidence of the independence of
the instability on the asymptotic behavior of the spacetime. And secondly, it
would shed a completely new light on the behavior of strongly coupled quantum
many-body systems modelled by gauge/gravity duality. The instability of the
quasinormal mode spectrum would indicate that also the excitation spectrum of
the dual quantum field theories would be unstable under perturbations. In
turn, this might rise [KL: raise] questions about the robustness of
conclusions for quantum many-body physics drawn from the behavior of
quasinormal modes. In particular, transport properties in these systems could
potentially be heavily model dependent, with a small deformation of the
original theory leading to very different transport coefficients. Let us add
here that despite the theoretically well established instability of the
quasinormal mode spectrum in the asymptotically flat case, extraction of
overtones from gravitational wave signals have been reported in Isi:2019aib ;
Giesler:2019uxc ; Capano:2021etf ; Capano:2022zqm .

To probe the stability of the quasinormal mode frequency (QNF) we follow the
example of Jaramillo:2020tuu  and turn to pseudospectrum analysis
Trefethen:2005 . Broadly speaking, pseudospectrum provides insight on how much
perturbations of a given size can displace the spectrum of an operator.
Accordingly, once the physically motivated notion of size is defined,
pseudospectra serve as a powerful tool to discuss spectral stability or lack
thereof.

In this work we initiate a systematic study of the stability of quasinormal
frequencies in AdS. We focus on two very simple cases: a real scalar field and
a transverse gauge field in a Schwarzschild AdS4+1 black brane (SAdS4+1). By
gauge/gravity duality these correspond to a scalar operator and a conserved
$R$-current of the strongly coupled maximally symmetric $\mathcal{N}=4$ Yang-
Mills theory Maldacena:1997re . Our choice of fluctuations is motivated by two
main reasons. Firstly, we consider the SAdS4+1 background because it is the
archetypical example for gauge/gravity duality. Secondly, we choose a real
scalar and a transverse gauge field for their simplicity and their lack of a
hydro mode. Hydro modes dominate the low energy spectrum of the dual QFT as
their corresponding QNF goes to zero as the momentum goes to zero.
Consequently, analyzing the stability of hydro modes is akin to studying the
stability of the hydrodynamic description of the dual quantum field theory
Kovtun:2005ev . As the study of hydrodynamic modes involves additional
technical complexities due to the presence of gauge symmetries and the
corresponding constraints, we will leave this topic for a follow-up work.

The organisation of this paper is as follows. Since the method of
pseudospectra is relatively new and perhaps little known in the context of
gauge/gravity duality we spend the first few sections reviewing some of the
underlying mathematics. Here we draw heavily from Trefethen:2005  and also
from a previous work on asymptotically flat space Jaramillo:2020tuu . In
particular we define the pseudospectrum and review its most important
properties in section 2.

In section 3 we focus specifically on the application of pseudospectra to
quasinormal modes. We point out the importance of choosing the proper
coordinates. Following Warnick:2013hba  we chose a coordinate slice which we
call regular. These coordinates have the important property that the equal-
time slices are spacelike outside the horizon but coincide with infalling
Eddington-Finkelstein coordinates at the horizon.111We note that infalling
coordinates are very well suited to impose outgoing boundary conditions. We
discuss the importance of the choice of norm and give a physical reasoning why
the chosen energy norm is the relevant one in the holographic context. Indeed,
the energy in regular coordinates is not conserved due to outflow on the
horizon. We conclude this section by discussing the holographic interpretation
of the instability of the quasinormal modes.

In section 4 we construct the specific differential operators that describe
the fluctuations of scalar and transverse gauge fields in an AdS black hole
background. We show explicitly that the non-Hermiticity is concentrated on the
horizon. As we discuss in detail this leads to a clear picture of the local
origin of the non-Hermiticity and non-normality.

In section 5 we first outline our numerical methods which are based on
pseudospectral methods Trefethen:2000 ; boyd . Next, we introduce the
selective pseudospectrum which tests the stability under random local
potential perturbations. And finally, to further explore the aforementioned
stability, we construct specific potentials to test some relevant regimes.

Section 6 is the core of the paper and contains the results of our
pseudospectrum calculations. We mostly concentrate on the results for the
energy norm but, to illustrate the norm dependence, we also present the result
for the $L^{2}$ norm in appendix C. Indeed it turns out that the
pseudospectrum for both the scalar field and the transverse vector show the
typical open contour lines indicating instability.

In section 7 we summarize our findings, present our conclusions and give an
outlook on further studies in particular the highly interesting case of
hydrodynamic modes.

Some of the more technical details on the calculation of the adjoint
differential operators and the numerical methods are presented in the
appendices A, B and E. The numerical values of the quasinormal frequencies are
presented in appendix D.

##  2 Pseudospectrum and Stability

In this section we collect relevant theorems and definitions of the
pseudospectrum of linear operators. We next specialize them to finite
dimensional operators (matrices) and illustrate some key notions through the
example of a $2\times 2$ matrix. For additional details, we refer the reader
to Trefethen:2005 .

Given a closed linear operator $L$ acting on a Hilbert space $\mathcal{H}$
with domain $\mathcal{D}(L)$, its spectrum $\sigma\left(L\right)$ is defined
as the set of points in the complex plane, $z$, where its resolvent:222Recall
that a closed linear operator $L$ is a linear operator such that if a sequence
$\\{x_{n}\\}$ in $\mathcal{H}$ converges to $x$, then the sequence
$\\{Lx_{n}\\}$ converges to $y$, where $Lx=y$. From now on, when talking about
operators, we implicitly assume that they are both closed and linear.

$\displaystyle\mathcal{R}(z;L)=(L-z)^{-1}\,,$ (1)

does not exist as a bounded operator on $\mathcal{H}$.

Within $\sigma\left(L\right)$, we define the eigenvalues $\\{\lambda_{i}\\}$
through the eigenvalue equation:

$Lu_{i}=\lambda_{i}u_{i}\,,$ (2)

where $u_{i}\in\mathcal{D}(L)$ is the eigenvector corresponding to the
eigenvalue $\lambda_{i}$. Note that with our definitions, the eigenvalues
correspond to isolated points in the spectrum (regions of dimension 0).
Nonetheless, it is important to stress that, in general, the spectrum also
contains regions of higher dimension, such as branch cuts (regions of
dimension 1). That said, Warnick:2013hba  proved that in the case of
Schwarzschild AdS black holes the spectrum of quasinormal frequencies indeed
consists of a countable discrete subset of the complex plane.

A particularly important property of eigenvalues is that as long as the
operator is self-adjoint (or, in more physical terms, the system is
conservative), the spectral theorem ensures that if we perturb the system with
a bounded operator of size $\varepsilon$ the eigenvalues of the perturbed
operator cannot suffer a displacement greater than $\varepsilon$
kato2013perturbation ; Courant:1989 .333The notion of size here is quite
nontrivial and we shall discuss it in greater detail later. In general, we
will consider that the size is given by the norm of the operator, which is
inherited from the function norm in $\mathcal{H}$. More generally, this
property holds for any normal operator $A$ satisfying
$\left[A,A^{\dagger}\right]=0$, with $A^{\dagger}$ the adjoint of $A$.
Physically, this ensures that any statement made relying on eigenvalue
analysis is guaranteed to be right up to corrections of the same size as the
error implicitly made when assuming a particular model.

However, many physical systems are inherently non-conservative. In these
systems, normality of the relevant operators is not guaranteed and small
perturbations could potentially alter the spectrum in a significant manner. As
any description of a physical system is always a model and consequently has
some inherent error associated with it, one can conclude Trefethen:2005  that
eigenvalue analysis is insufficient when dealing with non-conservative systems
defined through non-self-adjoint (and potentially non-normal) operators as
eigenvalues may not be stable.

In order to characterize the stability of eigenvalues we need to introduce the
notion of $\varepsilon$-pseudospectrum, which can be defined in three
mathematically equivalent ways Trefethen:2005 :

######  Def. 2.1 (Resolvent norm approach).

Given a closed linear operator $L$ acting on a Hilbert space $\mathcal{H}$
with domain $\mathcal{D}(L)$, and $\varepsilon>0$, the
$\varepsilon$-pseudospectrum $\sigma_{\varepsilon}(L)$ is

$\sigma_{\varepsilon}(L)=\\{z\in\mathbb{C}:\lVert\mathcal{R}(z;L)\rVert>1/\varepsilon\\}\,,$
(3)

with the convention $\lVert\mathcal{R}(z;L)\rVert=\infty$ for $z\in\sigma(L)$.

######  Def. 2.2 (Perturbative approach).

Given a closed linear operator $L$ acting on a Hilbert space $\mathcal{H}$
with domain $\mathcal{D}(L)$, and $\varepsilon>0$, the
$\varepsilon$-pseudospectrum $\sigma_{\varepsilon}(L)$ is

$\sigma_{\varepsilon}(L)=\\{z\in\mathbb{C},\exists V,\rVert
V\lVert<\varepsilon:z\in\sigma(L+V)\\}\,.$ (4)

######  Def. 2.3 (Pseudoeigenvalue approach).

Given a closed linear operator $L$ acting on a Hilbert space $\mathcal{H}$
with domain $\mathcal{D}(L)$, and $\varepsilon>0$, the
$\varepsilon$-pseudospectrum $\sigma_{\varepsilon}(L)$ is

$\sigma_{\varepsilon}(L)=\\{z\in\mathbb{C},\exists
u^{\varepsilon}\in\mathcal{D}(L):\;\rVert(L-z)u^{\varepsilon}\lVert<\varepsilon\lVert
u^{\varepsilon}\rVert\\}\,,$ (5)

where $u^{\varepsilon}$ is a $\varepsilon$-pseudoeigenvector with
$\varepsilon$-pseudoeigenvalue $z$.

Note that, contrary to the spectrum, the pseudospectrum depends on the
operator norm as it needs a notion of what constitutes a large or small
perturbation to quantify stability.

######  Def. 2.4 (Operator norm).

Given a bounded linear operator $V$ acting on a Hilbert space $\mathcal{H}$
equipped with a function norm $\lVert\cdot\rVert$, we define its norm $\lVert
V\rVert$ as:

$\norm{V}=\underset{u\in\mathcal{H}}{\max}\frac{\norm{Vu}}{\norm{u}}\,.$ (6)

Definition 2.2 corresponds to the physical intuition we were seeking: the
$\varepsilon$-pseudospectrum constitutes the maximal region containing all
possible displacements of the eigenvalues under perturbations of size
$\varepsilon$. It is then quite natural to represent the pseudospectrum as a
contour map indicating the boundaries of these regions for multiple values of
$\varepsilon$ (see e.g. figure 1 for one such a map corresponding to a
$2\times 2$ matrix).

On the other hand, definition 2.1, despite lacking such a clear physical
interpretation, is significantly more powerful as it allows to establish a few
relevant theorems whose proofs can be found in chapter 4 of Trefethen:2005 :

######  Thm. 2.1.

Given a closed linear operator $L$ acting on a Hilbert space $\mathcal{H}$
with domain $\mathcal{D}(L)$, the $\varepsilon$-pseudospectrum
$\sigma_{\varepsilon}(L)$ is a nonempty open set of $\,\mathbb{C}$, and any
bounded connected component of $\sigma_{\varepsilon}(L)$ has a nonempty
intersection with $\sigma(L)$.

######  Thm. 2.2.

Given a closed linear operator $L$ acting on a Hilbert space $\mathcal{H}$
with domain $\mathcal{D}(L)$, the pseudospectra are strictly nested supersets
of the spectrum: $\cap_{\varepsilon>0}\,\sigma_{\varepsilon}(L)=\sigma(L)$,
and conversely, for any $\delta>0$,
$\sigma_{\varepsilon+\delta}(L)\supseteq\sigma_{\varepsilon}(L)+\Delta_{\delta}$
with $\Delta_{\delta}$ the open disk of radius $\delta$.

######  Thm. 2.3.

Given a closed linear operator $L$ acting on a Hilbert space $\mathcal{H}$
with domain $\mathcal{D}(L)$, $z\in\mathbb{C}$ and $\varepsilon>0$, we have
$\norm{\mathcal{R}(\overline{z};L^{\dagger})}=\norm{\mathcal{R}(z;L)}$,
$\sigma(L^{\dagger})=\overline{\sigma(L)}$ and
$\sigma_{\varepsilon}(L^{\dagger})=\overline{\sigma_{\varepsilon}(L)}$.444Regarding
notation, we denote the complex conjugate with a bar and complex conjugate
transpose with an asterisk.

Definition 2.3 hints at the numerical difficulties arising when obtaining
eigenvalues and eigenvectors of non-normal operators. When operating at a
precision of $\varepsilon$ we are unable to distinguish between an
$\varepsilon$-pseudoeigenvalue and a genuine eigenvalue, which, by virtue of
definition 2.2, implies that we cannot obtain the spectrum of the desired
operator but instead that of a perturbed one. Thus, in the absence of
normality, one needs to be especially careful when employing numerical methods
to prevent loss of predictivity. One practical consequence is that even
numerical errors can lead to displacements and therefore it is often necessary
in numerical approximations to compute quasinormal frequencies with (much)
higher precision than the typical double precision floating-point machine
numbers.

Besides pseudospectra, the condition numbers $\\{\kappa_{i}\\}$ are another
useful tool to study the stability of eigenvalues. They quantify the effect of
perturbations of size $\varepsilon$ through the knowledge of the orthogonality
between the eigenvectors of the unperturbed operator and its adjoint (often
referred to as right- and left-eigenvectors, respectively). Heuristically, the
condition numbers exploit the lack of orthogonality between left- and right-
eigenvectors in non-normal operators to look for potential instabilities. This
is formalized in the following definition and theorem Trefethen:2005 :

######  Def. 2.5 (Condition numbers).

Given a closed linear operator $L$ acting on a Hilbert space $\mathcal{H}$ and
its adjoint $L^{\dagger}$, we define the condition number $\kappa_{i}$
associated with the eigenvalue $\lambda_{i}$ of $L$ as:

$\kappa_{i}=\frac{\norm{v_{i}}\norm{u_{i}}}{|\expectationvalue{v_{i},u_{i}}|}\,,$
(7)

where $\expectationvalue{\cdot,\cdot}$ is the inner product associated with
the norm $\norm{\cdot}$, $u_{i}$ is the right-eigenvector satisfying
$Lu_{i}=\lambda_{i}u_{i}$ and $v_{i}$ is the left-eigenvector satisfying
$L^{\dagger}v_{i}=\overline{\lambda_{i}}v_{i}$

######  Thm. 2.4 (Stability and condition numbers).

Given a closed linear operator $L$ acting on a Hilbert space $\mathcal{H}$,
whose spectrum contains a set of eigenvalues $\\{\lambda_{i}\\}$, and a
bounded perturbation $V$ with $\norm{V}=\varepsilon$ we have:

$|\lambda_{i}(\varepsilon)-\lambda_{i}|\leq\varepsilon\kappa_{i}\,,$ (8)

where $\\{\lambda_{i}(\varepsilon)\\}$ are the eigenvalues of the perturbed
operator $L(\varepsilon)=L+V$

###### Proof.

Let us consider a closed operator $L$, its adjoint $L^{\dagger}$ and $u_{i}$,
$v_{i}$ such that:

$Lu_{i}=\lambda_{i}u_{i}\,,\qquad
L^{\dagger}v_{i}=\overline{\lambda_{i}}v_{i}\,.$ (9)

Defining the perturbed operator $L(\varepsilon)=L+V$ with
$\norm{V}=\varepsilon$, the perturbed right-eigenvector
$u_{i}(\varepsilon)=u_{i}+\varepsilon\delta
u_{i}+\mathcal{O}\left(\varepsilon^{2}\right)$ satisfies:

$L(\varepsilon)u_{i}(\varepsilon)=\lambda_{i}(\varepsilon)u_{i}(\varepsilon)\,,$
(10)

with $\lambda_{i}(\varepsilon)$ the corresponding perturbed eigenvalue. Then,
keeping only the leading order in $\varepsilon$ we have:

$\displaystyle\lambda_{i}(\varepsilon)$
$\displaystyle=\frac{\expectationvalue{v_{i},L(\varepsilon)u_{i}(\varepsilon)}}{\expectationvalue{v_{i},u_{i}(\varepsilon)}}=\frac{\expectationvalue{v_{i},Lu_{i}}+\expectationvalue{v_{i},Vu_{i}}+\varepsilon\expectationvalue{v_{i},L\delta
u_{i}}}{\expectationvalue{v_{i},u_{i}}+\varepsilon\expectationvalue{v_{i},\delta
u_{i}}}$
$\displaystyle=\lambda_{i}+\varepsilon\frac{\expectationvalue{L^{\dagger}v_{i},\delta
u_{i}}}{\expectationvalue{v_{i},u_{i}}}-\varepsilon\lambda_{i}\frac{\expectationvalue{v_{i},\delta
u_{i}}}{\expectationvalue{v_{i},u_{i}}}+\frac{\expectationvalue{v_{i},Vu_{i}}}{\expectationvalue{v_{i},u_{i}}}$
$\displaystyle=\lambda_{i}+\varepsilon\lambda_{i}\frac{\expectationvalue{v_{i},\delta
u_{i}}}{\expectationvalue{v_{i},u_{i}}}-\varepsilon\lambda_{i}\frac{\expectationvalue{v_{i},\delta
u_{i}}}{\expectationvalue{v_{i},u_{i}}}+\frac{\expectationvalue{v_{i},Vu_{i}}}{\expectationvalue{v_{i},u_{i}}}=\lambda_{i}+\frac{\expectationvalue{v_{i},Vu_{i}}}{\expectationvalue{v_{i},u_{i}}}\,,$
(11)

from where we recover the result we wanted to prove:

$|\lambda_{i}(\varepsilon)-\lambda_{i}|=\left|\frac{\expectationvalue{v_{i},Vu_{i}}}{\expectationvalue{v_{i},u_{i}}}\right|\leq\left|\frac{\norm{v_{i}}\norm{Vu_{i}}}{\expectationvalue{v_{i},u_{i}}}\right|\leq\left|\frac{\norm{v_{i}}\norm{V}\norm{u_{i}}}{\expectationvalue{v_{i},u_{i}}}\right|=\varepsilon\kappa_{i}\,.$
(12)

∎

For a normal operator all eigenvalues are stable and have condition number
1.555For a normal operator $L$ and its adjoint $L^{\dagger}$ can be
simultaneously diagonalized and thus $v_{i}=u_{i}$ and
$\expectationvalue{v_{i},u_{i}}=\norm{v_{i}}\norm{u_{i}}$. Interestingly
enough, even for non-normal operators some eigenvalues may have condition
number 1 and thus be stable; these are referred to as normal eigenvalues.

Lastly, we note that throughout this section we have assumed $\varepsilon$ to
be small, but we have yet to give a formal description of smallness. As
indicated by definition 2.2, pseudospectrum analysis is, in spirit, equivalent
to a problem of perturbation theory in Quantum Mechanics and thus we decide to
inherit the latter’s definition of smallness:

######  Def. 2.6 (Smallness).

Given a closed linear operator of size $\varepsilon$, we say that
$\varepsilon$ is small if it is negligible when compared to the minimum
distance $d_{\min}$ between disconnected regions of the spectrum,
$\varepsilon/d_{\min}\ll 1$.

###  2.1 Matrix Pseudospectrum

Throughout this work, we will limit ourselves to studying pseudospectra of
matrices arising from the discretization of differential operators.
Consequently, in this subsection, we present some relevant results specialized
to matrices. Note that all previous theorems and definitions still hold since
$N\times N$ matrices are bounded linear operators acting on an $N$ dimensional
Hilbert space comprised of vectors with $N$ components. Again further details
on the following theorems and definitions can be found on Trefethen:2005 .

We begin by constructing a more practical definition of matrix norm using
definition 2.4:

######  Thm. 2.5 (Matrix norm).

Given an $N$-dimensional Hilbert space $\mathcal{V}$ equipped with a norm
$\norm{\cdot}$ and an $N\times N$ matrix $M$; the norm of $M$ is:

$\norm{M}=\sqrt{\rho\left(M^{\dagger}M\right)}\,,$ (13)

with $\rho$ the spectral radius:

$\rho(A)=\underset{\lambda\in\sigma\left(A\right)}{\max}\left\\{|\lambda|\right\\}\,.$
(14)

###### Proof.

We begin by pointing out that, for matrices, all elements in the spectrum are
eigenvalues. Then, by definition 2.4, we have:

$\norm{M}=\underset{u\in\mathcal{V}}{\max}\frac{\norm{Mu}}{\norm{u}}=\underset{u\in\mathcal{V}}{\max}\frac{\sqrt{\expectationvalue{Mu,Mu}}}{\norm{u}}=\underset{u\in\mathcal{V}}{\max}\frac{\sqrt{\expectationvalue{u,M^{\dagger}Mu}}}{\norm{u}}=\underset{u\in\mathcal{V}}{\max}\frac{\sqrt{\expectationvalue{M^{\dagger}Mu,u}}}{\norm{u}}\,,$
(15)

where $\expectationvalue{\cdot,\cdot}$ is the inner product associated with
the norm and $M^{\dagger}$ the adjoint with respect to the aforementioned
inner product.

As $M^{\dagger}M$ is non-negative definite and self-adjoint, we can find an
orthonormal basis $\\{e_{i}\\}$ satisfying
$M^{\dagger}Me_{i}=\lambda_{i}e_{i}$ with $\lambda_{i}$ real positive
eigenvalues. Choosing $\lambda_{n}=\max\\{\lambda_{i}\\}$, we then have:

$\expectationvalue{u,M^{\dagger}Mu}=\sum_{i,j}\overline{u^{i}}u^{j}\expectationvalue{e_{i},M^{\dagger}Me_{j}}=\sum_{i}\overline{u^{i}}u^{i}\lambda_{i}\leq\sum_{i}\lambda_{n}\,\overline{u^{i}}u^{i}=\lambda_{n}\norm{u}^{2}=|\lambda_{n}|\norm{u}^{2}\,,$
(16)

which implies:

$\norm{M}=\underset{u\in\mathcal{V}}{\max}\frac{\sqrt{\expectationvalue{u,M^{\dagger}Mu}}}{\norm{u}}=\sqrt{|\lambda_{n}|}=\sqrt{\underset{\lambda\in\sigma\left(M^{\dagger}M\right)}{\max}\left\\{|\lambda|\right\\}}=\sqrt{\rho\left(M^{\dagger}M\right)}\,.$
(17)

∎

Using this expression of the matrix norm, we can rewrite definition 2.1 in
terms of the minimum generalized singular value. This statement is formalized
in the following theorem:

######  Thm. 2.6.

Given an $N\times N$ matrix $M$ acting on an $N$-dimensional Hilbert space
$\mathcal{V}$ equipped with a norm $\norm{\cdot}$, the
$\varepsilon$-pseudospectrum $\sigma_{\varepsilon}(M)$ is

$\sigma_{\varepsilon}(M)=\\{z\in\mathbb{C}:s_{\min}(M-z)<\varepsilon\\}\,,$
(18)

where $s_{\min}$ is the smallest generalized singular value:

$s_{\min}(A)=\sqrt{\underset{\lambda\in\sigma(A^{\dagger}A)}{\min}{|\lambda|}}\,.$
(19)

###### Proof.

First, we note that given an invertible matrix $A$, we have:

$\underset{\lambda\in\sigma(A^{-1})}{\max}\\{|\lambda|\\}=\left(\underset{\lambda\in\sigma(A)}{\min}\\{|\lambda|\\}\right)^{-1}\,,$
(20)

which follows trivially from the diagonalization of $A$. Then, the norm of the
inverse can be written as follows:

$\norm{A^{-1}}=\left(\underset{\lambda\in\sigma(A)}{\min}\\{|\lambda|\\}\right)^{-1/2}=\left[s_{\min}(A)\right]^{-1}\,.$
(21)

Note that we can consistently define the norm of the inverse of a non-
invertible matrix to be $\infty$, thereby extending this proof to arbitrary
matrices. Now, specializing to the resolvent, we have:

$\norm{\mathcal{R}(z;M)}=\left[s_{\min}(M-z)\right]^{-1}\Longrightarrow\norm{\mathcal{R}(z;M)}<1/\varepsilon\Leftrightarrow
s_{\min}(M-z)>\varepsilon\,,$ (22)

which, applying definition 2.1, implies that the theorem holds. ∎

Having expressed the pseudospectrum in terms of the generalized singular
values, one can now establish the following corollary of theorem 2.2:

######  Cor. 2.1.

Given an $N\times N$ normal matrix $M$ acting on an $N$-dimensional Hilbert
space $\mathcal{V}$ equipped with a norm $\norm{\cdot}$, its
$\varepsilon$-pseudospectrum $\sigma_{\varepsilon}(M)$ is

$\sigma_{\varepsilon}(M)=\sigma(M)+\Delta_{\varepsilon}\,,$ (23)

with $\Delta_{\varepsilon}$ the open disk of radius $\varepsilon$.

###### Proof.

Let us consider a normal matrix $M$ and its adjoint $M^{\dagger}$. As $M$ and
$M^{\dagger}$ commute, they can be simultaneously diagonalized. Furthermore,
theorem 2.3 implies that the eigenvalues of $M$ and $M^{\dagger}$ are complex
conjugate to each other, i.e., $\lambda_{i}^{\dagger}=\overline{\lambda_{i}}$
for $\lambda_{i}\in\sigma(M)$ and
$\lambda_{i}^{\dagger}\in\sigma(M^{\dagger})$. Then, the eigenvalues of
$(M-z)^{\dagger}(M-z)$, $\left\\{\tilde{\lambda}\right\\}$ are given by:

$\left\\{\tilde{\lambda}\right\\}=\sigma\left(M^{\dagger}M-zM^{\dagger}+\overline{z}M+|z|^{2}\right)=\left\\{|\lambda|^{2}-z\overline{\lambda}-\overline{z}\lambda+|z|^{2}\right\\}=\left\\{|z-\lambda|^{2}\right\\}\,,$
(24)

where $\left\\{\lambda\right\\}=\sigma(M)$; and, consequently, $s_{\min}(M-z)$
is:

$s_{\min}(M-z)=|z-\lambda_{i}|\,,$ (25)

for some $\lambda_{i}\in\left\\{\lambda\right\\}$.

Lastly, recalling theorem 2.6, we conclude that the
$\varepsilon$-pseudospectrum is given by:

$\sigma_{\varepsilon}(M-z)=\\{z\in\mathbb{C},\lambda_{i}\in\sigma(M):|z-\lambda_{i}|<\varepsilon\\}=\sigma(M)+\Delta_{\varepsilon}\,.$
(26)

∎

Up to this point, we have been totally generic and have not assumed any
particular norm. However, when working with $N\times N$ matrices it is
particularly convenient to choose the $N\times N$ matrix $\ell^{2}$-norm
($\ell^{2}_{N}$-norm) defined as:

$\norm{u}_{2}=\sqrt{\sum_{a=1}^{N}|u^{a}|^{2}}\,.$ (27)

As we will discuss later, the $\ell^{2}_{N}$ norm allows for significant
optimization of the pseudospectrum algorithm, thus it is convenient to always
employ it. However, we have a preferred norm, the one induced by the
discretization of the original Hilbert space where the differential operator
acts. Then, to exploit the advantages of the $\ell^{2}_{N}$-norm, we need to
establish a connection between both norms.666Relating definitions in a generic
norm to the $\ell^{2}_{N}$-norm is identical to performing a change of basis
in the Hilbert space from a ”coordinate” basis with nontrivial metric to an
”orthogonal” basis with euclidean metric.

######  Thm. 2.7.

Given the $\ell^{2}_{N}$-norm $\norm{\cdot}_{2}$ and generic $G$-norm
$\norm{\cdot}_{G}$ such that:

$\expectationvalue{v,u}_{G}=G_{ij}\overline{v^{i}}u^{j}\,,$ (28)

with $G=F^{*}F$ a symmetric positive definite $N\times N$ matrix.

  * •

The $\varepsilon$-pseudospectrum of a matrix $M$ in the $G$-norm
$\sigma_{\varepsilon}^{G}(M)$ satisfies

$\sigma_{\varepsilon}^{G}(M)=\sigma_{\varepsilon}^{\ell^{2}_{N}}\left(FMF^{-1}\right)\,,$
(29)

with $\sigma_{\varepsilon}^{\ell^{2}_{N}}(M)$ the pseudospectrum in the
$\ell^{2}_{N}$-norm.

  * •

The condition number of the eigenvalue $\lambda_{i}$ of a matrix $M$ in the
G-norm $\kappa_{i}^{G}$ satisfies:

$\kappa_{i}^{G}=\frac{\norm{\tilde{v}_{i}}_{2}\norm{\tilde{u}_{i}}_{2}}{|\expectationvalue{\tilde{v}_{i},\tilde{u}_{i}}_{2}|}\,,$
(30)

where $\tilde{u}_{i}$ and $\tilde{v}_{i}$ fulfill:

$FMF^{-1}\tilde{u}_{i}=\lambda_{i}\tilde{u}_{i}\,,\qquad\left(FMF^{-1}\right)^{*}\tilde{v}_{i}=\overline{\lambda_{i}}\tilde{v}_{i}\,.$
(31)

###### Proof.

First, let us assume an $N$-dimensional Hilbert space $\mathcal{V}$ equipped
with the $G$-norm $\norm{\cdot}_{G}$ that induces the following inner product:

$\expectationvalue{v,u}_{G}=v^{*}Gu\,.$ (32)

For this to represent a well-defined inner product, the matrix $G$ has to be
symmetric and positive definite; thus it can always be expressed as $G=F^{*}F$
using a Cholesky decomposition (for a proof see, for instance, ch. 4 of
Horn:2017 ).777In the language of differential geometry, performing a Cholesky
decomposition is equivalent to choosing triangular vielbeins. For any metric,
this can always be achieved by exploiting the internal O($N$) symmetry of the
”orthogonal” coordinates.

  * •

Considering now an $N\times N$ matrix $A$, we have:

$\displaystyle\norm{A}_{G}^{2}$
$\displaystyle=\underset{u\in\mathcal{V}}{\max}\frac{\expectationvalue{Au,Au}_{G}}{\expectationvalue{u,u}_{G}}=\underset{u\in\mathcal{V}}{\max}\frac{\left(Au\right)^{*}GAu}{u^{*}Gu}=\underset{u\in\mathcal{V}}{\max}\frac{u^{*}A^{*}F^{*}FAu}{u^{*}F^{*}Fu}$
$\displaystyle=\underset{u\in\mathcal{V}}{\max}\frac{\expectationvalue{FAu,FAu}_{2}}{\expectationvalue{Fu,Fu}_{2}}=\underset{\tilde{u}\in\mathcal{V}}{\max}\frac{\expectationvalue{FAF^{-1}\tilde{u},FAF^{-1}\tilde{u}}_{2}}{\expectationvalue{\tilde{u},\tilde{u}}_{2}}=\norm{FAF^{-1}}_{2}^{2}\,,$
(33)

where we have introduced $\tilde{u}=Fu$. This relates the $G$-norm to the
$\ell^{2}_{N}$ norm for any generic matrix. In particular, for the resolvent
of a matrix $M$:

$\displaystyle\norm{\mathcal{R}(z;M)}_{G}$
$\displaystyle=\norm{F\mathcal{R}(z;M)F^{-1}}_{2}=\norm{F(M-z)^{-1}F^{-1}}_{2}$
$\displaystyle=\norm{(FMF^{-1}-z)^{-1}}_{2}=\norm{\mathcal{R}(z;FMF^{-1})}_{2}\,,$
(34)

and using definition 2.1 we trivially recover (29).

  * •

Now we assume an $N\times N$ matrix $M$ with left- and right-eigenvalues
$v_{i}$ and $u_{i}$:

$Mu_{i}=\lambda_{i}u_{i}\,,\qquad
M^{\dagger}v_{i}=\overline{\lambda_{i}}v_{i}\,.$ (35)

The condition number $\kappa_{i}^{G}$ corresponding to the eigenvalue
$\lambda_{i}$ is then given by:

$\kappa_{i}^{G}=\frac{\sqrt{v_{i}^{*}Gv_{i}}\sqrt{u_{i}^{*}Gu_{i}}}{|v_{i}^{*}Gu_{i}|}=\frac{\sqrt{(Fv_{i})^{*}Fv_{i}}\sqrt{(Fu_{i})^{*}Fu_{i}}}{|(Fv_{i})^{*}Fu_{i}|}=\frac{\norm{\tilde{v}_{i}}_{2}\norm{\tilde{u}_{i}}_{2}}{|\expectationvalue{\tilde{v}_{i},\tilde{u}_{i}}_{2}|}\,,$
(36)

where we have introduced $\tilde{u}_{i}=Fu_{i}$ and $\tilde{v}_{i}=Fv_{i}$.

Now, we note that $M^{\dagger}$ in the G norm is defined as:

$M^{\dagger}=\left(GMG^{-1}\right)^{*}=F^{-1}\left(FMF^{-1}\right)^{*}F\,,$
(37)

which trivially follows from the definition of adjoint:

$\expectationvalue{M^{\dagger}g,f}_{G}=\expectationvalue{g,Mf}_{G}\Rightarrow
g^{*}\left(M^{\dagger}\right)^{*}Gf=g^{*}GMf\Rightarrow\left(M^{\dagger}\right)^{*}G=GM\,,$
(38)

with $f$ and $g$ two arbitrary vectors.

Then, we can rewrite (35) in terms of $\tilde{u}_{i}$ and $\tilde{v}_{i}$ as:

$\displaystyle MF^{-1}\tilde{u}_{i}=F^{-1}\lambda_{i}\tilde{u}_{i}\Rightarrow
FMF^{-1}\tilde{u}_{i}=\lambda_{i}\tilde{u}_{i}\,,$ (39) $\displaystyle
F^{-1}\left(FMF^{-1}\right)^{*}FF^{-1}\tilde{v}_{i}=F^{-1}\overline{\lambda_{i}}\tilde{v}_{i}\Rightarrow\left(FMF^{-1}\right)^{*}\tilde{v}_{i}=F^{-1}\overline{\lambda_{i}}\tilde{v}_{i}\,,$
(40)

thus concluding the proof.

∎

Lastly, let us address a nontrivial subtlety of the discretization process.
Definition 2.2 of the $\varepsilon$-pseudospectrum assumes bounded
perturbations, or, in more physical terms, potential perturbations. However,
in a discretized setting all operators become bounded and one cannot
distinguish the discretized version of an originally bounded operator from
that of an originally unbounded one. Then, in many cases it will be more
interesting to probe a particular type of perturbations using definition 2.2,
either choosing a specific form of the perturbation matrix or generating
random matrices satisfying some physically motivated constraints.

####  2.1.1 Norm Dependence in a Simple Example

In order to make the formalism of this section a bit more transparent, we
present the pseudospectrum of a $2\times 2$ matrix $A$

$A=\begin{pmatrix}-1&0\\\ -50&-2\end{pmatrix}\,,$ (41)

with eigenvalues $\omega_{1}=-2$ and $\omega_{2}=-1$.

Furthermore, to make the norm dependence explicit, we present the analysis in
two norms, the $\ell^{2}_{2}$-norm and the $G$-norm:

$\norm{u}_{G}=G_{ij}\overline{u^{i}}u^{j}\,,$ (42)

with $G$

$G=\left(\begin{array}[]{cc}2\cdot 10^{4}&50\\\ 50&1\\\
\end{array}\right)=\left(\begin{array}[]{cc}100\sqrt{2}&0\\\\[6.0pt]
\frac{1}{2\sqrt{2}}&\sqrt{\frac{7}{8}}\end{array}\right)\left(\begin{array}[]{cc}100\sqrt{2}&\frac{1}{2\sqrt{2}}\\\\[6.0pt]
0&\sqrt{\frac{7}{8}}\end{array}\right)=F^{*}F\,.$ (43)

  * •

$\ell^{2}_{2}$-norm:

In this case, the adjoint of $A$ is its complex conjugate transpose

$A^{\dagger}=\begin{pmatrix}-1&-50\\\ 0&-2\end{pmatrix}\,,$ (44)

and $A$ is not normal:

$A^{\dagger}A-AA^{\dagger}=\left(\begin{array}[]{cc}-2500&-50\\\ -50&2500\\\
\end{array}\right)\neq 0\,.$ (45)

Then, one could potentially find instabilities associated with non-normality.
A suspicion confirmed by the condition numbers

$\kappa_{1}=\sqrt{2501}\approx 50\,,\qquad\kappa_{2}=\sqrt{2501}\approx 50\,,$
(46)

which differ greatly from 1, and the pseudospectrum of figure 1(a), which is
characterized by extended $\varepsilon$-pseudospectrum boundaries that
eventually engulf both eigenvalues for small values of $\varepsilon$.

  * •

$G$-norm:

In this case, $A$ is self-adjoint:

$A^{\dagger}=\overline{G^{-1}AG}=\begin{pmatrix}-1&0\\\
-50&-2\end{pmatrix}=A\,,$ (47)

and consequently its spectrum is stable.

Here, the condition numbers are

$\kappa_{1}=1\,,\qquad\kappa_{2}=1\,,$ (48)

and the pseudospectrum is characterized by concentric circles of radius
$\varepsilon$ around the eigenvalues (figure 1(b)), both denoting stability.

As expected, the choice of norm is nontrivial; the notion of stability depends
on the definition of small perturbations. It is also interesting to remark
that figures 1(a) and 1(b) show very common features of the pseudospectrum of
non-normal and normal operators, respectively. In general, the
$\varepsilon$-pseudospectrum of a normal operator is a set disks of radius
$\varepsilon$ centered on the eigenvalues (as indicated by corollary 2.1);
while for a non-normal operator it presents extended regions that may contain
multiple eigenvalues.888Note that for normal operators the
$\varepsilon$-pseudospectrum does engulf multiple eigenvalues if $\varepsilon$
is large. As expected, this behavior appears for $\varepsilon$ larger than
half the distance between the two eigenvalues.

(a) $\ell^{2}_{2}$-norm

(b) $G$-norm

Figure 1: Pseudospectrum of the $A$ matrix (41) in two different norms. In
both figures the red dots correspond to the eigenvalues, the white lines to
boundaries of different $\varepsilon$-pseudospectra and the heat map to the
logarithm of the inverse of the resolvent. In (a) the solid red and yellow
lines represent the boundaries of $\sigma(A)+\Delta_{0.04}$ and
$\sigma_{0.04}(A)$, respectively; note the great discrepancy associated with
the non-normality. On the other hand, in (b) the dashed red line and the solid
yellow line are superimposed and correspond to the boundaries of
$\sigma(A)+\Delta_{0.1}$ and $\sigma_{0.1}(A)$, respectively. The solid black
line in (b) represents the boundary of $\sigma_{0.5}(A)$.

##  3 Quasinormal Modes and Pseudospectra

Quasinormal modes (QNMs) in AdS are solutions to the linearized equations of
motion in a black hole background with well-defined quasinormal frequencies
(QNFs), subject to outgoing boundary conditions on the horizon and
normalizable boundary conditions on the AdS boundary Berti:2009kk .999Imposing
normalizable boundary conditions corresponds to choosing source-less
excitations in the dual QFT or, equivalently, to selecting solutions where the
leading term on the AdS boundary vanishes. Note that in some cases the
subleading mode can also be identified with the source; this corresponds to a
different choice of quantization scheme in the dual QFT Klebanov:1999tb . For
simplicity, in the present work we do not consider that possibility. QNMs are
small excitations of the black hole spacetime that eventually fall into the
horizon and die out, provided the background is stable.101010Here, stability
of the background entails that perturbations decay as opposed to growing
exponentially. From the dual QFT perspective, they correspond to excitations
over the thermal equilibrium that eventually thermalize. The decaying behavior
is reflected in the associated QNFs, which are complex as the system is
inherently non-conservative.

As argued in the previous section, the spectra of non-conservative systems are
potentially unstable. One thus needs to complement the study of QNFs with
pseudospectrum analysis. In order to achieve this, we need to cast the problem
of finding QNFs as an eigenvalue equation and define a physically motivated
norm. In this section we implement these two steps to determine the
pseudospectrum of fluctuations in asymptotic AdS black hole geometries.

Specifically, as a background we consider a static $(\text{d}+1)$-dimensional
AdS planar black hole (black brane) spacetime, whose dual is a d-dimensional
strongly coupled QFT in thermal equilibrium at temperature given by the
Hawking temperature of the black brane. The corresponding metric in Poincaré
coordinates is given by:

$ds^{2}=\frac{r^{2}}{l^{2}}\left(-f(r)dt^{2}+\delta_{ij}dx^{i}dx^{j}\right)+\frac{l^{2}}{r^{2}}\frac{dr^{2}}{f(r)}\,,$
(49)

where $l$ is the AdS radius and $f(r)$ is the blackening factor whose largest
zero we denote by $r_{h}$. In these coordinates, the Hawking temperature $T$
is given by:

$T=\frac{r_{h}^{2}}{4\pi l^{2}}\left.\frac{df(r)}{dr}\right|_{r=r_{h}}\,,$
(50)

and the AdS boundary and the horizon are located at $r=\infty$ and $r=r_{h}$,
respectively.

###  3.1 Construction of the Eigenvalue Problem

Our approach relies on utilizing regular coordinates, as introduced in
Warnick:2013hba ; Ficek:2023imn , to express the linearized equations of
motion as an eigenvalue problem while also implementing the outgoing boundary
conditions as regularity conditions on the horizon.111111Regular coordinates
are the AdS equivalent to the hyperboloidal ones used in asymptotically flat
spacetimes (see Schmidt:1993rcx ; Zenginoglu:2011jz ; Ansorg:2016ztf ;
PanossoMacedo:2018hab  for applications to QNMs and Bizon:2020qnd ;
PanossoMacedo:2023qzp  for a more recent overview).

For the sake of completeness, in this subsection we present the approach in a
general setting, only assuming a black brane background.121212Throughout this
subsection we will not address the normalizable boundary conditions. As we
will see later, these can be easily implemented as Dirichlet or Neumann
boundary conditions with an adequate redefinition of the QNMs. In 3.1.1, we
specialize this discussion to a scalar in an AdS4+1 Schwarzschild black brane
(SAdS4+1) background, hoping to shed more light on the nature of the approach.

In general, we consider linearized equations of motion of the form:

$\left[-\nabla_{M}\nabla^{M}+K^{M}\nabla_{M}+U\right]\Phi=0\,,$ (51)

with $\Phi$ a multi-component field, $K^{M}$ and $U$ functions and
$\nabla_{M}$ the covariant derivative.

Exploiting the symmetries of the black brane, we choose coordinates
$\\{\mathscr{t},\boldsymbol{\mathscr{x}},\mathscr{r}\\}$ such that the linear
operator of equation (51) only depends on $\mathscr{r}$. Furthermore, in order
to eliminate $\boldsymbol{\mathscr{x}}$ derivatives we express $\Phi$ as a
superposition of Fourier modes:

$\Phi(\mathscr{t},\boldsymbol{\mathscr{x}},\mathscr{r})=\int\frac{d^{\text{d}-1}k}{(2\pi)^{\text{d}-1}}\Phi(\mathscr{t},\mathbf{k},\mathscr{r})e^{i\mathbf{k}\boldsymbol{\mathscr{x}}}\,,$
(52)

and solve for $\Phi(\mathscr{t},\mathbf{k},\mathscr{r})$ instead. We can then
rewrite equation (51) as

$g^{\mathscr{t}\mathscr{t}}\partial_{\mathscr{t}}^{2}\Phi(\mathscr{t},\mathbf{k},\mathscr{r})=F_{1}\left[\partial_{\mathscr{r}}^{2},\partial_{\mathscr{r}};\mathbf{k},\mathscr{r}\right]\Phi(\mathscr{t},\mathbf{k},\mathscr{r})+F_{2}\left[\partial_{\mathscr{r}};\mathbf{k},\mathscr{r}\right]\partial_{\mathscr{t}}\Phi(\mathscr{t},\mathbf{k},\mathscr{r})\,,$
(53)

where
$F_{1}\left[\partial_{\mathscr{r}}^{2},\partial_{\mathscr{r}};\mathbf{k},\mathscr{r}\right]$
and $F_{2}\left[\partial_{\mathscr{r}};\mathbf{k},\mathscr{r}\right]$ are
differential operators on $\mathscr{r}$ of second and first order,
respectively. As long as the $\mathscr{t}$-hypersurfaces are not null
($g^{\mathscr{t}\mathscr{t}}\neq 0$), we can introduce the auxiliary variable
$\Psi=\partial_{\mathscr{t}}\Phi$ and rewrite the equations of motion (53)
as:131313For notational convenience, we refer to the $x^{M}=\text{const}$
hypersurface as $x^{M}$-hypersurface.

$\partial_{\mathscr{t}}\begin{pmatrix}\Phi(\mathscr{t},\mathbf{k},\mathscr{r})\\\
\Psi(\mathscr{t},\mathbf{k},\mathscr{r})\end{pmatrix}=\begin{pmatrix}0&1\\\
L_{1}\left[\partial_{\mathscr{r}}^{2},\partial_{\mathscr{r}};\mathbf{k},\mathscr{r}\right]&L_{2}\left[\partial_{\mathscr{r}};\mathbf{k},\mathscr{r}\right]\end{pmatrix}\begin{pmatrix}\Phi(\mathscr{t},\mathbf{k},\mathscr{r})\\\
\Psi(\mathscr{t},\mathbf{k},\mathscr{r})\end{pmatrix}\,,$ (54)

where we have defined
$L_{i}=\left(g^{\mathscr{t}\mathscr{t}}\right)^{-1}F_{i}$.

Recalling that a QNM $\Phi(\mathscr{t},\mathbf{k},\mathscr{r})$ has a well-
defined QNF $\omega$ or, equivalently, that it is an eigenfunction of the
Killing vector associated with stationarity,
$\mathfrak{t}=\partial_{\mathscr{t}}$:

$\mathfrak{t}\begin{pmatrix}\Phi(\mathscr{t},\mathbf{k},\mathscr{r})\\\
\Psi(\mathscr{t},\mathbf{k},\mathscr{r})\end{pmatrix}=\partial_{\mathscr{t}}\begin{pmatrix}\Phi(\mathscr{t},\mathbf{k},\mathscr{r})\\\
\Psi(\mathscr{t},\mathbf{k},\mathscr{r})\end{pmatrix}=-i\omega\begin{pmatrix}\Phi(\mathscr{t},\mathbf{k},\mathscr{r})\\\
\Psi(\mathscr{t},\mathbf{k},\mathscr{r})\end{pmatrix}\,,$ (55)

we can cast expression (54) as a standard eigenvalue problem:

$\omega\begin{pmatrix}\Phi(\mathscr{t},\mathbf{k},\mathscr{r})\\\
\Psi(\mathscr{t},\mathbf{k},\mathscr{r})\end{pmatrix}=i\begin{pmatrix}0&1\\\
L_{1}\left[\partial_{\mathscr{r}}^{2},\partial_{\mathscr{r}};\mathbf{k},\mathscr{r}\right]&L_{2}\left[\partial_{\mathscr{r}};\mathbf{k},\mathscr{r}\right]\end{pmatrix}\begin{pmatrix}\Phi(\mathscr{t},\mathbf{k},\mathscr{r})\\\
\Psi(\mathscr{t},\mathbf{k},\mathscr{r})\end{pmatrix}\,.$ (56)

Note that to write the equations of motion as an eigenvalue problem it is
fundamental that the $\mathscr{t}$-hypersurfaces are spacelike as opposed to
null. Otherwise, we would obtain a generalized eigenvalue equation of the
form:

$\omega
F_{2}\left[\partial_{\mathscr{r}};\mathbf{k},\mathscr{r}\right]\Phi(\mathscr{t},\mathbf{k},\mathscr{r})=-iF_{1}\left[\partial_{\mathscr{r}}^{2},\partial_{\mathscr{r}};\mathbf{k},\mathscr{r}\right]\Phi(\mathscr{t},\mathbf{k},\mathscr{r})\,,$
(57)

instead of the standard eigenvalue problem discussed in section
2.141414Technically, one can extend the formalism presented in section 2 to
also hold for generalized eigenvalue problems of the form $\left(L-\lambda
A\right)u=0$. However, the adequate notion of $\varepsilon$-pseudospectrum
becomes unclear as one has to choose whether to perturb $L$, $A$ or both (see
Chapter 45 of Trefethen:2005  and references therein for an extensive
discussion).

Regarding the boundary conditions, we seek to choose
$\\{\mathscr{t},\boldsymbol{\mathscr{x}},\mathscr{r}\\}$ such that demanding
regularity on the horizon is equivalent to imposing outgoing boundary
conditions. This is typically achieved working with infalling Eddington-
Finkelstein (IEF) coordinates $\\{u,\tilde{\mathbf{x}},\tilde{r}\\}$, which,
in terms of the Poincaré coordinates $\\{t,\mathbf{x},r\\}$, are given by:

$u=t+\int\frac{dr}{f(r)}\left(\frac{l}{r}\right)^{2}\,,\qquad\tilde{\mathbf{x}}=\mathbf{x}\,,\qquad\tilde{r}=r\,.$
(58)

However, writing the metric (49) in IEF coordinates:

$ds^{2}=\frac{\tilde{r}^{2}}{l^{2}}\left(-f(\tilde{r})du^{2}+\delta_{ij}dx^{i}dx^{j}\right)+2dud\tilde{r}\,,$
(59)

we conclude that the $u$-hypersurfaces are null. Thus IEF coordinates do not
result in an standard eigenvalue problem for the QNMs.

Nonetheless, we can define regular coordinates that resemble the IEF ones near
the horizon while having spacelike $\mathscr{t}$-hypersurfaces outside the
horizon. More concretely, we construct regular coordinates such that the
spacelike $\mathscr{t}$-hypersurfaces match the null IEF $u$-hypersurfaces
exactly on the horizon Warnick:2013hba . A particularly simple choice is given
by:

$\mathscr{t}=u-\frac{l^{2}}{r_{h}}\left(1-\frac{r_{h}}{\tilde{r}}\right)=t-\frac{l^{2}}{r_{h}}\left(1-\frac{r_{h}}{r}\right)+\int\frac{dr}{f(r)}\left(\frac{l}{r}\right)^{2}\,,\qquad\boldsymbol{\mathscr{x}}=\tilde{\mathbf{x}}=\mathbf{x}\,,\qquad\mathscr{r}=\tilde{r}=r\,,$
(60)

in terms of which the metric (49) is:

$ds^{2}=\frac{\mathscr{r}^{2}}{l^{2}}\left(-f(\mathscr{r})d\mathscr{t}^{2}+\delta_{ij}d\mathscr{x}^{i}d\mathscr{x}^{j}\right)+2\left(1-f(\mathscr{r})\right)d\mathscr{t}d\mathscr{r}+\frac{l^{2}}{\mathscr{r}^{2}}\left(2-f(\mathscr{r})\right)d\mathscr{r}^{2}\,.$
(61)

In practice, for most numerical applications it is more convenient to
introduce a dimensionless compactified radial coordinate $\rho$:

$\rho=1-\frac{r_{h}}{\mathscr{r}}\,,$ (62)

such that the horizon is mapped to $\rho=0$ and the AdS boundary to $\rho=1$.
In terms of these compactified regular coordinates
$\\{\mathscr{t},\boldsymbol{\mathscr{x}},\rho\\}$, the black brane metric (49)
takes the form:

$ds^{2}=\frac{l^{2}}{z_{h}^{2}(1-\rho)^{2}}\left(-\mathscr{f}(\rho)d\mathscr{t}^{2}+\delta_{ij}d\mathscr{x}^{i}d\mathscr{x}^{j}+2\left(1-\mathscr{f}(\rho)\right)z_{h}d\mathscr{t}d\rho+\left(2-\mathscr{f}(\rho)\right)z_{h}^{2}d\rho^{2}\right)\,.$
(63)

where $\mathscr{f}(\rho)=f\left(r\left(\rho\right)\right)$ and
$z_{h}=l^{2}/r_{h}$.

It is worth stressing that the components of the Killing vector
$\mathfrak{t}=\partial_{\mathscr{t}}=\partial_{t}$ are invariant under the
change of coordinates from regular to Poincaré coordinates. Consequently, the
QNFs defined as the eigenvalues of $\partial_{t}$ and $\partial_{\mathscr{t}}$
match. In gauge/gravity duality this is important as the holographic
dictionary identifies the poles of the retarded propagators with the
eigenvalues of $\partial_{t}$.151515Here, by the holographic dictionary we
refer to the standard one formulated in Poincaré coordinates. Note that one
can always formulate a new holographic dictionary in a different coordinate
set and identify the QNFs with the poles of the retarded propagators.

Before moving on to a concrete example, we summarize the core principles
behind the regular coordinates: Given an AdS${}_{\text{d}+1}$ black brane
spacetime whose metric in Poincaré coordinates is given by equation (49) a
regular coordinate set
$\\{\mathscr{t},\boldsymbol{\mathscr{x}},\mathscr{r}\\}$ is one that fulfils
the following three requirements: the metric components are functions only of
$\mathscr{r}$, $\mathscr{t}$-hypersurfaces are spacelike outside the horizon,
and those $\mathscr{t}$-hypersurfaces match the IEF $u$-hypersurfaces on the
horizon. Imposing outgoing boundary conditions in these coordinates
corresponds to demanding regularity on the horizon. Moreover, equations of
motion for the QNMs of the form (51) can be written as a standard eigenvalue
problem.

####  3.1.1 Example: Real Scalar in $\text{SAdS}_{4+1}$

To better clarify the role played by the regular coordinates in the
construction of the eigenvalue problem, here we specialize the discussion to a
real scalar field $\phi$ with action

$\mathcal{S}[\phi]=-\frac{1}{2}\int
d^{4+1}x\sqrt{-g}\,\,\left(\partial_{M}\phi\partial^{M}\phi+m^{2}\phi^{2}\right)\,,$
(64)

in a SAdS4+1 background. The corresponding linearized equations of motion are
given by:

$\left[-\nabla_{M}\nabla^{M}+m^{2}\right]\phi=0\,.$ (65)

For better emphasis on the key aspects of our approach, we present this
discussion in Poincaré, IEF and compactified regular coordinates.

Poincaré Coordinates: the SAdS4+1 metric in Poincaré coordinates is given by
(49) with $f(r)=1-\left(\frac{r_{h}}{r}\right)^{4}$ and $\\{i,j\\}=1,2,3$. In
terms of the Fourier modes $\phi(t,\mathbf{k},r)$, the linearized equation of
motion (65) for a QNM trivially produces the following eigenvalue problem:

$\omega^{2}\,\phi(t,\mathbf{k},r)=f(r)\left[\mathbf{k}^{2}+\frac{m^{2}r^{2}}{l^{2}}-\frac{\left(r^{5}f(r)\right)^{\prime}}{rl^{4}}\partial_{r}-\frac{r^{4}f(r)}{l^{4}}\partial_{r}^{2}\right]\phi(t,\mathbf{k},r)\,,$
(66)

where we have used that the QNM is an eigenfunction of $\partial_{t}$ with QNF
$\omega$.

In order to impose outgoing boundary conditions, we first analyze the behavior
near the horizon $r=r_{h}$. Solving the equations of motion with a Frobenius
series, we find:

$\phi(t,\mathbf{k},r)=e^{-i\omega t}\left\\{c_{1}\exp\left[-i\frac{\omega
l^{2}}{4r_{h}}\log(r-r_{h})\right]\left(1+...\right)+c_{2}\exp\left[i\frac{\omega
l^{2}}{4r_{h}}\log(r-r_{h})\right]\left(1+...\right)\right\\}\,.$ (67)

As the outgoing solution corresponds to fixing $c_{2}=0$, imposing outgoing
boundary conditions correspond to fixing a near-horizon behavior:

$\phi(t,\mathbf{k},r\rightarrow r_{h})\approx
c_{1}\exp\left[-i\omega\left(t+\frac{l^{2}}{4r_{h}}\log(r-r_{h})\right)\right]\,.$
(68)

This choice for the outgoing solution matches our physical interpretation of a
mode that falls into the horizon. Defining the tortoise coordinate $r_{*}$

$r_{*}=\frac{l^{2}}{2r_{h}}\left[\arctan\left(\frac{r}{r_{h}}\right)-\operatorname{arctanh}\left(\frac{r_{h}}{r}\right)\right]=\frac{l^{2}}{4r_{h}}\log(r-r_{h})+\mathcal{O}\left((r-r_{h})^{0}\right)\,,$
(69)

the near-horizon behavior can be written as two modes, one falling towards the
brane and one exiting it:

$\phi(t,\mathbf{k},r\rightarrow r_{h})\approx
c_{1}\exp\left[-i\omega\left(t+r_{*}\right)\right]+c_{2}\exp\left[-i\omega\left(t-r_{*}\right)\right]\,,$
(70)

and choosing a QNM that falls into the horizon corresponds to discarding the
ingoing mode, i.e., taking $c_{2}=0$.

IEF Coordinates: the SAdS4+1 metric in IEF coordinates is given by (59) with
$f(\tilde{r})=1-\left(\frac{r_{h}}{\tilde{r}}\right)^{4}$ and
$\\{i,j\\}=1,2,3$. The linearized equation of motion for a QNM (65) reads:

$\omega\left[\frac{3}{\tilde{r}}+2\partial_{\tilde{r}}\right]\phi(u,\mathbf{k},\tilde{r})=i\left[\frac{\mathbf{k}^{2}l^{2}}{\tilde{r}^{2}}+m^{2}-\frac{\left(\tilde{r}^{5}f(\tilde{r})\right)^{\prime}}{\tilde{r}^{3}l^{2}}\partial_{\tilde{r}}-\frac{\tilde{r}^{2}f(\tilde{r})}{l^{2}}\partial_{\tilde{r}}^{2}\right]\phi(u,\mathbf{k},\tilde{r})\,,$
(71)

and it cannot be cast as a standard eigenvalue problem due to the derivative
term. The near-horizon solution is given by:

$\phi(u,\mathbf{k},\tilde{r})=e^{-i\omega
u}\left\\{c_{1}\left(1+...\right)+c_{2}\exp\left[i\frac{\omega
l^{2}}{2r_{h}}\log(\tilde{r}-r_{h})\right]\left(1+...\right)\right\\}\,,$ (72)

with the outgoing solution corresponding to fixing $c_{2}=0$.161616This
follows from the Poincaré result noting that under the coordinate
transformation (58) we have:  $e^{-i\omega
u}=\exp\left[-i\omega\left(t+\frac{l^{2}}{4r_{h}}\log\left(r-r_{h}\right)\right)\right]+\mathcal{O}\left((r-r_{h})^{0}\right)\,.$
Consequently, imposing outgoing boundary conditions corresponds to demanding
regularity on the horizon:

$\phi(t,\mathbf{k},\tilde{r}\rightarrow r_{h})\approx c_{1}e^{-i\omega u}\,.$
(73)

Compactified Regular Coordinates: the SAdS4+1 metric in compactified regular
coordinates $\\{\mathscr{t},\boldsymbol{\mathscr{x}},\rho\\}$ is given by (63)
with $\mathscr{f}(\rho)=1-\left(1-\rho\right)^{4}$ and $\\{i,j\\}=1,2,3$.
Then, the linearized equation of motion for a QNM (65) is:

$\displaystyle
z_{h}^{2}\omega^{2}\left[\mathscr{f}(\rho)-2\right]\phi(\mathscr{t},\mathbf{k},\mathscr{r})=iz_{h}\omega\left[(1-\rho)^{3}\left(\frac{\mathscr{f}(\rho)-1}{(1-\rho)^{3}}\right)^{\prime}+2\left(\mathscr{f}(\rho)-1\right)\partial_{\rho}\right]\phi(\mathscr{t},\mathbf{k},\mathscr{r})$
$\displaystyle-\left[\frac{m^{2}l^{2}}{(1-\rho)^{2}}+z_{h}^{2}\mathbf{k}^{2}-(1-\rho)^{3}\left(\frac{\mathscr{f}(\rho)}{(1-\rho)^{3}}\right)^{\prime}\partial_{\rho}-\mathscr{f}(\rho)\partial_{\rho}^{2}\right]\phi(\mathscr{t},\mathbf{k},\mathscr{r})\,,$
(74)

and, introducing the dimensionless variables

$\mathfrak{w}=z_{h}\omega=\frac{\omega}{\pi
T}\,,\quad\boldsymbol{\mathfrak{q}}=z_{h}\mathbf{k}=\frac{\mathbf{k}}{\pi
T}\,,$ (75)

and the auxiliary field
$\psi(\mathscr{t},\mathbf{k},\mathscr{r})=z_{h}\partial_{\mathscr{t}}\phi(\mathscr{t},\mathbf{k},\mathscr{r})$,
equation (3.1.1) can be written as an eigenvalue problem 171717One could cast
the eigenvalue problem without introducing the dimensionless variables
$\mathfrak{w}$ and $\boldsymbol{\mathfrak{q}}$. However, it is convenient to
use them as they allows us to remove all $z_{h}$ dependence and thus reduce
the parameter space.

$\mathfrak{w}\begin{pmatrix}\phi(\mathscr{t},\mathbf{k},\rho)\\\
\psi(\mathscr{t},\mathbf{k},\rho)\end{pmatrix}=i\begin{pmatrix}0&1\\\
L_{1}\left[\partial_{\mathscr{\rho}}^{2},\partial_{\mathscr{\rho}};\boldsymbol{\mathfrak{q}},\rho\right]&L_{2}\left[\partial_{\rho};\boldsymbol{\mathfrak{q}},\rho\right]\end{pmatrix}\begin{pmatrix}\phi(\mathscr{t},\mathbf{k},\rho)\\\
\psi(\mathscr{t},\mathbf{k},\rho)\end{pmatrix}\,,$ (76)

where the differential operators $L_{1}$ and $L_{2}$ take the form:

$\displaystyle
L_{1}\left[\partial_{\rho}^{2},\partial_{\rho};\boldsymbol{\mathfrak{q}},\rho\right]$
$\displaystyle=\left[\mathscr{f}(\rho)-2\right]^{-1}\left[\frac{m^{2}l^{2}}{(1-\rho)^{2}}+\boldsymbol{\mathfrak{q}}^{2}-(1-\rho)^{3}\left(\frac{\mathscr{f}(\rho)}{(1-\rho)^{3}}\right)^{\prime}\partial_{\rho}-\mathscr{f}(\rho)\partial_{\rho}^{2}\right]\,,$
(77) $\displaystyle
L_{2}\left[\partial_{\rho};\boldsymbol{\mathfrak{q}},\rho\right]$
$\displaystyle=\left[\mathscr{f}(\rho)-2\right]^{-1}\left[(1-\rho)^{3}\left(\frac{\mathscr{f}(\rho)-1}{(1-\rho)^{3}}\right)^{\prime}+2\left(\mathscr{f}(\rho)-1\right)\partial_{\rho}\right]\,.$
(78)

With regard to the near-horizon behavior, we have:

$\phi(\mathscr{t},\mathbf{k},\rho)=e^{-i\omega\mathscr{t}}\left\\{c_{1}\left(1+...\right)+c_{2}\exp\left[i\frac{z_{h}\omega}{2}\log(\rho)\right]\left(1+...\right)\right\\}\,,$
(79)

with the outgoing solution corresponding to fixing $c_{2}=0$.181818Once again,
this follows from the Poincaré result as under the coordinate transformation
(60) we have:
$e^{-i\omega\mathscr{t}}=\exp\left[-i\omega\left(t+\frac{l^{2}}{4r_{h}}\log\left(r-r_{h}\right)\right)\right]+\mathcal{O}\left((r-r_{h})^{0}\right)\,.$
Then, imposing outgoing boundary conditions corresponds to demanding
regularity on the horizon:

$\phi(\mathscr{t},\mathbf{k},\rho\rightarrow 0)\approx
c_{1}e^{-i\omega\mathscr{t}}\,.$ (80)

Therefore, it is clear that regular coordinates constitute a sort of middle
ground between IEF and Poincaré coordinates. They exploit the regular behavior
on the horizon of the IEF coordinates while maintaining the spacelike
character of the constant-time hypersurfaces found in Poincaré coordinates. In
fact, with our particular choice of regular coordinates, the
$\mathscr{t}$-hypersurfaces can be thought of as interpolating between $u$\-
and $t$-hypersurfaces (see figure 2).

Figure 2: Penrose Diagram of SAdS4+1. The AdS boundary is denoted by
$\mathcal{J}$, $\mathcal{H}^{+}$ ($\mathcal{H}^{-}$) represents the future
(past) horizon and $i^{+}$ ($i^{-}$) denotes the future (past) time-like
infinity. The dashed blue lines correspond to regular
$\mathscr{t}$-hypersurfaces and the red and cyan dashed lines represent
Poincaré $t$ and IEF $u$-hypersurfaces, respectively.

###  3.2 Choice of Norm

As discussed in section 2, pseudospectra are norm dependent. Moreover, we
concluded that the spectrum of an operator may be unstable in one norm but not
in a different one (as illustrated in figure 1), which shows the nontrivial
role of the norm. Ideally, we want to find a physically motivated norm such
that our physical notion of smallness matches the mathematical one.

Implicitly, we have been assuming that the QNMs are fluctuations small enough
not to affect the background, or equivalently, that their backreaction is
negligible.191919This is the underlying assumption behind studying their
linearized equations of motion on a fixed background. Consequently, the
physically motivated notion of size for the QNMs is their contribution to the
energy-momentum tensor, which is directly related to their backreaction. Thus,
the adequate norm is the energy norm $\norm{\cdot}_{E}$, introduced in
Jaramillo:2020tuu ; Gasperin:2021kfv  in the case of quasinormal modes and
generally advocated for in Trefethen:2005 202020For an application of
pseudosepctra with energy norm in hydrodynamics see reddy1993pseudospectra .
It defines the norm of a QNM $\Phi$ as its energy in a constant-time
hypersurface:

$\norm{\Phi}_{E}^{2}=-\int_{\mathscr{t}=\text{const}}\star
J\left[\Phi\right]\,,$ (81)

where $\star J\left[\Phi\right]$ is the Hodge dual of the current 1-form

$J\left[\Phi\right]=\mathfrak{t}^{M}T_{MN}\left[\Phi\right]dx^{N}\,,$ (82)

with $T_{MN}\left[\Phi\right]$ the QNM’s leading order contribution to the
energy-momentum tensor.

Therefore, the energy norm is generically given by:212121Note that the leading
order contribution of the fluctuations to the energy momentum tensor is
quadratic. This occurs because the linear contribution vanishes as it is
proportional to the equations of motion of the background.

$\displaystyle\norm{\Phi}_{E}^{2}=\int_{\mathscr{t}=\text{const}}d^{\text{d}}x\,\,\Phi^{*}(\mathscr{t},\boldsymbol{\mathscr{x}},\mathscr{r})$
$\displaystyle\left[\overleftarrow{\partial}_{M}H_{1}(\mathscr{r})\overrightarrow{\partial}^{M}+H_{M}(\mathscr{r})\overrightarrow{\partial}^{M}+\right.$
$\displaystyle\left.+\overleftarrow{\partial}^{M}H_{M}^{*}(\mathscr{r})+H_{2}(\mathscr{r})\right]\Phi(\mathscr{t},\boldsymbol{\mathscr{x}},\mathscr{r})\,,$
(83)

where $H_{1}(\mathscr{r}),H_{2}(\mathscr{r}),H_{M}(\mathscr{r})$ are functions
respecting the black brane symmetries. Expression (3.2) can be further
simplified decomposing into Fourier modes and introducing the auxiliary field
$\Psi(\mathscr{t},\mathbf{k},\mathscr{r})=\partial_{\mathscr{t}}\Phi(\mathscr{t},\mathbf{k},\mathscr{r})$,
obtaining:

$\norm{\Phi}_{E}^{2}=\int_{\mathscr{t}=\text{const}}\frac{d^{\text{d}-1}{k}}{(2\pi)^{\text{d}-1}}\,d\mathscr{r}\,\,\begin{pmatrix}\Phi^{*}(\mathscr{t},\mathbf{k},\mathscr{r})&\Psi^{*}(\mathscr{t},\mathbf{k},\mathscr{r})\end{pmatrix}\mathcal{G}\left[\overleftarrow{\partial}_{\mathscr{r}},\overrightarrow{\partial}_{\mathscr{r}};\mathbf{k},\mathscr{r}\right]\begin{pmatrix}\Phi(\mathscr{t},\mathbf{k},\mathscr{r})\\\\[6.0pt]
\Psi(\mathscr{t},\mathbf{k},\mathscr{r})\end{pmatrix}\,,$ (84)

where we have defined

$\mathcal{G}\left[\overleftarrow{\partial}_{\mathscr{r}},\overrightarrow{\partial}_{\mathscr{r}};\mathbf{k},\mathscr{r}\right]=\begin{pmatrix}\mathcal{G}_{11}\left[\overleftarrow{\partial}_{\mathscr{r}},\overrightarrow{\partial}_{\mathscr{r}};\mathbf{k},\mathscr{r}\right]&\mathcal{G}_{12}\left[\overleftarrow{\partial}_{\mathscr{r}};\mathbf{k},\mathscr{r}\right]\\\\[6.0pt]
\mathcal{G}_{21}\left[\overrightarrow{\partial}_{\mathscr{r}};\mathbf{k},\mathscr{r}\right]&\mathcal{G}_{11}[\mathscr{r}]\end{pmatrix}\,,$
(85)

with $\mathcal{G}_{ab}$ differential operators on $\mathscr{r}$.

As discussed in the previous section, we present the eigenvalue problem for
the Fourier modes, i.e., we work out the spectral stability of the operator
acting on the subspace of QNMs with well-defined momentum $\mathbf{k}$.
Consequently, as we do not consider mixing between subspaces, we can drop the
integral over $\mathbf{k}$ and introduce a constant prefactor $C$ instead.
Furthermore, as we are mainly interested in the operator norm (which is
independent of constant prefactors by virtue of definition 2.4), we can set
$C=1$ without losing any generality. With all these considerations, the
relevant inner product can be written as:

$\expectationvalue{\Phi_{1},\Phi_{2}}_{E}=\int_{\mathscr{t}=\text{const}}d\mathscr{r}\,\,\begin{pmatrix}\Phi_{1}^{*}(\mathscr{t},\mathbf{k},\mathscr{r})&\Psi_{1}^{*}(\mathscr{t},\mathbf{k},\mathscr{r})\end{pmatrix}\mathcal{G}\left[\overleftarrow{\partial}_{\mathscr{r}},\overrightarrow{\partial}_{\mathscr{r}};\mathbf{k},\mathscr{r}\right]\begin{pmatrix}\Phi_{2}(\mathscr{t},\mathbf{k},\mathscr{r})\\\\[6.0pt]
\Psi_{2}(\mathscr{t},\mathbf{k},\mathscr{r})\end{pmatrix}\,.$ (86)

Lastly, note that one needs to explicitly check that the energy norm is
positive definite; otherwise, it would not be a well-defined norm. Whether
this condition is satisfied depends on the particular choice of Hilbert space.
Consequently, one must select a function space that enforces the desired
asymptotic behavior for the QNMs while also ensuring that all functions
residing in the aforementioned space have positive energy. In section 4 we
discuss this topic in detail for our particular model.

####  3.2.1 The Nature of the Energy Norm

In this subsection we address the coordinate dependence of the energy norm.
Definition (81) depends on the selection of constant-time hypersurfaces, and
consequently, one could wonder if regular $\mathscr{t}$-hypersurfaces actually
represent a well-motivated choice.

A physically well-defined energy should decay over time, as the fluctuations
fall into the black brane. In order to showcase that the energy norm (81)
satisfies this property, we consider the simple model introduced in subsection
3.1.1. In particular, to stress the importance of choosing regular
$\mathscr{t}$-hypersurfaces, we compute the time derivative of energy for the
real scalar field in Poincaré, IEF and compactified regular coordinates.
Completely analogous considerations hold in the case of the transverse gauge
field.

Poincaré Coordinates: the energy of a QNM in Poincaré coordinates is given by

$E[\phi]=\frac{1}{2l^{5}}\int_{t=\text{const}}d^{3}\mathbf{x}\,dr\,\,r^{3}\,\left[m^{2}l^{2}\phi^{2}+\frac{l^{4}}{r^{2}}\left(\partial_{\mathbf{x}}\phi\right)^{2}+r^{2}f(r)\left(\partial_{r}\phi\right)^{2}+\frac{l^{4}}{r^{2}f(r)}\left(\partial_{t}\phi\right)^{2}\right]\,,$
(87)

where, for notational convenience, we omit the dependence of the field $\phi$
on the coordinates. Using the conservation equation of the current (82):

$d\star J=0\,,$ (88)

we can conclude that the time derivative of the energy vanishes

$\frac{d}{dt}E[\phi]=\frac{r_{h}^{5}}{l^{5}}f(r_{h})\int_{r=r_{h}}d^{3}\mathbf{x}\,\,\partial_{t}\phi\partial_{r}\phi=0\,$
(89)

as it is proportional to the blackening factor evaluated at the
horizon.222222In general, we would also have a term arising from the AdS
boundary. However, in all the coordinates studied in this section, this term
vanishes due to the normalizable boundary conditions. Physically this behavior
was to be expected, as the AdS boundary behaves as a wall for the QNMs. A
quasinormal mode however is an exponentially decaying solution and one would
expect the energy to decay as well. The resolution to this puzzle is that in
Poincaré coordinates a QNM is actually singular on the horizon and thus the
energy is not well defined. We also note that the term proportional to the
time derivatives in the energy is singular on the horizon. For these reasons
the Poincaré coordinates are not suited to define a norm on which the
pseudospectrum analysis could be based.

IEF Coordinates: in IEF coordinates, the energy of a QNM is

$E[\phi]=\frac{1}{2l^{5}}\int_{u=\text{const}}d^{3}\tilde{\mathbb{x}}\,d\tilde{r}\,\,\tilde{r}^{3}\,\left[m^{2}l^{2}\phi^{2}+\frac{l^{4}}{\tilde{r}^{2}}\left(\partial_{\tilde{\mathbb{x}}}\phi\right)^{2}+\tilde{r}^{2}f(\tilde{r})\left(\partial_{\tilde{r}}\phi\right)^{2}\right]\,,$
(90)

and its time derivative is:

$\frac{d}{du}E[\phi]=-\frac{r_{h}^{5}}{l^{5}}\int_{\tilde{r}=r_{h}}d^{3}\tilde{\mathbf{x}}\,\,(\partial_{u}\phi)^{2}<0\,.$
(91)

Thus, in IEF the decay rate of the energy is a function of the time derivative
of the QNMs. Indeed, this was the picture we expected, proving that the IEF
coordinates are well suited to define a physically motivated energy. However,
IEF coordinates do not generate a standard eigenvalue problem; and thus we do
not use them.

Compactified Regular Coordinates: in compactified regular coordinates, the
energy of a QNM associated with a real scalar with action (64) is

$\displaystyle
E[\phi]=\frac{l^{3}}{2z_{h}^{3}}\int_{\mathscr{t}=\text{const}}d^{3}\boldsymbol{\mathscr{x}}\,\frac{d\rho}{(1-\rho)^{3}}\,\,\Biggl{[}$
$\displaystyle\left(\frac{m^{2}l^{2}}{z_{h}(1-\rho)^{2}}\phi^{2}+z_{h}\left(\partial_{\boldsymbol{\mathscr{x}}}\phi\right)^{2}\right)+\frac{\mathscr{f}(\rho)}{z_{h}}\left(\partial_{\rho}\phi\right)^{2}$
$\displaystyle+z_{h}\left(2-\mathscr{f}(\rho)\right)\left(\partial_{\mathscr{t}}\phi\right)^{2}\Biggr{]}\,,$
(92)

and its time derivative is

$\frac{d}{d\mathscr{t}}E[\phi]=-\frac{l^{3}}{z_{h}^{3}}\int_{\rho=0}d^{3}\boldsymbol{\mathscr{x}}\,\,(\partial_{\mathscr{t}}\phi)^{2}<0\,,$
(93)

which, as IEF coordinates, does not vanish. Since the compactified regular
coordinates also allow to compute the QNMs from a standard eigenvalue problem
it is the best suited coordinate system for the pseudospectrum analysis.

####  3.2.2 Time Independence of the Operator Norm

Since we have found that the energy norm in regular coordinates is time
dependent, one could ask what $\mathscr{t}$-hypersurface we should choose to
define it. This is however not a relevant question for the pseudospectrum
analysis. Pseudospectra depend on the norm (3.2) only trough the operator norm
(6), which, by virtue of its definition, is time independent.

For functions belonging to the relevant function space defined by the boundary
conditions, the energy norm indeed defines a positive definite inner product
and thus promotes the function space to a Hilbert space. Functions change over
time, moving around in the Hilbert space. However, the Hilbert space itself
does not change, as it contains all possible function configurations.
Consequently, noting that the operator norm is defined as the maximum value of
the norm of the operator acting on a normalized member of the aforementioned
Hilbert space, we can trivially conclude that the operator norm is time
independent. Thus, we need not worry about choosing one particular
$\mathscr{t}$-hypersurface or another.

###  3.3 Spectral and Background Instability

Before moving on, this is a good point to discuss in detail the holographic
perspective on the (in)stability probed through pseudospectra. From the
gravitational side, pseudospectra probe how much the QNFs change if the
background is slightly modified in some way (e.g., by a change of the geometry
and/or the background value of the fields). Consequently, in the dual quantum
field theory, pseudospectra help us estimate how much the poles of the
retarded Green’s functions might change if the theory is slightly perturbed.
In both cases, these perturbations should be understood as perturbations to
the Lagrangian, leading to a change in the spectrum of excitations. To stress
this point, in this section we will explicitly refer to the (in)stability
associated with these perturbations as spectral (in)stability.

Spectral instability is in nature very different from the usual background
instability which occurs when a linear fluctuation grows indefinitely. The
latter signifies an instability of the background or equivalently the ground
state of the dual field theory. In contrast, spectral instability denotes
instability of the spectrum of the model under slight deformations.
Heuristically, spectral instability implies that a small perturbation to the
Lagrangian can have a large effect on the energy spectrum.

Hence, spectral instability implies that the frequencies and damping rates
computed through an idealized model will differ from those observed in nature.
In the context of holography this suggests that the transport properties of a
quantum many-body system are likely to differ significantly from those
predicted by the holographic model.

It is worth mentioning that in some cases spectral instability may be
connected to background instability. A model with no QNFs in the upper half of
the complex plane, i.e., with a stable background, could exhibit enough
spectral instability so that a small perturbation could drive a QNF into the
upper half of the complex plane. This would imply that when modelling a
quantum many-body system through the aforementioned model, one may not even be
able to accurately predict the equilibrium state. The small error made in the
idealization of the system could lead to very different conclusions.232323An
example of this is the pseudospectrum of the Orr-Sommerfeld operator for
Couette flow in hydrodynamics reddy1993pseudospectra ;
trefethen1993hydrodynamic . In holography this might become relevant in models
with second order phase transitions such as the holographic superconductor
Amado:2009ts ; Gubser:2008px ; Hartnoll:2008vx . There, spectral instability
might be of particular importance in the presence of superflow and for the
validity of the Landau criterion Amado:2013aea ; Gouteraux:2022kpo . It is
important to note, however, that we have not found evidence of this kind of
instability for the scalar and transverse gauge fields studied in this paper.

##  4 Holographic Model

After presenting the general formalism, we now proceed to construct the
specific fluctuations whose stability we are going to study. We seek to probe
the spectral stability of a real scalar $\phi$ and of the transverse
components of a $U(1)$ gauge field $A_{M}$ in a SAdS4+1 background with
compactified regular coordinates
$\\{\mathscr{t},\boldsymbol{\mathscr{x}},\rho\\}$ and metric (63). Their
respective actions are given by:

$\mathcal{S}[\phi]=-\frac{1}{2}\int
d^{4+1}x\,\,\sqrt{-g}\left(\partial_{M}\phi\partial^{M}\phi+m^{2}\phi^{2}\right)\,,$
(94) $\mathcal{S}[A]=-\frac{1}{4}\int
d^{4+1}x\,\,\sqrt{-g}\left(F_{MN}F^{MN}\right)\,,$ (95)

where $F_{MN}=\partial_{M}A_{N}-\partial_{N}A_{M}$ is the field strength
tensor and $m$ is the scalar mass, which we assume to be above the
Breitenlohner-Freedman (BF) bound ($m^{2}l^{2}>-4$) Breitenlohner:1982bm .
This latter requirement ensures that the QNM is normalizable and that the dual
QFT is unitary Aharony:1999ti .

###  4.1 Real Scalar

The eigenvalue problem for a real scalar field with action (94) has already
been presented in equations (76)-(78) of example 3.1.1; concluding that the
relevant operator $L$, whose spectral stability we seek to study is:

$L=i\begin{pmatrix}0&1\\\
L_{1}\left[\partial_{\mathscr{\rho}}^{2},\partial_{\mathscr{\rho}};\boldsymbol{\mathfrak{q}},\rho\right]&L_{2}\left[\partial_{\rho};\boldsymbol{\mathfrak{q}},\rho\right]\end{pmatrix}\,.$
(96)

With respect to the boundary conditions, we note that, as indicated in section
3.1, the outgoing behavior is automatically satisfied as long as we demand
regularity on the horizon. On the other hand, normalizable boundary conditions
have to be imposed manually, discarding the leading mode on the AdS boundary.
In order to identify the aforementioned mode, we solve the equation of motion
(65) near the AdS boundary using a Frobenius series, obtaining the following
solution:

$\phi(\mathscr{t},\mathbf{k},\rho)=e^{-i\omega\mathscr{t}}\left\\{c_{-}(\rho-1)^{\Delta_{-}}(1+...)+c_{+}(\rho-1)^{\Delta_{+}}(1+...)\right\\}\,,$
(97)

with $\Delta_{\pm}=2\pm\sqrt{4+m^{2}l^{2}}$. Imposing normalizable boundary
conditions corresponds to fixing $c_{-}=0$.

Regarding the energy norm, we have the following expression for the inner
product:242424We have eliminated all constant prefactors from the energy norm
(98) as they play no role in the operator norm.

$\expectationvalue{\phi_{1},\phi_{2}}_{E}=\int_{\mathscr{t}=\text{const}}d\rho\,\,\begin{pmatrix}\overline{\phi_{1}}(\mathscr{t},\mathbf{k},\rho)&\overline{\psi_{1}}(\mathscr{t},\mathbf{k},\rho)\end{pmatrix}\mathcal{G}\left[\overleftarrow{\partial}_{\rho},\overrightarrow{\partial}_{\rho};\boldsymbol{\mathfrak{q}},\rho\right]\begin{pmatrix}\phi_{2}(\mathscr{t},\mathbf{k},\rho)\\\\[6.0pt]
\psi_{2}(\mathscr{t},\mathbf{k},\rho)\end{pmatrix}\,,$ (98)
$\mathcal{G}\left[\overleftarrow{\partial}_{\rho},\overrightarrow{\partial}_{\rho};\boldsymbol{\mathfrak{q}},\rho\right]=\begin{pmatrix}\frac{m^{2}l^{2}}{(1-\rho)^{5}}+\frac{\boldsymbol{\mathfrak{q}}^{2}}{(1-\rho)^{3}}+\overleftarrow{\partial}_{\rho}\frac{\mathscr{f}(\rho)}{(1-\rho)^{3}}\overrightarrow{\partial}_{\rho}&0\\\\[6.0pt]
0&\frac{2-\mathscr{f}(\rho)}{(1-\rho)^{3}}\end{pmatrix}\,,$ (99)

derived from the subsequent energy-momentum tensor:

$T_{MN}=\nabla_{M}\phi\nabla_{N}\phi-\frac{1}{2}g_{MN}\left[\nabla_{S}\phi\nabla^{S}\phi+m^{2}\phi^{2}\right]\,.$
(100)

In order to ensure that the energy norm (98) is positive definite, while also
enforcing the adequate asymptotic behavior for the QNMs, we chose to work in
the space of regular functions with the following behavior on the AdS
boundary:

$\phi(\mathscr{t},\mathbf{k},\rho\rightarrow
1)=A(\mathscr{t},\mathbf{k})\,(\rho-1)^{n}\,,\qquad n>2\,.$ (101)

With this choice, for masses above the BF bound, the leading mode in (97) is
automatically discarded and the energy norm is positive definite
Breitenlohner:1982bm . It is worth pointing out that the chosen function
space, besides being mathematically convenient, is also the physically
relevant one as it contains all possible asymptotic behaviors for scalar QNMs
with masses above the BF bound in the standard quantization.

In view of the above, in order to select functions belonging to the chosen
Hilbert space, it is convenient to work with the rescaled fields

$\begin{pmatrix}\tilde{\phi}\\\
\tilde{\psi}\end{pmatrix}=(1-\rho)^{-2}\begin{pmatrix}\phi\\\
\psi\end{pmatrix}\,,$ (102)

for which imposing the asymptotic behavior (101) amounts to fixing Dirichlet
boundary conditions on the AdS boundary.

Lastly, note that, given the expression of the adjoint of $L$ with respect to
the energy norm $L^{\dagger}$:252525The explicit computation of $L^{\dagger}$
can be found in appendix A.1.

$L^{\dagger}=L+\begin{pmatrix}0&0\\\\[6.0pt]
0&-2i\,\delta(\rho)\,\frac{\mathscr{f}(\rho)-1}{\mathscr{f}(\rho)-2}\end{pmatrix}\,,$
(103)

we can conclude that $L$ is non-normal and that this non-normality is
associated with the existence of a horizon. This matches our initial physical
intuition, the system is non-conservative because the fluctuations eventually
fall into the horizon and die out.

###  4.2 Transverse Gauge Field

The linearized equations of motion for a gauge field with action (95) are:

$\nabla^{M}F_{MN}=0\,.$ (104)

Decomposing into Fourier modes $A_{\mu}(\mathscr{t},\mathbf{k},\rho)$ and
assuming the momentum $\mathbf{k}$ to be oriented in the $\mathscr{x}_{3}$
direction, the equations of motion decouple into two sectors: transverse
$\\{A_{1},A_{2}\\}$ and longitudinal
$\\{A_{\mathscr{t}},A_{3},A_{\mathscr{r}}\\}$. The members of the transverse
channel transform as vectors under the unbroken $O(2)$ symmetry, while the
components of the longitudinal sector are invariants Kovtun:2005ev .

In the present work, we only study the transverse sector, i.e., we consider
QNMs $a(\mathscr{t},\mathbf{k},\rho)$ satisfying the following eigenvalue
problem:

$\displaystyle\mathfrak{w}\begin{pmatrix}a(\mathscr{t},\mathbf{k},\rho)\\\
\alpha(\mathscr{t},\mathbf{k},\rho)\end{pmatrix}$
$\displaystyle=i\begin{pmatrix}0&1\\\
L_{1}\left[\partial_{\mathscr{\rho}}^{2},\partial_{\mathscr{\rho}};\boldsymbol{\mathfrak{q}},\rho\right]&L_{2}\left[\partial_{\rho};\boldsymbol{\mathfrak{q}},\rho\right]\end{pmatrix}\begin{pmatrix}a(\mathscr{t},\mathbf{k},\rho)\\\
\alpha(\mathscr{t},\mathbf{k},\rho)\end{pmatrix}\,,$ (105) $\displaystyle
L_{1}\left[\partial_{\rho}^{2},\partial_{\rho};\boldsymbol{\mathfrak{q}},\rho\right]$
$\displaystyle=\left[\mathscr{f}(\rho)-2\right]^{-1}\left[\boldsymbol{\mathfrak{q}}^{2}-(1-\rho)\left(\frac{\mathscr{f}(\rho)}{1-\rho}\right)^{\prime}\partial_{\rho}-\mathscr{f}(\rho)\partial_{\rho}^{2}\right]\,,$
(106) $\displaystyle
L_{2}\left[\partial_{\rho};\boldsymbol{\mathfrak{q}},\rho\right]$
$\displaystyle=\left[\mathscr{f}(\rho)-2\right]^{-1}\left[(1-\rho)\left(\frac{\mathscr{f}(\rho)-1}{(1-\rho)}\right)^{\prime}+2\left(\mathscr{f}(\rho)-1\right)\partial_{\rho}\right]\,,$
(107)

where we used again the dimensionless variables $\mathfrak{w}=z_{h}\omega$,
$\boldsymbol{\mathfrak{q}}=z_{h}\mathbf{k}$ and introduced the auxiliary field
$\alpha(\mathscr{t},\mathbf{k},\rho)=z_{h}\partial_{\mathscr{t}}a(\mathscr{t},\mathbf{k},\rho)=-i\mathfrak{w}a(\mathscr{t},\mathbf{k},\rho)$.

It is important to stress that the transverse sector is gauge invariant, as we
can eliminate the $\\{\mathscr{x}_{1},\mathscr{x}_{2}\\}$ dependence of all
functions given the symmetries of the problem. Explicitly, we have that the
gauge transformation of the transverse gauge field is the following:

$a\rightarrow a+d\chi=a$ (108)

where, as indicated, we disregard any possible dependence of $\mathscr{\chi}$
on $\\{\mathscr{x}_{1},\mathscr{x}_{2}\\}$. This ensures that the
pseudospectrum and the condition numbers only probe stability under gauge-
invariant perturbations. In fact, the gauge invariant electric field is simply
$E_{1,2}=i\mathfrak{w}A_{1,2}$.

As with the scalar, outgoing boundary conditions are automatically satisfied
demanding regularity on the horizon, while normalizable boundary conditions
have to be imposed explicitly. Solving the equations of motion (104) with a
Frobenius series, we obtain the following near-AdS boundary behavior:

$a(\mathscr{t},\mathbf{k},\rho)=e^{-i\omega\mathscr{t}}\left\\{c_{-}(1+...)+c_{+}(1-\rho)^{2}(1+...)\right\\}\,,$
(109)

and thus, imposing normalizable boundary conditions corresponds to fixing
$c_{-}=0$.

With regard to the energy norm, we recall that the energy-momentum tensor for
a gauge field with action (95) is given by:

$T_{MN}=F_{MS}\tensor{F}{{}_{N}^{S}}-\frac{1}{4}g_{MN}F_{SP}\tensor{F}{{}^{S}{}^{P}}\,,$
(110)

which yields the following inner product:

$\expectationvalue{a_{1},a_{2}}_{E}=\int_{\mathscr{t}=\text{const}}d\rho\,\,\begin{pmatrix}\overline{a_{1}}(\mathscr{t},\mathbf{k},\rho)&\overline{\alpha_{1}}(\mathscr{t},\mathbf{k},\rho)\end{pmatrix}\mathcal{G}\left[\overleftarrow{\partial}_{\rho},\overrightarrow{\partial}_{\rho};\boldsymbol{\mathfrak{q}},\rho\right]\begin{pmatrix}a_{2}(\mathscr{t},\mathbf{k},\rho)\\\\[6.0pt]
\alpha_{2}(\mathscr{t},\mathbf{k},\rho)\end{pmatrix}\,,$ (111)
$\mathcal{G}\left[\overleftarrow{\partial}_{\rho},\overrightarrow{\partial}_{\rho};\boldsymbol{\mathfrak{q}},\rho\right]=\begin{pmatrix}\frac{\boldsymbol{\mathfrak{q}}^{2}}{1-\rho}+\overleftarrow{\partial}_{\rho}\frac{\mathscr{f}(\rho)}{1-\rho}\overrightarrow{\partial}_{\rho}&0\\\\[6.0pt]
0&\frac{2-\mathscr{f}(\rho)}{1-\rho}\end{pmatrix}\,.$ (112)

In this case, as we have no mass term, the energy norm is always positive
definite. Thus, we consider the function space comprised of functions
satisfying outgoing boundary conditions on the horizon with the following
behavior on the AdS boundary:

$a(\mathscr{t},\mathbf{k},\rho\rightarrow
1)=A(\mathscr{t},\mathbf{k})(\rho-1)^{n}\,,\qquad n>1\,.$ (113)

This ensures the desired asymptotic behavior for the QNMs while also
guaranteeing that the norm is non-divergent. Similarly to the case of the
scalar, we find it more convenient to work with the rescaled fields:

$\begin{pmatrix}\tilde{a}\\\
\tilde{\alpha}\end{pmatrix}=(1-\rho)^{-1}\begin{pmatrix}a\\\
\alpha\end{pmatrix}\,,$ (114)

for which imposing the asymptotic behavior (113) amounts to fixing Dirichlet
boundary conditions on the AdS boundary.

To conclude, we note that, identically to what we found for the scalar field,
$L$ is non-normal with respect to the energy norm. Its adjoint with respect to
the aforementioned norm is given by

$L^{\dagger}=L+\begin{pmatrix}0&0\\\\[6.0pt]
0&-2i\delta(\rho)\frac{\mathscr{f}(\rho)-1}{\mathscr{f}(\rho)-2}\end{pmatrix}\,,$
(115)

which once again, relates the non-normality to the existence of a
horizon.262626The explicit computation of $L^{\dagger}$ can be found in
appendix A.2.

##  5 Numerical method

We approach the stability analysis numerically; discretizing the radial
coordinate $\rho$ in a Chebyshev grid with points:

$\rho_{j}=\frac{1}{2}\left[1-\cos\left(\frac{j\pi}{N}\right)\right]\,,\qquad
j=0,1,...,N\,,$ (116)

which in turn allows us to discretize the differential operators using the
corresponding Chebyshev differentiation matrices Trefethen:2000 . This choice
is equivalent to approximating the QNMs by a series of Chebyshev polynomials
$T_{n}$ with $n=\\{0,1,...,N\\}$. As a rule of thumb, this allows us to
approximate at most $N/2$ QNFs Trefethen:2000 ; boyd .

Regarding the energy norm we proceed as indicated in Jaramillo:2020tuu  and
construct a $G_{E}$ matrix defined as the discretized version of the original
energy norm. Labelling by $u$ the rescaled scalar doublet
$\left(\tilde{\phi},\,\tilde{\psi}\right)^{T}$ and the rescaled gauge field
doublet $\left(\tilde{a},\,\tilde{\alpha}\right)^{T}$, we construct $G_{E}$
such that:

$\displaystyle\lim_{N\rightarrow\infty}u_{N}^{*}G_{E}u_{N}=\expectationvalue{u,u}_{E}\,,$
(117)

where $u_{N}$ is the vector arising form the discretization of $u$.272727There
is an important subtlety in the construction of $G_{E}$ that should be
addressed. One needs to construct $G_{E}$ in a grid of at least twice the size
of the original one and, at the end, interpolate back. This ensures
consistency with the discretization process: when working in a grid with $N+1$
points the maximum resolution is given by polynomials of degree $N+1$. Thus,
for the discretized norm to be exact for polynomials of such degree, one needs
to construct it on a grid with at least $2(N+1)$ points. A detailed discussion
on the construction of the $G_{E}$ matrix can be found in appendix B.

To conclude the discretization process, we need to numerically impose the
function spaces introduced in the previous section. Note that, in any
numerical method, we can only get regular solutions. Consequently, the
discretization immediately selects the space of regular functions, and thus we
only need to ensure the adequate behavior on the AdS boundary. This
corresponds to imposing Dirichlet boundary conditions for the rescaled scalar
and gauge fields, which can be achieved by removing the rows and columns
corresponding to the AdS boundary from all discretized operators, including
the $G_{E}$ matrix.282828This is equivalent to reducing the space in which the
matrices act to that of vectors vanishing on the AdS boundary.

With the original system fully discretized, we can proceed to study the
spectral stability. We use Wolfram Engine to compute condition numbers and
pseudospectra as indicated in theorem 2.7. It is particularly convenient to
rewrite the problem in terms of the matrix $\ell^{2}$-norm as it reduces the
computation of the pseudospectrum to obtaining the smallest eigenvalue of a
Hermitian matrix; which we locate using Arnoldi iteration (see e.g. ch. 28 of
Trefethen:2005 ). With this procedure, we achieve
$\mathcal{O}\left(N^{2}\right)$ runtime for each point $z$ in the complex
plane where we compute the norm of the resolvent.

In order to gain greater insight into the nature of the (in)stability, we also
explore the selective pseudospectra associated with local potential
perturbations to the original equations of motion. Concretely, we consider
perturbations to equations (65) and (104) of the form:

$\left[-\nabla_{M}\nabla^{M}+m^{2}+\frac{V(\rho)}{l^{2}}\right]\phi=0\,,\qquad\nabla^{M}F_{MN}-\frac{V(\rho)}{l^{2}}A_{N}=0\,,$
(118)

where, in order to preserve the asymptotic behavior on the AdS boundary, we
choose potentials vanishing on the aforementioned boundary $V(1)=0$. Note that
the potential term added to the gauge field dynamics does, in general, break
gauge invariance. Nonetheless, as a transverse gauge field $A_{M}$ is gauge
invariant, the chosen potential term preserves the gauge symmetry.

The main difference between full and selective pseudospectra lies in the kind
of stability each of them probes. The former explores the stability under
generic bounded perturbations, while the latter only considers local
perturbations. Physically, these local perturbations can be interpreted as
effective interactions arising from small deviations from the perfect SAdS4+1
background. Therefore, probing the stability under them is akin to determining
the underlying model dependence of the system.292929In general, we consider
that corrections arising from deformations of the perfect SAdS4+1 background
(associated, for instance, with the existence of extra fields or non-minimal
couplings in the action) appear as non-local interactions. Nonetheless, we
assume that the leading order behavior can be modelled with a local potential
(a very common picture in effective field theories). Under this
interpretation, and recalling that we have no control over the nature of the
perturbations entering the full pseudospectrum, we consider the selective
pseudospectrum to be slightly more physically motivated.

The selective pseudospectrum is computed using definition 2.2 with randomly
generated potential perturbations constructed as diagonal matrices with random
entries and normalized to a given size. This analysis is complemented with the
computation of the QNFs for the perturbed system with the following
deterministic potentials:

$\displaystyle V_{1}(\rho)$ $\displaystyle=A_{1}(1-\rho)\cos(2\pi\rho)\,,$
(119a) $\displaystyle V_{2}(\rho)$
$\displaystyle=A_{2}(1-\rho)\cos(90\pi\rho)\,,$ (119b) $\displaystyle
V_{3}(\rho)$
$\displaystyle=A_{3}(1-\rho)\left\\{1-\tanh\left[20\rho\right]\right\\}\,,$
(119c) $\displaystyle V_{4}(\rho)$
$\displaystyle=A_{4}(1-\rho)\left\\{1-\tanh\left[20(1-\rho)\right]\right\\}\,,$
(119d)

which shed light on a few interesting regimes. With $V_{1}$ and $V_{2}$, we
probe the effect of long and short $\rho$-wavelength (wavelength in the $\rho$
direction) perturbations; while with $V_{3}$ and $V_{4}$, we analyze the
stability under localized perturbations near the horizon (IR of the dual QFT)
and the boundary (UV of the dual QFT). The $\\{A_{i}\\}$ are normalization
constants to fix the magnitude of the perturbation. We plot these potentials
in figure 3.

Figure 3: Deterministic potentials (119) with $A_{i}=1$. Recall that the
horizon is at $\rho=0$ and the boundary at $\rho=1$.

The grid establishes a cutoff in the $\rho$-wavelength of the perturbations.
Thus, as a byproduct of the discretization, we are only truly sensitive to the
stability under perturbations above the cutoff.

Finally, to stress the nontrivial nature of the norm, we also consider the
pseudospectrum in the $L^{2}$-norm

$\norm{u}_{L^{2}}=\int
d\rho\,\,u^{*}(\mathscr{t},\mathbf{k},\rho)u(\mathscr{t},\mathbf{k},\rho)\,.$
(120)

Although the physical relevance of this norm for the present problem is less
clear, it might be informative to compare the general features of the
pseudospectra in the energy norm to those in the $L^{2}$-norm. We direct the
interested reader to Appendix C where we present the pseudospectrum analysis
in the $L^{2}$-norm.

##  6 Results

Here, we present the results of the analysis described in the previous
sections. Our numerical simulations are performed in a grid of 120 points with
a precision $5\times$MachinePrecision.

(a) Scalar with $m^{2}l^{2}=0$ and $\boldsymbol{\mathfrak{q}}^{2}=0$.

(b) Transverse Gauge Field with $\boldsymbol{\mathfrak{q}}^{2}=0$.

Figure 4: Convergence test for QNFs computed with $5\times$MachinePrecision.
The error is defined as
$\left|1-\frac{|\mathfrak{w}|}{|\mathfrak{w}_{\text{ref}}|}\right|$, with
$\mathfrak{w}_{\text{ref}}$ the reference value obtained with
$10\times$MachinePrecision in grid of $400$ points. The missing points are
such that their error is below the numerical accuracy.

(a) Scalar with $m^{2}l^{2}=0$ and $\boldsymbol{\mathfrak{q}}^{2}=0$.

(b) Transverse Gauge Field with $\boldsymbol{\mathfrak{q}}^{2}=0$.

Figure 5: Convergence test for QNFs computed with MachinePrecision. The error
is defined as
$\left|1-\frac{|\mathfrak{w}|}{|\mathfrak{w}_{\text{ref}}|}\right|$, with
$\mathfrak{w}_{\text{ref}}$ the reference value obtained with
$10\times$MachinePrecision in grid of $400$ points. The large error of the
fourth QNF hints towards spectral instability.

We require large grids and high precision to ensure that the numerical
procedure does not heavily affect our results. Physically, the discretization
and numerical round-offs can be understood as perturbations to the original
problem and, as we want to analyze the stability of the latter, we need to be
especially careful and avoid undesired effects associated with them. As we are
mainly interested in the stability of the three or four lowest-lying QNFs, in
figure 4 we test the validity of our numerics by computing those QNFs in our
setup and comparing them to the ones obtained using $10\times$MachinePrecision
in a grid of 400 points.303030Henceforth, we abuse notation and denote by QNF
the pair $\\{\mathfrak{w},-\mathfrak{w}^{*}\\}$. As the equations of motion
are invariant under simultaneous conjugation and multiplication by $-1$, the
spectrum (and pseudospectrum) is symmetric under reflections along the
imaginary axis and consequently the eigenvalues always appear in pairs
$\\{\mathfrak{w},-\mathfrak{w}^{*}\\}$. Furthermore, we choose to order the
QNFs by their imaginary part such the first QNF is the one with the largest
imaginary part. Thus, we conclude that with our choice of parameters we are
ensuring that the effect of the numerics on the first four QNFs is smaller
than $10^{-50}$ in all cases. It is interesting to note that the need for such
large grids and precisions is a good indicator of the underlying instability
of the problem. Small perturbations associated with the numerics have effects
much larger than their typical scale (see figure 5).

###  6.1 Real Scalar in $\text{SAdS}_{4+1}$

(a) $m^{2}l^{2}=0$, $\mathfrak{q}=0$.

(b) $m^{2}l^{2}=0$, $\mathfrak{q}=10$.

(c) $m^{2}l^{2}=-3$, $\mathfrak{q}=0$.

(d) $m^{2}l^{2}=-3$, $\mathfrak{q}=10$.

Figure 6: Close-up of the scalar pseudospectrum in the energy norm around the
first QNF for different values of $\mathfrak{q}$ and $m^{2}l^{2}$. The red dot
corresponds to the QNF, the white lines represent the boundaries of various
full $\varepsilon$-pseudospectra, and the dashed blue circle symbolizes a
circle with a radius of $10^{-1}$ centered on the QNF. The heat map
corresponds to the logarithm in base 10 of the inverse of the resolvent, while
the blue and cyan dots indicate selective $\varepsilon$-pseudospectra computed
with random local potential perturbations of size $10^{-1}$ and $10^{-3}$,
respectively. Remarkably, only in (a) do we observe instability of the first
QNF under local potential perturbations.

In figures 6 and 7 we present full and selective pseudospectra and the
corresponding condition numbers in the energy norm for different values of
$m^{2}l^{2}$ and $\mathfrak{q}=\sqrt{\boldsymbol{\mathfrak{q}}^{2}}$. Recall
that throughout this section we are working in units of $z_{h}=(\pi\,T)^{-1}$.
The $\varepsilon$-pseudospectra exhibit extended open regions denoting
instability.313131Note that all QNFs are contained within some region of the
$\varepsilon$-pseudospectra. However, due to limited resolution in the plots,
not all regions of the $\varepsilon$-pseudospectrum are observable. For
example, in the full pseudospectrum shown in figure 7(a), we cannot appreciate
the $10^{-7}$-pseudospectrum around the first QNF observable in figure 6(a).
The same limitation applies to the selective pseudospectra, where small
regions are covered by the red dots representing the QNFs (e.g., the selective
$10^{-3}$-pseudospectrum observed around the first QNF in figure 6(a) is
covered by a red dot in figure 7(a)). As observed in asymptotically flat
Jaramillo:2020tuu ; Destounis:2021lum  and de Sitter spacetimes Sarkar:2023rhp
, the instability increases the further away the QNFs are from the real axis.
In the dual quantum field theory this implies that excitations are
increasingly spectrally unstable the more short-lived they are.

Interestingly, we find that one needs perturbations of size $\sim
10^{-0.4}\sim 0.4$ to drive the QNFs to the upper half of the complex plane,
as indicated by the full pseudospectra in figure 7. In our units the typical
distance between eigenvalues is $\sim 2$, this implies that all backgrounds
arising as a small deviation from SAdS4+1 are stable, i.e., they do not have
exponentially growing QNMs.

(a) $m^{2}l^{2}=0$, $\mathfrak{q}=0$.

(b) $m^{2}l^{2}=0$, $\mathfrak{q}=10$.

(c) $m^{2}l^{2}=-3$, $\mathfrak{q}=0$.

(d) $m^{2}l^{2}=-3$, $\mathfrak{q}=10$.

Figure 7: Scalar pseudospectrum in the energy norm for different values of
$\mathfrak{q}$ and $m^{2}l^{2}$. In the lower panels, we present selective and
full pseudospectra. The red dots represent the (unperturbed) QNFs in the
typical “Christmas Tree” configuration. The white lines denote the boundaries
of different full $\varepsilon$-pseudospectra. The heat map corresponds to the
logarithm in base 10 of the inverse of the resolvent, while the blue, cyan,
green, and yellow dots indicate different selective
$\varepsilon$-pseudospectra computed with random local potential perturbations
of size $10^{-1}$, $10^{-3}$, $10^{-5}$, and $10^{-7}$; respectively. In the
upper panels, we represent the condition numbers. Most notably, for small
values of $\varepsilon$, the full $\varepsilon$-pseudospectra present open
regions containing multiple QNFs, which signals spectral instability.

(a) Imaginary part of the QNFs for $\mathfrak{q}=0$.

(b) Condition numbers for $\mathfrak{q}=0$.

(c) Imaginary part of the QNFs for $\mathfrak{q}=10$.

(d) Condition numbers for $\mathfrak{q}=10$.

Figure 8: Mass dependence in the energy norm for the scalar fluctuations.

(a) Imaginary part of the QNFs for $m^{2}l^{2}=0$.

(b) Condition numbers for $m^{2}l^{2}=0$.

(c) Imaginary part of the QNFs for $m^{2}l^{2}=-3$.

(d) Condition numbers for $m^{2}l^{2}=-3$.

Figure 9: Momentum dependence in the energy norm for the scalar fluctuations.

Focusing now on the selective pseudospectrum, we note that its structure seems
to resemble the structure of the full pseudospectrum. Nonetheless, the scales
clearly do not match; the selective pseudospectrum shows significantly more
stability. Hence, we can conclude that most of the instability arises from
perturbations which cannot be interpreted as local potential perturbations.

Regarding the mass and momentum dependence, stability decreases with mass and
increases with momentum (see figures 8 and 9). Such behavior can be understood
in the context of the instability being related to the distance to the real
axis: the imaginary part of the QNFs is reduced as mass (momentum) decreases
(increases), and thus we observe more stability. For the QFT this implies that
high-momentum fluctuations of operators with small scaling dimension are more
stable.

(a) $m^{2}l^{2}=0$, $\mathfrak{q}=0$.

(b) $m^{2}l^{2}=0$, $\mathfrak{q}=10$.

(c) $m^{2}l^{2}=-3$, $\mathfrak{q}=0$.

(d) $m^{2}l^{2}=-3$, $\mathfrak{q}=10$.

Figure 10: Effect on the first scalar QNF of the deterministic perturbations
(119) with size $\norm{V_{i}}_{E}=10^{-1}$. The unperturbed QNF is shown in
red, while the perturbed QNFs are depicted in blue. The dashed blue line
represents the circle of radius $10^{-1}$ centered in the unperturbed QNF.
Remarkably, only for $\mathfrak{q}=0$ do we observe instability under near-
horizon perturbations ($V_{3}$).

(a) $m^{2}l^{2}=0$, $\mathfrak{q}=0$.

(b) $m^{2}l^{2}=0$, $\mathfrak{q}=10$.

(c) $m^{2}l^{2}=-3$, $\mathfrak{q}=0$.

(d) $m^{2}l^{2}=-3$, $\mathfrak{q}=10$.

Figure 11: Effect on the spectrum of the scalar of the deterministic
perturbations (119) with size $\norm{V_{i}}_{E}=10^{-1}$. In the lower panels
we present the spectra and in the upper ones the condition numbers for the
lowest QNFs. The unperturbed QNFs are shown in red, while the perturbed ones
are depicted in blue. The plotted region of the spectra is stable under long
$\rho$-wavelength perturbations ($V_{1}$) and near-boundary perturbations
($V_{4}$). Note also that, as indicated by the condition numbers, the
perturbed QNFs associated with $V_{2}$ and $V_{3}$ are more stable than the
unperturbed QNFs.

We now focus on the first QNF (see figure 6), which has the smallest imaginary
part and thus dominates the low energy spectrum of the QFT. Remarkably, we
find that for $m^{2}l^{2}=0$ at $\mathfrak{q}=0$ that QNF is unstable under
local potential perturbations as the selective $\varepsilon$-pseudospectra
extends beyond the circle of radius $\varepsilon$ centered around the QNF. To
further explore the nature of this instability, in figure 10 we analyze the
effect of the deterministic potential perturbations (119) on the first QNF. We
observe that the instability is only present for near-horizon perturbations,
concluding that it is mainly associated with deformations of the IR of the
dual QFT. Moreover, we also find that, although we did not have statistical
evidence in the selective pseudospectra, the first QNF is actually unstable
under local potential perturbations for $m^{2}l^{2}=-3$ (figure 10(c)). At
$\mathfrak{q}=10$ we find that the QNM is stable, as expected by the enhanced
stability at large $\mathfrak{q}$.

Continuing our exploration of the effect of the deterministic potentials, in
figure 11 we present the spectrum and the condition numbers for the first 3
perturbed QNFs. We observe that large $\rho$-wavelength ($V_{1}$) and near-
boundary ($V_{4}$) perturbations have very little effect on the spectrum. On
the other hand, short $\rho$-wavelength ($V_{2}$) and near-horizon
perturbations ($V_{3}$) displace the QNFs significantly, generating new
branching structures. Remarkably, these perturbed QNFs are more stable than
the unperturbed ones as indicated by the condition numbers plotted in the
upper panels of figure 11. Thus, the spectra of poles of Green’s functions in
strongly coupled holographic field theories might have some preference for
branching structures over the unperturbed “Christmas Tree” configuration.

We end this section by briefly addressing the results for the pseudospectrum
of the scalar field in the $L^{2}$-norm. We refer the reader to Appendix C and
here we simply point out that the pseudospectra are qualitatively similar to
those in the energy norm. However, in the $L^{2}$-norm we find that the
stability is enhanced under local potential perturbations and reduced under
generic perturbations.

###  6.2 Transverse Gauge Field in $\text{SAdS}_{4+1}$

In this section we present the pseudospectrum of the fluctuations of the
transverse gauge field introduced in section 4.2. As we detail in the
following, we observe the same generic features as for the real scalar. In
particular we find the QNFs of the gauge field to be unstable.

(a) $\mathfrak{q}=0$.

(b) $\mathfrak{q}=10$.

Figure 12: Close-up of the transverse gauge field pseudospectrum in the energy
norm around the first QNF for different values of $\mathfrak{q}$. The red dot
corresponds to the QNF, the white lines represent the boundaries of various
full $\varepsilon$-pseudospectra, and the dashed blue circle symbolizes a
circle with a radius of $10^{-1}$ centered on the QNF. The heat map
corresponds to the logarithm in base 10 of the inverse of the resolvent, while
the blue dots indicate selective $\varepsilon$-pseudospectra computed with
random local potential perturbations.

(a) $\mathfrak{q}=0$.

(b) $\mathfrak{q}=10$.

Figure 13: Transverse gauge field pseudospectrum in the energy norm for
different values of $\mathfrak{q}$. In the lower panels, we present selective
and full pseudospectra. The red dots represent the QNFs, and the white lines
denote the boundaries of different full $\varepsilon$-pseudospectra. The heat
map corresponds to the logarithm in base 10 of the inverse of the resolvent,
while the blue, cyan, green, and yellow dots indicate different selective
$\varepsilon$-pseudospectra computed with random local potential perturbations
of size $10^{-1}$, $10^{-3}$, $10^{-5}$, and $10^{-7}$; respectively. In the
upper panels, we represent the condition numbers.

(a) $\mathfrak{q}=0$.

(b) $\mathfrak{q}=10$.

(c) $\mathfrak{q}=0$.

(d) $\mathfrak{q}=10$.

Figure 14: Effect on the transverse gauge field spectrum of the deterministic
perturbations (119) with size $\norm{V_{i}}_{E}=10^{-1}$. The unperturbed QNF
is shown in red, while the perturbed QNFs are depicted in blue. In (a) and (b)
the dashed blue line represents the circle of radius $10^{-1}$ centered in the
unperturbed QNF. In the upper panels of (c) and (d) we present the condition
numbers for the lowest QNFs.

(a) Imaginary part of the QNFs.

(b) Condition numbers.

Figure 15: Momentum dependence in the energy norm for the transverse gauge
field.

First, in figures 12 and 13 we display the pseudospectrum of the transverse
gauge field. The presence of open contour lines around the QNFs indicates that
they are unstable under generic perturbations. However, for the first QNF we
have not found statistical evidence of its instability under random local
potential perturbations (see the selective pseudospectra in figure 12(a)). It
is interesting to note that this is similar to what we found for the scalar
field of mass $m^{2}l^{2}=-3$ (see figure 6(c)). From the conformal field
theory point of view both of these fields correspond to operators of conformal
dimension three.

In order to further characterize the nature of the instability, in figure 14
we show the effect of the deterministic perturbations (119) on the QNF
spectrum. At zero momentum we find that the spectrum is unstable under the
perturbations (119b) and (119c). In particular in figure 14(a) the first QNF
is shown to be unstable under near-horizon perturbations (119c). Moreover,
from studying the evolution of the condition numbers shown in that figure, one
can conclude that near-horizon and short $\rho$-wavelength perturbations
result in less unstable QNFs.

As in the previous section, it is also interesting to analyze the dependence
of the stability of the QNFs on the momentum of the fluctuation. Therefore, in
figure 15 we plot the momentum dependence of the QNFs and their condition
numbers. The QNFs get closer to the real axis as momentum increases and, as
for the real scalar, we observe that their condition numbers decrease in the
process (although staying above one for the range of momentum studied).

Despite the qualitative similarities, it is worth noting that the actual
values of the pseudospectra do not match the cases studied in the previous
subsection. Specifically, we observe that the gauge field is more stable than
the massless scalar and less stable than the scalar with mass $m^{2}l^{2}=-3$.
This agrees with the instability being related to the imaginary part of the
QNFs. The gauge field QNFs are closer to the real axis than those of the
massless scalar and further away than those of the scalar with mass
$m^{2}l^{2}=-3$.

Finally, it is illustrative to analyze the pseudospectrum in the $L^{2}$-norm.
We present our results in appendix C. Again, we observe that in the
$L^{2}$-norm the stability is enhanced under local potential perturbations and
reduced under generic perturbations.

##  7 Conclusions

In this work we have initiated the study of the spectral stability of
quasinormal modes in asymptotically AdS black hole geometries. Gauge/gravity
duality maps the frequencies of those modes (QNFs) to the spectrum of
collective excitations of a strongly coupled quantum many-body system.
Therefore, probing the spectral stability of QNFs is akin to probing the
stability of the spectra of strongly coupled holographic field theories.

Following previous works in asymptotically flat Jaramillo:2020tuu ;
Destounis:2021lum  and de Sitter spacetimes Sarkar:2023rhp , we probed
spectral stability through pseudospectra and condition numbers. To construct
the pseudospectra we needed to choose a physically motivated norm and cast the
equations of motion for the quasinormal modes in the form of an eigenvalue
problem. Consequently, in section 3 we proposed a generic approach valid for
static black brane backgrounds, a family of spacetimes very prominent in the
literature. Our approach relied on two cornerstones. First, the regular
coordinates, which enabled us to write the equations of motion as an
eigenvalue problem while also translating outgoing boundary conditions to
regularity conditions on the horizon. And second, the energy norm, which
represents the physically relevant notion of size for the QNMs.

After discussing the general approach, we proceeded to study the stability of
a real scalar field and a transverse gauge field in SAdS4+1. As previously
observed in asymptotically flat Jaramillo:2020tuu ; Destounis:2021lum  and dS
spacetimes Sarkar:2023rhp , we found that the QNFs are unstable and that the
instability is associated with the existence of a horizon. Moreover, we
observed that QNFs are less unstable the closer they are to the real axis, a
feature also present in the aforementioned spacetimes.

Regarding the scalar, we found an increase in stability at high momenta and
small masses. This implies that in the dual QFT the instability is milder for
high-momentum fluctuations of operators with small mass dimension. In order to
further study the nature of the instability, we probed the spectrum of QNFs
with local potential perturbations which could arise as effective interactions
representing deviations from the exact SAdS4+1 background. Remarkably, we
inferred that only a small portion of the total instability is associated with
these perturbations. Moreover, we concluded that small perturbations cannot
drive the QNFs to the upper half of the complex plane. We also found that for
small momenta, the first QNF is unstable, with the instability arising from
localized perturbations near the horizon. Thus, small deformations of the IR
can destabilize even the lowest QNF.

With respect to the gauge field, we concluded that its pseudospectrum is
qualitatively very similar to that of the scalar field. This is even more
pronounced for the scalar of mass $m^{2}l^{2}=-3$. In that case both the
scalar and the transverse gauge field correspond to operators with the same
conformal dimension. The resemblance between the pseudospectra of the scalar
and the gauge field might be due to the fact that the spectra of both fields
are very similar. Exploring this potential connection would indicate a degree
of universality, which would be intriguing to investigate in future research.

We have proven that AdS spacetimes also present spectral instability
associated with the existence of an event horizon. This supports the physical
intuition that the observed instability is a direct consequence of the non-
conservative character of black holes and thus is independent of the
asymptotic behavior of the spacetime.

With regard to gauge/gravity duality, our results imply that thermal
excitation spectra of strongly coupled quantum field theories dual to SAdS4+1
(or some small deviation from it) are unstable, i.e., small perturbations to
the original model have a large effect on the spectra. This suggests that a
mathematical model would not be able to accurately capture the actual spectra
of the real physical system. That said, it is important to remark that the
spectral instability decreases with the decay width of the fluctuations. We
expect the imprint of the aforementioned instability to be less pronounced for
modes near the real axis. In a practical setting it implies that, at most, one
might hope to be able to capture the leading order behavior (dominated by the
longest-lived excitation). On the other hand, it might not be possible to
model the subleading features of decay trustfully as the higher QNFs are
increasingly unstable.

In the realm of gauge/gravity duality, pseudospectral analysis in asymptotic
AdS geometries has many applications with potentially deep implications. We
list some of them in the following.

  * •

Study the stability of hydrodynamic modes. As we have seen, the closer a QNF
is to the real axis, the less unstable it becomes. In particular that might
indicate that hydrodynamic modes enjoy privileged stability properties in the
pseudospectrum. On the other hand, if this effect is not sufficiently
pronounced it would be interesting to see whether small local perturbations at
zero momentum could drive the hydro mode’s QNF to the upper half plane
indicating instabilities of flow patterns. In fact such instabilities have
been argued to arise in Navier-Stokes equations and explain the early onset of
turbulence trefethen1993hydrodynamic .

  * •

Study the stability of the collision between the hydro mode’s QNF and the
first non-hydro mode’s QNF. This would shed a new perspective on the validity
of the hydrodynamic approximation. Currently, the hydrodynamic expansion is
postulated to be valid up to the energy where the aforementioned collision
happens Withers:2018srf ; GrozdanovKovtun:2019 ; Grozdanov:2019 .
Consequently, finding instability would indicate that the limit of validity
itself is unstable. A similar comment applies to pole skipping points. These
are points at which the residue of a pole of a holographic Green’s function
vanishes Amado:2007yr ; Amado:2008ji . Such points are of particular interest
in relation to quantum chaos Grozdanov:2017ajz ; Blake:2017ris ; Blake:2018leo
.

  * •

Study the spectral stability at the phase transition in the holographic
superconductor model Hartnoll:2008vx . It would be interesting to see whether
the critical temperature is stable under the generic perturbations captured by
the pseudospectrum.

  * •

Study the spectral stability of AdS Reissner-Nordström black branes. The
pseudospectrum analysis of charged black branes might have important
consequences for the dual description of quantum critical phases of matter
zaanen2015holographic ; Hartnoll:2018xxg .

  * •

Study the spectral stability of other solutions frequently used in holography.
Beyond completing the current literature of QNFs in AdS spacetimes, it would
be interesting to explore the existence of some degree of universality in the
pseudospectrum. In particular, as spectra in AdS typically share the
"Christmas Tree" structure observed in the current work (see, for instance,
Kovtun:2005ev ; Amado:2009ts ), this study would answer whether the
overarching features observed here are fundamental to the structure of the
spectrum or depend heavily on the particular setup.

  * •

It has recently been suggested that quasinormal modes can serve as a signature
for having successfully simulated (quantum) gravity on a (quantum) computer
Biggs:2023sqw ; Maldacena:2023acv . In view of our results, this idea suffers
from similar problems as the idea of black hole spectroscopy from
gravitational wave signals. If such (quantum) simulations can be realized, the
task of identifying the quasinormal modes will be subject to similar
difficulties as for astrophysical black holes (see Jaramillo:2021tmt  for a
relevant discussion of this issue).

  * •

A different way of studying the poles of holographic Green’s functions is to
look for complex solutions for the momenta upon fixing the frequency to be a
real number Amado:2007pv . This is particularly important for the diffusive
mode and has implications for causality Landsteiner:2012gn  (see also the
recent preprint Gavassino:2023mad ). It would be interesting to study the
pseudospectra of these complex momentum modes.

## Acknowledgements

We thank J.L. Jaramillo and V. Boyanov for valuable discussions, and the
anonymous referee for the constructive comments which led to the inclusion of
section 3.3. The work of D.A and K.L. is supported through the grants
CEX2020-001007-S and PID2021-123017NB-100, PID2021-127726NB-I00 funded by
MCIN/AEI/10.13039/501100011033 and by ERDF “A way of making Europe”. The work
of D.G.F. is supported by JAEIntroICU-2022-IFT-02.

##  Appendix A Computation of $L^{\dagger}$

Here we present a detailed computation of the adjoint operator $L^{\dagger}$
for the real scalar field (equation (103)) and transverse gauge field
(equation (115)).

###  A.1 Real Scalar Field

For the real scalar field, the differential operator $L$ is given by

$L=i\begin{pmatrix}0&1\\\
L_{1}\left[\partial_{\mathscr{\rho}}^{2},\partial_{\mathscr{\rho}};\boldsymbol{\mathfrak{q}},\rho\right]&L_{2}\left[\partial_{\rho};\boldsymbol{\mathfrak{q}},\rho\right]\end{pmatrix}\,,$
(121)

with

$\displaystyle
L_{1}\left[\partial_{\rho}^{2},\partial_{\rho};\boldsymbol{\mathfrak{q}},\rho\right]$
$\displaystyle=\left[\mathscr{f}(\rho)-2\right]^{-1}\left[\frac{m^{2}l^{2}}{(1-\rho)^{2}}+\boldsymbol{\mathfrak{q}}^{2}-(1-\rho)^{3}\left(\frac{\mathscr{f}(\rho)}{(1-\rho)^{3}}\right)^{\prime}\partial_{\rho}-\mathscr{f}(\rho)\partial_{\rho}^{2}\right]\,,$
(122) $\displaystyle
L_{2}\left[\partial_{\rho};\boldsymbol{\mathfrak{q}},\rho\right]$
$\displaystyle=\left[\mathscr{f}(\rho)-2\right]^{-1}\left[(1-\rho)^{3}\left(\frac{\mathscr{f}(\rho)-1}{(1-\rho)^{3}}\right)^{\prime}+2\left(\mathscr{f}(\rho)-1\right)\partial_{\rho}\right]\,,$
(123)

and the inner product induced by the energy norm is:

$\displaystyle\expectationvalue{\phi_{1},\phi_{2}}_{E}=$
$\displaystyle\int_{0}^{1}\frac{d\rho}{(1-\rho)^{3}}\Biggl{[}\left(\frac{m^{2}l^{2}}{(1-\rho)^{2}}+\boldsymbol{\mathfrak{q}}^{2}\right)\overline{\phi_{1}}(\mathscr{t},\mathbf{k},\rho)\phi_{2}(\mathscr{t},\mathbf{k},\rho)+$
$\displaystyle+\mathscr{f}(\rho)\partial_{\rho}\overline{\phi_{1}}(\mathscr{t},\mathbf{k},\rho)\partial_{\rho}\phi_{2}(\mathscr{t},\mathbf{k},\rho)-\left(\mathscr{f}(\rho)-2\right)\overline{\psi_{1}}(\mathscr{t},\mathbf{k},\rho)\psi_{2}(\mathscr{t},\mathbf{k},\rho)\Biggr{]}\,.$
(124)

For notational convenience, we drop the explicit dependence of all functions
and denote $\rho$ derivatives with a prime.

Now, recalling the equality:

$\expectationvalue{L^{\dagger}\left[\phi_{1}\right],\phi_{2}}_{E}=\expectationvalue{\phi_{1},L\left[\phi_{2}\right]}_{E}\,,$
(125)

we can compute $L^{\dagger}$ integrating by parts

$\displaystyle\expectationvalue{\phi_{1},L\left[\phi_{2}\right]}_{E}=$
$\displaystyle\,i\int_{0}^{1}\frac{d\rho}{(1-\rho)^{3}}\Biggl{[}\left(\frac{m^{2}l^{2}}{(1-\rho)^{2}}+\boldsymbol{\mathfrak{q}}^{2}\right)\overline{\phi_{1}}\psi_{2}+\overline{\phi_{1}}^{\prime}\mathscr{f}\psi_{2}^{\prime}-\overline{\psi_{1}}\left\\{\left(\frac{m^{2}l^{2}}{(1-\rho)^{2}}+\boldsymbol{\mathfrak{q}}^{2}\right)\phi_{2}\right.$
$\displaystyle\left.-(1-\rho)^{3}\left(\frac{\mathscr{f}\phi_{2}^{\prime}}{(1-\rho)^{3}}\right)^{\prime}+(1-\rho)^{3}\left(\frac{(\mathscr{f}-1)\psi_{2}}{(1-\rho)^{3}}\right)^{\prime}+\left(\mathscr{f}-1\right)\psi_{2}^{\prime}\right\\}\Biggr{]}$
$\displaystyle=$
$\displaystyle-i\int_{0}^{1}\frac{d\rho}{(1-\rho)^{3}}\Biggl{[}\left(\frac{m^{2}l^{2}}{(1-\rho)^{2}}+\boldsymbol{\mathfrak{q}}^{2}\right)\phi_{2}\overline{\psi_{1}}+\phi_{2}^{\prime}\mathscr{f}\overline{\psi_{1}}^{\prime}-\psi_{2}\left\\{\left(\frac{m^{2}l^{2}}{(1-\rho)^{2}}+\boldsymbol{\mathfrak{q}}^{2}\right)\overline{\phi_{1}}\right.$
$\displaystyle\left.-(1-\rho)^{3}\left(\frac{\mathscr{f}\overline{\phi_{1}}^{\prime}}{(1-\rho)^{3}}\right)^{\prime}+(1-\rho)^{3}\left(\frac{(\mathscr{f}-1)\overline{\psi_{1}}}{(1-\rho)^{3}}\right)^{\prime}+\left(\mathscr{f}-1\right)\overline{\psi_{1}}^{\prime}\right\\}\Biggr{]}$
$\displaystyle+\left.i\frac{\mathscr{f}\overline{\phi_{1}}^{\prime}\psi_{2}}{(1-\rho)^{3}}\right|_{\rho=0}^{\rho=1}+\left.i\frac{\mathscr{f}\overline{\psi_{1}}\phi_{2}^{\prime}}{(1-\rho)^{3}}\right|_{\rho=0}^{\rho=1}-\left.2i\frac{(\mathscr{f}-1)\overline{\psi_{1}}\psi_{2}}{(1-\rho)^{3}}\right|_{\rho=0}^{\rho=1}$
$\displaystyle=$
$\displaystyle\expectationvalue{L\left[\phi_{1}\right],\phi_{2}}_{E}+\left.i\frac{\mathscr{f}\overline{\phi_{1}}^{\prime}\psi_{2}}{(1-\rho)^{3}}\right|_{\rho=0}^{\rho=1}+\left.i\frac{\mathscr{f}\overline{\psi_{1}}\phi_{2}^{\prime}}{(1-\rho)^{3}}\right|_{\rho=0}^{\rho=1}-\left.2i\frac{(\mathscr{f}-1)\overline{\psi_{1}}\psi_{2}}{(1-\rho)^{3}}\right|_{\rho=0}^{\rho=1}$
$\displaystyle=$
$\displaystyle\expectationvalue{L^{\dagger}\left[\phi_{1}\right],\phi_{2}}_{E}\,.$
(126)

The boundary terms vanish for $\rho=1$ as the chosen function space satisfies
the boundary condition (101). However at $\rho=0$, we have $\mathscr{f}(0)=0$,
and we get a non-zero contribution arising from the last term. Then, the final
expression for (A.1) is given by:

$\displaystyle\expectationvalue{\phi_{1},L\left[\phi_{2}\right]}_{E}=$
$\displaystyle\expectationvalue{L\left[\phi_{1}\right],\phi_{2}}_{E}-\left.2i\frac{(\mathscr{f}-1)\overline{\psi_{1}}\psi_{2}}{(1-\rho)^{3}}\right|_{\rho=0}$
$\displaystyle=$
$\displaystyle\expectationvalue{L\left[\phi_{1}\right],\phi_{2}}_{E}-i\int_{0}^{1}\frac{d\rho}{(1-\rho)^{3}}(\mathscr{f}-2)\left[2\delta(\rho)\frac{\mathscr{f}(\rho)-1}{\mathscr{f}(\rho)-2}\right]\overline{\psi_{1}}\psi_{2}$
$\displaystyle=$
$\displaystyle\expectationvalue{L\left[\phi_{1}\right],\phi_{2}}_{E}+\expectationvalue{\delta
L\left[\phi_{1}\right],\phi_{2}}_{E}=\expectationvalue{L^{\dagger}\left[\phi_{1}\right],\phi_{2}}_{E}\,,$
(127)

from where we recover the expression for the adjoint given in equation (103):

$L^{\dagger}=L+\delta L=L+\begin{pmatrix}0&0\\\\[6.0pt]
0&-2i\delta(\rho)\frac{\mathscr{f}(\rho)-1}{\mathscr{f}(\rho)-2}\end{pmatrix}\,.$
(128)

###  A.2 Transverse Gauge Field

For the transverse gauge field, the differential operator $L$ is given by

$L=i\begin{pmatrix}0&1\\\
L_{1}\left[\partial_{\mathscr{\rho}}^{2},\partial_{\mathscr{\rho}};\boldsymbol{\mathfrak{q}},\rho\right]&L_{2}\left[\partial_{\rho};\boldsymbol{\mathfrak{q}},\rho\right]\end{pmatrix}\,,$
(129)

with

$\displaystyle
L_{1}\left[\partial_{\rho}^{2},\partial_{\rho};\boldsymbol{\mathfrak{q}},\rho\right]$
$\displaystyle=\left[\mathscr{f}(\rho)-2\right]^{-1}\left[\boldsymbol{\mathfrak{q}}^{2}-(1-\rho)\left(\frac{\mathscr{f}(\rho)}{1-\rho}\right)^{\prime}\partial_{\rho}-\mathscr{f}(\rho)\partial_{\rho}^{2}\right]\,,$
(130) $\displaystyle
L_{2}\left[\partial_{\rho};\boldsymbol{\mathfrak{q}},\rho\right]$
$\displaystyle=\left[\mathscr{f}(\rho)-2\right]^{-1}\left[(1-\rho)\left(\frac{\mathscr{f}(\rho)-1}{(1-\rho)}\right)^{\prime}+2\left(\mathscr{f}(\rho)-1\right)\partial_{\rho}\right]\,,$
(131)

and the inner product induced by the energy norm is:

$\displaystyle\expectationvalue{a_{1},a_{2}}_{E}=$
$\displaystyle\int_{0}^{1}\frac{d\rho}{(1-\rho)}\Bigl{[}\boldsymbol{\mathfrak{q}}^{2}\overline{a_{1}}(\mathscr{t},\mathbf{k},\rho)a_{2}(\mathscr{t},\mathbf{k},\rho)+\mathscr{f}(\rho)\partial_{\rho}\overline{a_{1}}(\mathscr{t},\mathbf{k},\rho)\partial_{\rho}a_{2}(\mathscr{t},\mathbf{k},\rho)$
$\displaystyle-\left(\mathscr{f}(\rho)-2\right)\overline{\alpha_{1}}(\mathscr{t},\mathbf{k},\rho)\alpha_{2}(\mathscr{t},\mathbf{k},\rho)\Bigr{]}\,.$
(132)

As with the scalar, we drop the explicit dependence of all functions, denote
$\rho$ derivatives with a prime and compute $L^{\dagger}$ integrating by
parts:

$\displaystyle\expectationvalue{a_{1},L\left[a_{2}\right]}_{E}=$
$\displaystyle
i\int_{0}^{1}\frac{d\rho}{(1-\rho)^{3}}\Biggl{[}\boldsymbol{\mathfrak{q}}^{2}\overline{a_{1}}\alpha_{2}+\overline{a_{1}}^{\prime}\mathscr{f}\alpha_{2}^{\prime}-\overline{\alpha_{1}}\left\\{\boldsymbol{\mathfrak{q}}^{2}a_{2}-(1-\rho)\left(\frac{\mathscr{f}a_{2}^{\prime}}{(1-\rho)}\right)^{\prime}\right.$
$\displaystyle\left.+(1-\rho)\left(\frac{(\mathscr{f}-1)\alpha_{2}}{(1-\rho)}\right)^{\prime}+\left(\mathscr{f}-1\right)\alpha_{2}^{\prime}\right\\}\Biggr{]}$
$\displaystyle=$
$\displaystyle-i\int_{0}^{1}\frac{d\rho}{(1-\rho)^{3}}\Biggl{[}\boldsymbol{\mathfrak{q}}^{2}a_{2}\overline{\alpha_{1}}+a_{2}^{\prime}\mathscr{f}\alpha_{1}^{\prime}-\alpha_{2}\left\\{\boldsymbol{\mathfrak{q}}^{2}\overline{a_{1}}-(1-\rho)\left(\frac{\mathscr{f}\overline{a_{1}}^{\prime}}{(1-\rho)}\right)^{\prime}\right.$
$\displaystyle\left.+(1-\rho)\left(\frac{(\mathscr{f}-1)\overline{\alpha_{1}}}{(1-\rho)}\right)^{\prime}+\left(\mathscr{f}-1\right)\overline{\alpha_{1}}^{\prime}\right\\}\Biggr{]}$
$\displaystyle+\left.i\frac{\mathscr{f}\overline{a_{1}}^{\prime}\alpha_{2}}{(1-\rho)}\right|_{\rho=0}^{\rho=1}+\left.i\frac{\mathscr{f}\overline{\alpha_{1}}a_{2}^{\prime}}{(1-\rho)}\right|_{\rho=0}^{\rho=1}-\left.2i\frac{(\mathscr{f}-1)\overline{\alpha_{1}}\alpha_{2}}{(1-\rho)}\right|_{\rho=0}^{\rho=1}$
$\displaystyle=$
$\displaystyle\expectationvalue{L\left[a_{1}\right],a_{2}}_{E}+\left.i\frac{\mathscr{f}\overline{a_{1}}^{\prime}\alpha_{2}}{(1-\rho)}\right|_{\rho=0}^{\rho=1}+\left.i\frac{\mathscr{f}\overline{\alpha_{1}}a_{2}^{\prime}}{(1-\rho)}\right|_{\rho=0}^{\rho=1}-\left.2i\frac{(\mathscr{f}-1)\overline{\alpha_{1}}\alpha_{2}}{(1-\rho)}\right|_{\rho=0}^{\rho=1}$
$\displaystyle=$
$\displaystyle\expectationvalue{L^{\dagger}\left[a_{1}\right],a_{2}}_{E}\,.$
(133)

The boundary terms vanish for $\rho=1$ as the chosen function space satisfies
the boundary condition (113). Nonetheless, at $\rho=0$, we have
$\mathscr{f}(0)=0$, and we get a non-zero contribution arising from the last
term. Then, the final expression for (A.2) is given by:

$\displaystyle\expectationvalue{a_{1},L\left[a_{2}\right]}_{E}=$
$\displaystyle\expectationvalue{L\left[a_{1}\right],a_{2}}_{E}-i\int_{0}^{1}\frac{d\rho}{(1-\rho)}(\mathscr{f}-2)\left[2\delta(\rho)\frac{\mathscr{f}(\rho)-1}{\mathscr{f}(\rho)-2}\right]\overline{\alpha_{1}}\alpha_{2}$
$\displaystyle=$
$\displaystyle\expectationvalue{L\left[a_{1}\right],a_{2}}_{E}+\expectationvalue{\delta
L\left[a_{1}\right],a_{2}}_{E}=\expectationvalue{L^{\dagger}\left[a_{1}\right],a_{2}}_{E}\,,$
(134)

from where we recover the expression for the adjoint given in equation (115):

$L^{\dagger}=L+\delta L=L+\begin{pmatrix}0&0\\\\[6.0pt]
0&-2i\delta(\rho)\frac{\mathscr{f}(\rho)-1}{\mathscr{f}(\rho)-2}\end{pmatrix}\,.$
(135)

##  Appendix B Discretization in the Chebyshev grid

In this section we give some more details regarding the discretization
process. This discussion is based on the appendix C of Jaramillo:2020tuu .

###  B.1 Collocation Method

Working in the Chebyshev grid (116) with $N+1$ points is equivalent to
approximating a function $f(\rho)$ by a power series

$f(\rho)\approx\sum_{n=0}^{N}c_{n}T_{n}(2\rho-1)\,,$ (136)

where $T_{n}$ are the Chebyshev polynomials

$T_{n}(x)=\cos\left(n\arccos\left(x\right)\right)\,.$ (137)