Hunting for fermionic instabilities in charged AdS black holes

Fermions scattering on a black hole background cannot develop an instability sourced by superradiance. However, in a global (or planar) AdS$_4$-Reissner-Nordstr\"om background fermions can violate the AdS$_2$ fermionic Breitenl\"ohner-Freedman bound as measured by a near horizon observer at zero temperature. This suggests that AdS-Reissner-Nordstr\"om black holes might still be unstable to Dirac perturbations. Motivated by this observation we search for linear mode instabilities of Dirac fields in these backgrounds but find none. This is in contrast with the scalar field case, where a violation of the near-horizon BF bound in the AdS-Reissner-Nordstr\"om background triggers the already known scalar condensation near-horizon linear instability (in the infinite horizon radius, a.k.a. planar, limit this is Gubser's instability that initiated the holographic superconductor programme). We consider both the standard and alternative AdS/CFT quantizations. These are reflective boundary conditions that have vanishing energy flux at the asymptotic boundary.


Introduction
It is a well known fact that bosonic waves impinging charged or rotating black holes can be amplified via superradiant scattering (see e.g. the review [1] and references therein). It follows that black holes perturbed by bosonic fields in the presence of a gravitational potential well − provided by, for example, an asymptotic anti-de Sitter (AdS) potential, a physical cavity at finite radius or by the mass of the field − can develop a superradiant instability. However, this is not the only instability that can be present in such systems.
Indeed, the family of superradiant black holes always has a configuration with zero temperature. Typically, such extreme black holes have a near-horizon geometry that is the direct product (or a fibration) of a base space (e.g. a sphere) and an AdS 2 space [2]. These extreme (and near-extreme) black holes, when confined in a gravitational wall, are unstable if the 'effective mass' of the perturbation (as seen by an AdS 2 observer) violates the 2-dimensional Breitenlöhner-Freedman (BF) bound [3,4], i.e. a bound that is required to have finite energy perturbations. The superradiant and near-horizon instabilities have a different physical nature. But, in a non-extremal black hole they are usually entangled. However, if L is the typical dimension of the gravitational well (e.g. the radius of the cavity or the AdS radius), they disentangle for small dimensionless horizon radius r + /L since the near-horizon instability is suppressed for r + /L 1, while the superradiant instability is still present [5][6][7].
Not less well known is the fact that fermionic waves, unlike bosons, cannot suffer from superradiant scattering amplification [8][9][10][11][12]1]. Perhaps less familiar is the fact that fermions, like bosons [13][14][15][16], can also violate the 2-dimensional BF bound of the nearhorizon geometry of near-extremal black holes, if the fermion charge is high enough. In particular, this can happen in a Reissner-Nordström (RN) black hole in an asymptotically planar AdS background [15,[17][18][19][20][21][22][23] or in an asymptotically global AdS geometry (section 4 below). This suggests that such a RN might be unstable to condensation of a fermionic cloud around the horizon. Moreover, if the features from the bosonic field extend to the fermion case, then this might be a linear instability.
Motivated by these considerations, in the present manuscript we will search for linear mode instabilities of Dirac fields in a global AdS 4 RN black hole. We will not find any such instabilities. The absence of linear mode instabilities in the global AdS RN background is consistent with the fact that they are also not present in the planar AdS limit, r + /L → ∞, as found previously in [15,[17][18][19][20][21][22][23]. In view of these 'no-go' findings, in the conclusion remarks of section 5, we will argue that the absence of linear mode instabilities in the AdS RN background still leaves room (but not necessarily) for the following possibility: for a large number of fermions and in the semiclassical limit, the violation of the AdS 2 BF bound might signal a non-linear instability of the system.
In this manuscript, we also take the opportunity to improve/clarify the understanding of the near-horizon condensation instability of scalar fields. In particular, following a similar analysis for scalar fields in a Minkowsky cavity [7], we will explicitly show that the BF bound criterion for instability in AdS-RN is quantitatively sharp (section 2.4). This scalar field analysis will fit smoothly in our presentation since reviewing its details will also allow to make later a direct comparison with the Dirac field case and pinpoint major differences between (half-)integer spins. We further take the opportunity to show that the unstable modes belong to a family of near-extremal modes that are connected to the normal modes of AdS when the dimensionless horizon radius shrinks to zero. (This is an interesting observation because in de Sitter black holes the 'normal mode de Sitter family' is distinct from the 'near-extremal family of modes', where we are using the nomenclature of [24][25][26]).
Necessarily, we will also correct some statements/conclusions or clarify analyses/inter-pretations that were presented in previous literature. Namely, it has been claimed that the AdS/CFT correspondence no-source boundary conditions (i.e. those that do not deform the boundary conformal field theory) − often denoted as 'standard' or, if allowed, 'alternative' quantizations − do not yield Dirac field solutions with vanishing energy flux at the asymptotic boundary of AdS [27]. First, in sections 2.1-2.2, we will review the Dirac equation in a AdS-RN background and do the necessary field redefinitions in the (physical) Dirac field that allow it to separate and even decouple. Then, in section 2.3 we will analyse in detail the AdS/CFT standard and alternative quantizations of a Dirac field. In particular, we will show that the requirement that the source vanishes implies that the energy flux at the conformal boundary also vanishes (as it could not be otherwise). In this sense, these boundary conditions can also be denoted as reflective boundary conditions. For a massive Dirac field these no-source boundary conditions translate into Dirichlet or Neumann boundary conditions in the auxiliary separated/decoupled Dirac radial fields. However, for a massless fermion (a Weyl field) these no-source boundary conditions − which are Dirichlet or Neumann conditions on appropriate projections of the original physical Dirac field − translate into mixed (Robin) boundary conditions for the auxiliary separated/decoupled Dirac radial fields. This seems to have been the origin of the incorrect claim of [27]. The misleading focus on the boundary condition for the auxiliary field (and associated consequences) is similar to the one taken on the boundary conditions of the Regge-Wheeler−Zerilli master fields (aka Kodama-Ishibashi fields) in the case of gravitational perturbations of AdS black holes [28,29]. 1 In this context it is also important to clarify that in [31,32] massless Dirac quasinormal modes in Schwarschild-AdS were computed imposing Dirichlet boundary conditions in the auxiliary separated/decoupled fields. Unless otherwise justified, these boundary conditions are deployed of any known physical motivation. In particular, as discussed above, they are not the AdS/CFT boundary conditions and allow energy flux to leak through the asymptotic boundary.
Finally, in section 4.4 we will describe our strategy to search (unsuccessfully) for linear mode instabilities (sourced by the 2-dimensional BF bound violation) of Dirac fields in global AdS 4 -RN black holes. We consider both the standard and alternative quantizations and we will highlight the differences between the scalar and fermion systems. For a reader interested on a future detailed analysis of the frequency spectra, we also compute (analytically) the normal modes of massive and massless Dirac field in global AdS.

AdS-RN black holes and an orthogonal vierbein
The gravitational g µν and Maxwell A µ fields of the AdS-RN BH are described by 2 where dΩ 2 2 is the line element of a unit radius 2-sphere, M and Q are the mass and charge parameters. We will find convenient to replace M and Q by the event horizon radius r + (where f (r + ) = 0) and chemical potential µ. The relation between these two pairs of parameters is In (2.1), C is an arbitrary integration constant and fixing it amounts to choose a particular gauge. One common gauge choice is C = 0 where one has A| ∞ = 0 and the chemical potential is µ = −A| r + = Q/r + (we will typically use this one when presenting our results).
The components of any tensor in the coordinate basis {dx µ } can be obtained from the components on the tetrad basis using the projectors e (a) and e (a) . For example, T ν µ = 2 This is a solution of the Lagrangian L = √ −g R − 2Λ − 1 2 F 2 with F = dA. Note that if we rescale A → κA then the charge q of the perturbation field (to be discussed in later sections) rescales as q → κ −1 q so that qA, and thus the gauge covariant derivative, remain invariant.
σ T (a) (b) . This example also illustrates that we can have mixedindex tensors with mixed components in the coordinate and tetrad bases.
The spin connection of non-coordinate based differential geometry can be introduced in terms of the affine connection Γ β µν of coordinate based differential geometry as γ (c)(a)(b) = e µ (c) ∂ ν e (a)µ − Γ β µν e (a)β e ν (b) . Equivalently, one can define the spin connection as with λ (a)(b)(c) = e α (a) ∂ β e (b)α − ∂ α e (b)β e β (c) , which allows to compute the spin connections without the use of the affine connection.
For a multi-index tensor with tetrad and coordinate indices the mixed-index covariant derivative is defined as (2.7) Onwards, we take the affine connection Γ ν µα to be given by the Christoffel symbols that covariantly conserve the metric, ∇ α g µν = 0. It follows that the spin coefficients ω µ(a) (b) defined in (2.6) are such that the vierbein is also covariantly conserved, ∇ α e ν . Further note that the spin coefficient is anti-symmetric in the tetrad pair of indices, ω µ(a)(b) = −ω µ(b)(a) .

Dirac equation in the AdS-RN background
One of our main goals will be to compare scalar (spin-0) and fermionic (spin 1/2) perturbations on the AdS-RN background. Therefore, below we briefly review the equations that these fields have to obey in a curved background, e.g. in the AdS-RN spacetime (2.1). For more detailed discussions see [33][34][35][36].
To start, consider spin-s fields in Minkowski spacetime. The spin of a field can be identified looking into how the field transforms under a Lorentz transformation, x µ → x µ = Λ µ ν x ν . Let M αβ = −M βα be the generators of Lorentz transformations (i.e. a basis of six 4 × 4 antisymmetric matrices obeying the Lorentz Lie algebra). A finite Lorentz transformation is described by Λ = exp 1 2 Ω αβ M αβ where Ω αβ = −Ω βα are six parameters describing the particular transformation Λ (boost, rotations) of the Lorentz group SO (3,1).
Under  1 2 Ω αβ S αβ . 4 In Minkowski spacetime, under Lorentz transformations the derivative of a spin-s field transforms as ∂ µ Ψ → ∂ µΨ (x) = Λ ν µ S[Λ]∂ ν Ψ(Λ −1 x). If we want to couple the spin-s field to a curved background, while preserving general covariance, one needs to promote the partial derivative ∂ µ to a covariant derivative D µ . This promotion is chosen such that any function of Ψ and D µ Ψ that is a scalar under Lorentz transformations in Minkowski spacetime remains a scalar − under general coordinate transformations and local changes in the vierbein − in the curved background. This is the case if, under an arbitrary Lorentz transformation, the covariant derivative D µ still transforms as a derivative of a spin-s field: It follows that the covariant derivative of a spin-s field that preserves Lorentz invariance in a curved background is where we took the opportunity to allow the spin-s field to have a charge q (that couples to the Maxwell background field A µ ), and Γ µ is a covariant spin connection which depends on the spin of the field it acts on since the generators S (a)(b) depend on whether we are looking into, for example, the scalar or spinor representation of the Lorentz group. In more detail, for a scalar (spin-0) field Ψ ≡ Φ the Lorentz group generator is simply S (a)(b) = 0. Therefore S[Λ] = 1 and the scalar covariant derivative (2.13) is simply D µ = ∂ µ − iqA µ . The action for a massive charged complex scalar field is given by S Φ = M d 4 x √ −g D µ ΦD µ Φ * + m 2 ΦΦ * where * stands for complex conjugation. The factor of √ −g is introduced to ensure that the Lagrangian L Φ is a scalar density and thus the action S Φ = M d 4 xL Φ is a scalar. Varying this action w.r.t. Φ * one gets the Klein-Gordon equation for the scalar field D µ D µ Φ − m 2 Φ = 0 (2.15) and similarly for Φ * . On the other hand, for a (spin-1 2 ) Dirac 4-spinor field 5 Ψ ≡ ψ, out of the gamma matrices γ (a) (2.8) that satisfy the covariant Clifford algebra (2.9) one can build the commutator that satisfies the Lorentz Lie algebra. 6 S (a)(b) is the generator of the Lorentz group in the spinor representation and replacing this (2.16) into (2.14) and then the latter Γ µ into (2.12) one gets the spinor covariant derivative D µ , namely (2.13), that acts on the Dirac spinor ψ.
The action that is Lorentz invariant and describes the coupling of a spin-1 2 fermion field ψ to a curved background M is where we have introduced the Dirac adjointψ = ψ † γ (0) with ψ † = (ψ * ) T being the Hermitian adjoint of the multi-component field ψ. One needs to work with the Dirac adjoint because the Fermi bilinearsψψ andψγ µ ψ transform covariantly (as a scalar and as a vector, respectively) under the Lorentz group (while the Hermitian partner objects do not).
Varying the action S D w.r.t.ψ and ψ, respectively, one gets the Dirac equations The Dirac spinor is a 4-component field ψ A with complex components A = 1, 2, 3, 4. In our study we will typically omit the spinorial indices and simply write ψ ≡ ψ A , γ (a) ≡ (γ (a) ) A B and S[Λ] ≡ S[Λ] A B . 6 More precisely, Dirac 4-spinor fields are invariant under internal local Lorentz transformations of the spinor representation of the SU (4) group. There are 15 Dirac matrices that provide a 4 × 4 fundamental representation of the SU(4) group [36]. These are the four vectors γ (a) introduced in (2.8), the six tensors S (a)(b) defined (2.16), the pseudoscalar γ (5) introduced in (2.10), and four axial vectors γ (5) γ (a) .
To find solutions of (2.18) it is advantageous to write the Dirac 4-spinor ψ in terms of the left-handed and right-handed 2-spinors Ψ − and Ψ + , respectively, as The chiral 2-spinors Ψ ± emerge naturally when we note that the pseudoscalar γ (5) defined in (2.10) obeys (γ (5) ) 2 = 1. Therefore we can introduce the Lorentz invariant projection operators P ± that project the Dirac 4-spinor ψ into the chiral spinors: and such that P 2 ± = P ± and P + P − = 0. 7 Moreover, the task of finding solutions of the Dirac equation in the AdS-RN black hole gets considerably simplified by the fact that under the separation anstaz [37,38]: the Dirac equations (2.18) reduce to a set of equations where the radial and angular functions of the fermion field are decoupled. This separation ansatz exploits the fact that ∂ t and ∂ φ are Killing vectors of the background AdS-RN solution. This allows to do a Fourier decomposition in these directions which introduces the frequency ω and azimuthal angular momentum m φ of the fermionic wavefunction. 8 Concretely, the radial functions R 1 (r), R 2 (r) obey the coupled system of first order ODEs where λ is a separation constant, while the angular functions S 1 (θ), S 2 (θ) satisfy the coupled system of first order ODEs  8 This frequency ω is measured in the gauge A|∞ = 0 − see (2.1) − and we will work preferentially on this gauge unless otherwise stated. In particular, all our numerical results will be using it. Note that in the alternative gauge A|r + = 0 (also often used) the associated frequency isω = ω − qµ.
Furthermore, the coupled pair of first order radial equations (2.22) can be decoupled in a pair of second order ODEs, one for R 1 (r) and the other for R 2 (r). For that we solve the first (second) equation in (2.22) w.r.t. R 2 (R 1 ) and replace it in the second (first) equation. We end up with two decoupled second order ODEs for R 1 and R 2 , where * denotes complex conjugation and we have defined Of course, we are only interested in solutions of (2.24) that also solve the original first order system (2.22). The requirement that (2.22) is solved imposes extra constraints on solutions of (2.24). This is best illustrated if we consider the Taylor expansion about the boundaries of the integration domain: the ODE pair (2.24) has four integration constants about each boundary but only two of them are independent when we further require that the solution solves the two first order ODEs (2.22); see discussion of (2.27) below. Similarly, the coupled pair of first order ODEs for S 1,2 (θ) can be written as a decoupled set of two second order ODEs for S 1 (θ) and S 2 (θ). They are hypergeometric equations and S 1,2 (θ) are the spin-1 2 weighted spherical harmonics. Regularity at θ = 0 and θ = π quantizes the angular separation constant as ( is a harmonic number related to the number of zeros of the wavefunction) with the azimuthal number being constrained as m φ ≤ . Unfortunately, the radial ODEs cannot be solved analytically 9 . We can however do a Frobenius analysis about the asymptotic boundary r → ∞ to find the asymptotic behaviours of R 1 (r) and R 2 (r). One finds that (for m = 0, 1 2 ) 10 10 For m = 1/2 one of the two independent solutions decays asymptotically as a power law in r and the other as a power law multiplied by a log r. We do not discuss further this special case (see [17] for more details).
where we used (−gf ) − 1 4 | r→∞ ∼ L 1/2 r −3/2 and, anticipating the AdS/CFT discussion below, we have introduced the conformal dimensions As expected for a coupled system (2.22) of two first order ODEs, there are two independent arbitrary constants (α 1 , β 1 ) in the asymptotic decay (2.27), that is to say, the decays of R 2 are fixed by the equations of motion as a function of (α 1 , β 1 ). The dots in (2.27) represent subleading terms that depend only on α 1 (in the ∆ − contribution) or β 1 (in the ∆ + terms). Note that if we take m → −m, we simply trade the role of the ∆ ± contributions. Therefore onwards we assume m ≥ 0 in our discussion. The m = 0 case is special since ∆ + = ∆ − (in this sense we can call this the fermionic BF bound). Choosing the gauge A| ∞ = 0, a Frobenius analysis of the first order equations of motion about the asymptotic boundary yields for this case 11 i.e. we can take the two independent integration constants associated to the coupled pair of first order ODEs to be α 1 and β 1 and the equations of motion then fix the decay of R 2 as a function of α 1 and β 1 .

Boundary conditions for the Dirac spinor in AdS-RN
To find the solution of the Dirac spinor field ψ and its Dirac adjointψ in the AdS-RN background we have to solve a system of two equations that are first order, namely (2.22), subject to boundary conditions imposed at the event horizon r = r + and at the asymptotic boundary r → ∞. Before imposing boundary conditions, such a system of two first order differential equations necessarily has two independent constants at the horizon boundary and another two independent constants at the asymptotic boundary (namely, α 1 and β 1 in (2.27)), which can be identified doing a Frobenius analysis at these two boundaries. To have a well posed formulation of the elliptic problem one should impose two boundary conditions that fix two of the independent constants and solve the equations of motion to find the other two. We certainly want the Dirac solutions to be regular at the event horizon: this 11 For a Dirac field (or scalar field) with phase ϕ, ψ = |ψ|e iϕ , U (1) gauge transformations with gauge parameter χ leave the action and equations of motion invariant and transform the Dirac (scalar) and Maxwell fields as ϕ →φ = ϕ + q χ , At →Ãt = At + ∇tχ. Thus, if in the gauge A|∞ = 0 (i.e. C = 0) we denote the frequency of the Dirac (scalar) field by ω then a transformation with gauge parameter χ = µ t into the gaugeÃ|∞ = µ (C = µ) changes the frequency intoω = ω − qµ. boundary condition fixes one of the constants 12 . One should then fix one of the asymptotic constants α 1 or β 1 (or a relation between them) with an appropriate boundary condition. But we certainly cannot fix both asymptotic independent constants: once the first is fixed, the second one must be found by solving the equations of motion in the bulk subject to the two aforementioned boundary conditions. 13 This poses the question: how do we choose a boundary condition at the asymptotic boundary that is physically relevant? Dirac spinor fields ψ are intrinsically quantum fields. The dynamics of such fields can be naturally described by a path integral formulation whereby one sums over all possible field configurations in configuration space to get the transition amplitude between two states. In particular, the partition function Z (i.e. the generating functional of correlation functions between operators) can also be naturally computed using the path integral formulation. Schematically one has, S cl [ψ,ψ] is the action evaluated on a solution of the classical equations of motion, that follow from the variation δS = 0 subject to the boundary conditions. As emphasised in [39][40][41][42]19] this statement that the action must be stationary when evaluated on a classical solution severely constrains the type of boundary conditions we can impose on the field ψ. Indeed, if δS = 0 then it is not necessarily true that δ(S + B) = 0 where B is the boundary term describing the desired boundary conditions (i.e. a total derivative term that does not change the equations of motion). That is to say, the physical choice we make for the boundary conditions must be such that δB = 0. In particular, in the context of the AdS/CFT correspondence, this condition fixes the form of the boundary term that must be added to the standard Dirac action (2.17) to have stationary solutions. Vice-versa, this boundary term fixes the boundary field theory. So which boundary term B should be added to the (2.17) to ensure that S is stationary in the class of histories contributing to the variational principle? To answer this question one first notes that the "radial" Dirac gamma matrix γ (1) defined in (2.8) satisfies (γ (1) ) 2 = 12 Alternatively, since we have a ODE system, we could use two boundary conditions to fix the two asymptotic independent constants and solving the equations of motion would yield the behaviour of the Dirac fields at the event horizon. However, this is not a good strategy because in general these solutions would not be smooth at the event horizon. 13 Similarly, for the hyperbolic problem we just need to give ψ (which automatically fixes ψ † ) as initial data − but not the time derivative of the 4-spinor − to determine the full time evolution. Again, this is a consequence of the fact that the Dirac equations are first order which implies that the momentum conjugate to the 4-spinor ψ is not its time derivative but simply ψ † . It is also this property that is ultimately responsible for the fact that the number of degrees of freedom per spatial point of the system is 4 (2 spins, up and down, for the fermion and another two for the anti-fermion). Indeed, the 4 spinor ψ − and thus ψ † − have 4 complex components and thus 8 real components. The number of degrees of freedom is half the dimension of the phase space, i.e. 4. If the Dirac equations were to be second order, then the momentum would be proportional to the time derivative of ψ and thus contribute with 4 extra degrees of freedom. I 4 and γ (1) = γ (1) † . It follows that we can decompose the Dirac spinor as where ψ ± (ψ ± ) are 4-eigenspinors of γ (1) with eigenvalue ±1 (∓1). 14 Using this property, including the associated properties listed in footnote 14, one finds that the terms in the Dirac action (2.17) that contain radial derivatives of the spinor are where we used the fact that D r = f 1/2 ∂ r . It follows that if we vary the Dirac action (2.17) w.r.t. ψ + and ψ − one gets, after integration by parts, where the bulk terms describe a contribution that vanishes when the equations of motion − which are equivalent to (2.18) − are satisfied and δS bdry is a boundary term resulting from integrating by parts the radial derivative terms (2.32) given by As discussed above, to have a well-posed boundary value problem, after requiring that the solution is regular at the event horizon we no longer have the freedom to fix both ψ + and ψ − at the asymptotic boundary (these are the two independent asymptotic constants of our pair of first order ODEs). Instead, we can either fix ψ + at the asymptotic boundary (in which case δψ + = 0) or fix the asymptotic value of ψ − (in which case δψ − = 0 at the boundary).
Suppose we want to fix ψ + at the asymptotic boundary (a similar analysis would apply if we wanted to fix ψ − ). In order to have a well-defined variational problem one should add a boundary term that cancels the contributionψ + δψ − in (2.34). Adding the boundary term [39][40][41][42]19] produces the desired effect since the total on-shell action becomes which indeed vanishes when δψ + = 0 (and thus δψ + = 0). We can also compute the momentum conjugate to ψ + andψ + by varying S tot w.r.t. ψ + andψ + , respectively, yielding In more detail, γ (1) ψ± = ±ψ± and γ (1)ψ ± = ∓ψ± and thus ψ± = 1 2 I4 ± γ (1) ψ andψ± = 1 2ψ I4 ∓ γ (1) . A few properties follow that are useful. For example,ψγ (1) In terms of the functions R 1 (r), R 2 (r) and S 1 (θ), S 2 (θ) introduced in the separation ansatz (2.21), the 4-spinors ψ ± are given by From the asymptotic decays of R 1,2 in (2.27) (valid for m = 0, 1 2 ) or in (2.29) (valid for m = 0) one finds that ψ ± decay as if 0 < mL < 1 2 ; (2.39) where α 1 , β 1 are the free constants introduced in (2.27) or (2.29) and the constants a(α 1 ), b(β 1 ),ã(α 1 ) andb(α 1 ) are fixed as functions of α 1 or β 1 (as described by their argument) by the equations of motion (details are irrelevant for our aim). The asymptotic decays of the Dirac adjointsψ ± follow straightforwardly from (2.39) with the exchange α 1 →ᾱ 1 , β 1 →β 1 , etc. For mL > 1 2 the only normalizable mode (i.e. with finite energy) is ψ + [43,15,17,21]. In the context of the AdS/CFT correspondence, the leading term of the asymptotic expansion lim r→∞ r ∆ − ψ + = 2α 1 is then identified with the source of the dual operator O which has mass dimension ∆ + . We have a well-posed boundary value problem if we impose smoothness of ψ + at the event horizon and a Dirichlet boundary condition for α 1 at the asymptotic boundary. In particular, if we do not want to deform the boundary field theory we impose the no-source Dirichlet boundary condition: α 1 = 0. We have no freedom left to fix asymptotically ψ − i.e. β − . Instead, β − and thus ψ − | ∞ is determined by solving the Dirac equations subject to the above boundary conditions. The expectation value Ō of the dual operator is given by the conjugate momentum Π + defined in (2.37): On the other hand for 0 ≤ mL < 1 2 both modes ψ ± are normalizable [43,15,17,21]. Thus we can still impose the standard quantization where we identify the lim r→∞ r ∆ − ψ + ≡ ψ (0) + as the source of the dual operatorŌ. In particular the no-source standard boundary condition for all possible masses: As mentioned above, for this range of masses both modes are normalisable, therefore we can also impose the so-called alternative quantization; where we identify the lim r→∞ r ∆ + ψ − ≡ ψ (0) − as the source of the dual operator O with mass dimension ∆ − . In particular, if we do not want to deform the boundary field theory we impose the no-source alternative boundary condition: The two quantizations (2.42) and (2.43) yield two distinct boundary conformal field theories [43,15,17]. For m = 0, note that the Dirichlet boundary condition on ψ + , ψ (0) + = 0, implies the Neumann condition in ψ − (i.e. the next-to-leading order term in the expansion for ψ − vanishes) and vice-versa. This follows straightforwardly from an inspection of (2.40).
It is important to emphasize that the no-source standard and alternative boundary conditions (2.42)-(2.43) that do not deform the boundary theory imply that the energy flux and fermion particle flux vanish at the asymptotic boundary . In this sense we can regard these as 'reflective' boundary conditions. The Dirac action (2.17) (and (2.36)) is left invariant if we rotate the phase of the Dirac spinor, ψ → e −iα ψ. The Dirac current associated to this symmetry is j µ =ψγ µ ψ and one can check that it is conserved, ∇ µ j µ = 0 after using the first order equations of motion (2.18). This is an internal vector symmetry since ψ ± transform in the same way under this symmetry. This current gives the charge flux or particle number flux of fermions. The associated conserved charge is γ ab is the associated induced metric, and ξ = ∂ t is the Killing vector describing time translations. In particular, j r | r→∞ gives the radial flux of particles at the asymptotic spacelike boundary Σ. One can also show that the energy flux across a spacelike boundary is proportional to the Dirac current. In particular, the energy flux across the asymptotic boundary Φ ∂t | ∞ is proportional to the particle flux j r | ∞ and is given by 15 (2.44) Inserting the asymptotic decays (2.27) for R 1,2 this yields That is, the energy flux at the asymptotic boundary vanishes if we impose the above discussed no-source Dirichlet boundary conditions α 1 = 0 or, for the alternative quantization, β 1 = 0 which do not deform the boundary conformal field theory. 15 Let again ξ = ∂t be the Killing vector field conjugate to the energy. The energy-momentum tensor for the Dirac field is Tµν = i 2 ψ γ (µ D ν) ψ − D (µψ γ ν) ψ and it is conserved ∇µT µν = 0. This conservation law together with the Killing equation, ∇ (µ ξ ν) = 0, imply that the 1-form Jµ = −Tµν ξ ν is conserved, d J = 0, where is the Hodge dual. We can then define the energy flux across the asymptotic hypersurface Σ (like the asymptotic boundary) as Φ ξ ≡ − Σ J = − Σ dVΣ Tµν ξ µ n ν where n ν is the unit normal vector to Σ and dVΣ is the induced volume on Σ.

Near-horizon geometry of the extreme AdS-RN black hole
The near-horizon geometry of the extremal AdS-RN black hole will play an important role in our discussions in sections 3 and 4. Therefore, we review it here. The limiting procedure described below was first presented in [2].
The extremal AdS-RN black hole is given by (2.1) with µ = µ ext given by (2.4). To obtain the near-horizon geometry, it is convenient to work in the gauge A| r + = 0 (C = µ; otherwise we can do a gauge transformation in the end). One first zooms in around the horizon region by making the coordinate transformations: where L AdS 2 is the AdS 2 radius (to be defined below). Now the near-horizon geometry is obtained by taking ε → 0 which yields solves the 2-dimensional Einstein-AdS equations, R µν = −L −2 AdS 2 g µν , if L AdS 2 is identified as a function of the AdS 4 radius L and the horizon radius r + as indicated in the first line of (2.48).

Scalar fields in a AdS-RN background and their instabilities
Scalar fields confined inside the gravitational potential (like the AdS potential or a box in an asymptotically Minkowski background) of a black hole can condense creating nearhorizon linear instabilities [13-16, 44, 5, 6] (for planar AdS, this instability triggered the holographic superconductor programme [13][14][15]). Essentially this happens because we can have scalar fields that obey the asymptotically AdS 4 UV Breitenlöhner-Freedman (BF) bound but violate the 2-dimensional BF bound associated to the AdS 2 × S 2 near-horizon geometry of the extremal black hole of the system. As we shall discuss in section 4.1, a similar violation of the 2-dimensional BF bound can occur for Dirac fields. In spite of this, as we will find in section 4.4, it turns out that Dirac fields are not linearly unstable to the near-horizon condensation mechanism. Therefore, before we discuss the fermionic case, it is important to revisit the scalar field case. This will allow to: 1) motivate the search of linear instabilities due to Dirac fields done in this manuscript, 2) eventually identify differences between the two spins that could help in understanding the opposite outcomes. We also take the opportunity to demonstrate: i) how remarkably sharp the near-horizon instability bound (3.7) is by comparing it with the numerical solutions of the Klein-Gordon equation, and ii) that the unstable modes are both peaked near the horizon but also connected to the AdS normal modes (that is to say, in the language of [24][25][26] the AdS and near-extremal families of modes coincide and describe the unstable modes).
Using the fact that the AdS-RN background (2.1) is static and spherically symmetric we can consider a separation ansatz for the scalar field (with mass m and charge q) with the Fourier decomposition where Y (θ) are the familiar (spin-0) spherical harmonics − which are regular when the separation constant of the system is quantized as λ = ( +1), = 0, 1, 2, · · · − and |m φ | ≤ is the azimuthal quantum number. The Klein-Gordon equation then yields the following equation for the radial function φ(r): A Taylor expansion around the asymptotic boundary yields the two independent solutions being the conformal dimensions of the field. Such a scalar field in AdS 4 is normalizable (i.e. it has finite energy) as long as its mass obeys the AdS 4 Breitenlöhner and Freedman (BF) bound, m 2 ≥ m 2 BF ≡ − 9 Such a scalar field that is stable in the UV region can however be unstable in the IR region. This is best understood if we take the near-horizon limit of (3.2). Concretely, applying the near-horizon coordinate transformation (2.47) together with the near-horizon frequency transformationω → ω ε/L 2 AdS 2 (so that e −iωt → e −i ωτ ) followed by the nearhorizon limit ε → 0 yields the radial Klein-Gordon equation in the near-horizon geometry (2.48): This is nothing else but the Klein-Gordon equation for a scalar field around AdS 2 with an electromagnetic potential A τ = α ρ. A Frobenius analysis of (3.4) yields 5) which determines the effective mass of the scalar field from the perspective of a near-horizon observer, Now, a scalar field with mass (3.6) in AdS 2 has unstable modes if it violates the . It follows that extremal AdS-RN 4 black holes should be unstable whenever the charge of the scalar field obeys (near-horizon instability bound) .
(3.7) Note that scalar fields can also induce instabilities due to another mechanism that is known as supperradiance. Unlike the near-horizon instability − which is suppressed in the limit r + /L → 0; indeed (3.7) goes as q 2 L 2 ≥ L 2 8r 2 + + O(1) − the superradiant instability is present for small r + /L 1 black holes. For example, for m = 0, from the perturbative results of [6] one finds that the superradiant instability in extremal AdS-RN 4 is present for scalar charges 17 Next, we solve the Klein-Gordon equation numerically to confirm that the near-horizon and superradiant instabilities are indeed present and to find how sharp the instability bounds (3.7) and (3.8) are. We present results for scalar masses above the unitarity bound m 2 BF + 1 = −5/4 so asymptotically we impose the Dirichlet boundary condition a = 0; see (3.3). 18 On the other hand at the horizon we require that the solution is regular in the future 17 This bound can be obtained from the expression for the frequency obtained in section III.D of [6].
Namely, the onset charge (3.8) of the superradiant instability is obtained by settingω = 0 and µ = µext in equation (55) of [6] and solving for the charge q. For further discussions between the entanglement of the superradiant and near-horizon instabilities and their different nature we ask the reader to see [6] and [7]. 18 For m 2 BF < m 2 < m 2 BF + 1 both modes are normalizable and thus we could also impose the Neumann boundary condition b = 0 (the so-called alternative quantization in the context of AdS/CFT). horizon which discards outgoing modes. To present the results, it is important to note that our system has a scaling symmetry [6] which means that the physical dimensionless quantities that are relevant for the problem are (this effectively sets L ≡ 1) First, we are interested on finding the onset of the instabilities namely, the scalar field charge q onset (µ, r + /L; mL, ) above which the system is unstable. This onset occurs when the frequency satisfiesω = ω − qµ = 0. The Klein-Gordon equation (3.2) is then solved as an eigenvalue problem for q = q onset . For concreteness, we fix = 1 (we need m φ ≥ 1 to have an instability). In the left plot of Fig. 1 we set m = 0 and we plot the dimensionless onset charge q onset L as a function of the dimensionless horizon radius r + /L for different values of the chemical potential µ = µ ext (1−10 −x ) that increasingly approaches the extremal value. From top to the bottom the green numerical curves describe chemical potentials with x = 2, 3, 6, 15. We see that as we get closer to extremality these onset curves increasingly approach (for values of r + /L larger than ∼ 0.25) the red dashed curve which describes the near-horizon bound (3.7). This strongly suggests that the instability, for large values of the horizon radius and near extremality, can be understood as due to the violation of the AdS 2 BF and that the associated near-horizon bound (3.7) is sharp (i.e. it is attained at extremality). On the other hand, as pointed out before, the near-horizon red dashed curve diverges as r + /L → 0. However, Fig. 1 shows that q onset L is finite for small r + /L. Actually, in this regime the numerical onset curves are well described by the superradiant bound (3.8) (blue dashed curve with negative slope). This suggests that for Figure 2. Scalar field frequency as a function of the dimensionless scalar charge qL for a AdS-RN black hole with µ = 0.99µ ext and r + /L = 0.5 (mL = 0 and = 1). Left panel: Imaginary part of the dimensionless frequency, Im(ωL). The system becomes unstable for q > q where q L ∼ 1.863. Right panel: Real part of the dimensionless frequency, Re(ωL), measured with respect to qµL. This quantity changes sign at q = q , i.e. when Im(ωL) changes sign. In both plots, the black dashed curves describe the analytic prediction of the asymptotic matching expansion (A.10). We find very good agreement with the numerical results (blue curves) for qL < 1.1 (say). This is a further justification of the crude assumption that we should match with 0 in the overlapping region.
However we find that these modes connect with the AdS normal modes as r + → 0.
small horizon radius and near extremality the instability has a superradiant nature and the superradiant bound (3.8) is sharp. For finite values of r + /L, i.e. away from the r + /L → 0 and r + /L → ∞ regions, the superradiant and near-horizon instabilities are entangled. These features are not unique to the massless case. For example, the onset charge plot for a scalar mass of mL = 2 is shown in the right panel of Fig. 1. Again, as extremality is approached the numerical green curves increasingly approach the near-horizon onset bound (3.7) (in this case we do not show the curve corresponding to the perturbative superradiant curve because it was not computed in [6] but we see that the behaviour of the onset curves for r + /L 1 is similar to the massless case).
To compare with what happens in the Dirac field case, it is enlightening to do the following exercise whose results are summarized in Fig. 2. We pick a particular AdS-RN background with chemical potential µ = 0.99µ ext and horizon radius r + /L = 0.5. We also fix the scalar mass to be mL = 0 and the scalar field harmonic number = 1. Then we solve the Klein-Gordon equation to find the imaginary and real parts of the frequency ωL as a function of the dimensionless scalar field charge qL: these are shown in the left and right panels, respectively, of Fig. 2. From the left panel we see that, in accordance with the conclusions of Fig. 1, for small qL the system is stable (since Im(ωL) < 0) but there is a critical charge q L ∼ 1.863 (vertical dotted line) above which the system becomes unstable. Precisely at this critical onset charge one has Re(ω) − qµ = 0 and this quantity is negative (positive) for q < q (q > q ). The inset plot of Fig. 1 zooms-in the region q < q . In section 4.4 we will find that the partner plot for Dirac fields is substantially different.
We take also the opportunity to understand better the frequency spectrum of scalar fields in AdS-RN. For global AdS RN black holes there are two quasinormal mode families [45][46][47]: one whose imaginary part grows negative without bound as the horizon radius r + /L decreases, and another whose imaginary part vanishes as r + /L → 0 and whose real part approaches the normal modes of AdS. The unstable modes are found in this second family. This could well be the complete story. However, in de Sitter black holes there is a third family of quasinormal modes − called the near-extremal family − whose wavefunctions are spatially peaked near the horizon and that is distinct from the de Sitter family (as the name suggests, the latter is connected to the normal modes of de Sitter when the black hole shrinks). This naturally raises the question: could it be that in AdS one also has a near-extremal family of quasinormal modes that is not connected to the AdS family? If so, do the near-horizon unstable modes with bound (3.7) fit in this near-extremal family? We find a negative answer to these questions: the unstable modes belong to the AdS family of modes and the near-extremal family coincides with the AdS family. To arrive to this conclusion we first use a matching asymptotic expansion similar to the one used in de Sitter [24][25][26]48] to find the frequency spectrum of the near-extremal family of quasinormal modes. This is done in Appendix A and here we just quote the final result. Namely, near-extremality and for small scalar field charge one finds that near-extremal modes have the frequency (for the lowest radial overtone p = 0) where R + = r + /L, e = qL and σ = r + −r − r + measures the distance away from extremality with r − (r + , µ, L) being the inner (Cauchy) horizon for which f (r − ) = 0. In Fig. 2, this analytical near-extremal frequency (with = 1) is described by the dashed black curve. We find that it matches quite well the numerical result for small scalar charge. This indicates that the unstable modes fit into the near-extremal family of modes. But they also fit into the AdS family of normal modes. That is to say, unlike in the de Sitter case, in AdS the near-extremal and AdS family of modes coincide. To see this is indeed the case we pick two solutions in Fig. 2 that have qL = 2.5 (orange diamond) and qL = 2 and (keeping µ = 0.99µ ext fixed) we follow this family of unstable modes as r + /L decreases to zero. 19 This is done in Fig. 3 for qL = 2.5 and Fig. 4 for qL = 2. In both cases we find that, as r + /L → 0, Im(ωL) → 0 and Re(ωL) → 3, which is indeed the normal mode frequency of AdS with = 1 (and lowest radial overtone). 19 Note that qL = 2.5 is well above the onset curves of Fig. 1 for any r+/L while the qL = 2 line is above the onset curves only above a certain horizon radius. So, for the latter charge, the system is unstable only above a critical value of r+/L, as shown in Fig. 4. Figure 3. Scalar field frequency as a function of the dimensionless horizon radius r + /L for a AdS-RN black hole with µ = 0.99µ ext and qL = 2.5 that is always above the near-horizon bound (3.7) (mL = 0 and = 1). Left panel: Imaginary part of the dimensionless frequency, Im(ωL). Right panel: Real part of the dimensionless frequency, Re(ωL). In both plots, the magenta dashed curves describe the analytic prediction of the perturbative expansion in r + /L about AdS. The unstable modes are thus connected to the AdS normal modes when r + → 0 (brown disk). For reference, the blue disk with r + /L = 0.5 is shown (which makes contact with Fig. 2). In Fig. 3 and Fig. 4 the magenta dashed lines departing from the normal mode of AdS describe the frequency that one obtains when we consider a perturbative expansion in r + /L and near-extremality about global AdS (and = 1, m = 0). This result is taken from [6] (we already mentioned it to get the bound (3.8)): where e = qL. So we see that not only the unstable modes approach the normal modes of AdS but they also do it at the expected rate in an expansion in r + /L. That is to say, the matching of our numerical results with the perturbative results (3.10) and (3.11) represents a non-trivial check of our results and illustrates the regime of validity of the perturbative results.
Now that we have highlighted the key features of the near-horizon (and superradiant) instabilities due to scalar perturbations in AdS-RN, we can proceed to the study of perturbations of Dirac fields in AdS-RN.

Searching for an instability of Dirac fields in the AdS-RN background
In section 3 we have seen that scalar fluctuations in the AdS-RN background give rise to the near-horizon scalar condensation instability. Moreover, we have seen that this instability is closely associated to the violation of the AdS 2 scalar BF bound. So much that the associated bound (3.7) for the onset of the instability is sharp. This naturally invites the questions: in the fermionic case can we also have a range of parameters where the AdS 2 fermionic BF bound is violated? If so what is the equivalent bound to (3.7) for the onset of the instability?
In this section we will address these questions. We will find that a near-horizon analysis of the Dirac equation indeed indicates that the AdS 2 fermionic BF bound can be violated near-extremality if the charge of the fermion is above a critical value (subsection 4.1). Encouraged by this result we will do a numerical analysis that will search for unstable modes in the region of parameters of interest (subsection 4.3). However, we will find no trace of instabilities, unlike in the scalar field case.

Argument for a near-horizon instability of Dirac fields
The Dirac equation in the near-horizon geometry (2.48) of the extreme AdS-RN black hole can be obtained taking the near-horizon limit of section 2.4 directly on the Dirac equation (2.24) for the extreme AdS-RN black hole. Concretely, applying the near-horizon coordinate transformation (2.47) together with the near-horizon frequency transformatioñ ω → ω ε/L 2 AdS 2 (so that e −iωt → e −i ωτ ) followed by the near-horizon limit ε → 0 yields the Dirac equation in the near-horizon geometry (2.48): 20 where the AdS 2 radius L AdS 2 and the Maxwell near-horizon parameter α are defined in (2.48) and λ is the angular eigenvalue quantized as in (2.26). Also, recall that m ≥ 0 and q ≥ 0 are the mass and charge of the fermionic field. Figure 5. Predicted Dirac field charge for the onset of an instability as a function of the horizon radius for an extremal AdS-RN black hole (µ = µ ext ). In both plots the red dashed curve is the lower bound of (4.4). The plot in the left (right) panel is for fermion mass mL = 0 (mL = 4) and harmonic number = 1/2. The near-horizon analysis of the 2-dimensional BF bound violation leading to (4.4) predicts that region B should be unstable while region A should be stable (at least with respect to the BF bound mechanism). Note that for small r + /L the system is not unstable because there is no superradiance for fermions and the predicted near-horizon instability is also suppressed. These Dirac figures can be (qualitatively) compared with Fig. 1 for the scalar field.
Asymptotically, as ρ → ∞, a Frobenius analysis of (4.1) finds that the solution R 1 (ρ) decays as where α 1 , β 1 are two arbitrary constants and we have introduced the AdS 2 conformal dimensions The 2-dimensional fermionic BF bound is obtained equating ∆ + = ∆ − : thus it is obeyed if m 2 eff ≥ 0, i.e. ∆ ± must be real. It follows that we can have situations where the Dirac field obeys the 4-dimensional fermionic BF bound, m 2 ≥ 0, but violates the 2-dimensional BF bound. When this happens, i.e. when m 2 eff < 0, one should expect an instability. This condition can be rewritten: the 2-dimensional BF bound is violated if the charge of the fermion is larger than The equality applies strictly to the extremal case; as we move away from extremality, by continuity the instability should still be present but a higher fermion charge is needed to trigger it. Fig. 5 illustrates the regions where (4.4) predicts instability/stability. At this level, we see that the near-horizon analysis of the possible violation of the AdS 2 BF bound for a Dirac field parallels very much the partner analysis done for a scalar field in section 3, with the minimum value for the charge (3.7) for the scalar case just replaced by the fermionic minimum value (4.4). In the scalar field case, we found (through a numerical study of linear perturbations in the AdS-RN background) that the violation of the 2-dimensional BF bound translates into the existence of a linear scalar condensation instability. Moreover, the near-horizon scalar bound (3.7) turns out to be very sharp, as best illustrated in Fig. 1. This scalar condensation linear instability indicates that nonlinearly the AdS-RN black hole, when perturbed by a scalar field evolves towards a new configuration − a hairy black hole (with a scalar condensate floating above the horizon) − with the same UV asymptotics (since the 4-dimensional BF bound is satisfied) but with a different near-horizon geometry where the 2-dimensional BF bound is no longer violated [49,5,50,6].
These considerations motivate the study done in this manuscript for a Dirac field. In this case the AdS 2 BF bound can also be violated: at extremality this occurs for a fermion charge that saturates (4.4). From the lessons learned in the scalar field case one might well expect that the AdS-RN black hole, when perturbed by a Dirac field, is linearly unstable. To confirm whether this is the case, in the rest of this section we will solve numerically the Dirac equation in the AdS-RN background to hunt for a signature of the near-horizon linear instability. However, unlike the scalar field case, we will not find any evidence of a linear instability.

Dirac normal modes of global AdS
Before looking for potential instabilities (or frequency spectrum of damped oscillations) of Dirac modes in the global AdS-RN black hole it is convenient to first compute the normal mode spectrum of Dirac fields in global AdS. Indeed, some families of AdS-RN perturbations must reduce to these in the limit where the horizon shrinks to zero. Massive (section 4.2.1) and massless (section 4.2.2) Dirac fields require a distinct analysis.

Massive normal modes
For massive fermions in global AdS, it is not easy to solve directly the Dirac equations to get the radial functions R 1,2 . Perhaps for this reason, the massive normal modes of global AdS have not yet been computed. There is however an appropriate combination of R 1,2 that yields equations of motion that are explicit hypergeometric equations. The linear combination for R 1,2 that we use below is motivated by a similar analysis done to compute the massive normal modes of fermions for de Sitter in [51].
For m = 0 and in global AdS, we introduce the new radial variable y = −ir/L and make the following field redefinitions where f 1,2 (y) are functions to be determined. In these conditions, the coupled system of Dirac equations (2.22) yields Adding and subtracting these two ODEs yields This pair of coupled first order ODEs can be straightforwardly rewritten as a decoupled pair of second order ODEs for f 1 and f 2 . Moreover, if we introduce the new radial coordinate z = y 2 and the field redefinitions each of the ODEs becomes a hypergeometric ODE with the standard form with parameters a i , b i and c i given by The most general solutions of (4.8) are [52] is the Gaussian (ordinary) hypergeometric function and A 1,2 , B 1,2 are arbitrary amplitudes. We can now plug (4.11) into (4.7) and then into (4.5) to get the most general solution for R 1 (r) and R 2 (r). Finally, but not least important, we can insert this most general solution for R 1,2 (r) into (2.38) to get the most general solution for the Dirac fields ψ ± (r). These are the physical fields that have to be regular everywhere and this constrains some of the amplitudes A 1,2 and B 1,2 and the frequencies. Namely, at the origin, r = 0, one finds that both ψ ± have two divergent terms of the form B 2 /r +3/2 and (2A 2 − B 2 )/r +1/2 . Regularity at the origin thus requires that one sets A 2 = 0 and B 2 = 0 and the other two amplitudes A 1 and B 1 are left arbitrary. It follows that the regular normal eigenmodes are We have not yet imposed the asymptotic boundary condition. A Frobenius analysis of (4.12) near the conformal boundary together with the use of (2.38) finds that ψ ± behaves as (2.39) or (2.41) with For m > 0 (m = 1/2) the no-source standard boundary condition (2.42) requires α 1 = 0. Using Γ[−p] = ∞ for p = 0, 1, 2, · · · this quantizes the frequency as ωL = + 2 + mL + 2p or ωL = −( + 1 + mL + 2p) , (standard quantization) (4.14) For 0 < mL < 1 2 we can also impose the alternative quantization (2.43), i.e. β 1 = 0. This quantizes the frequency spectrum as (also with radial overtone p = 0, 1, 2, · · · )

Massless normal modes
In this section we find the normal modes in global AdS for a massless fermionic field. These have been previously discussed in [53,27] but these references have not identified the full spectra of frequencies.
For m = 0 in global AdS, introducing the change of coordinates and field redefinition the radial equation (2.24) can be rewritten as a hypergeometric ODE in the standard form Its most general solution is [52] Introducing this into (4.16) one gets R 1,2 (r) (note that R 2 = R * 1 as discussed in the next section). Plugging this into (2.38) one finds the most general solution for the Dirac fields ψ ± (r). At the origin, r = 0, these ψ ± have a divergent term proportional to C 2 r − −3/2 . Regularity at the origin thus requires that we set C 2 = 0. Now we need to impose the asymptotic boundary condition. One finds that asymptotically ψ ± decays as (2.40) with As explained previously, for m = 0 we can impose either the standard or alternative boundary conditions. The no-source standard boundary condition (2.42), α 1 (λ + ωL) − iβ 1 = 0, quantizes the frequency spectrum as ωL = + 2 + 2p and ωL = − ( + 1 + 2p) , p = 0, 1, 2, . . . , (Standard quantization) (4.20) On the other hand, for the no-source alternative quantization (2.43), α 1 (λ − ωL) + iβ 1 = 0, the normal mode frequencies of a massless Dirac field in global AdS are: ωL = +1+2p and ωL = − ( + 2 + 2p) , p = 0, 1, 2, . . . , (Alternative quantization) (4.21) The positive frequencies in (4.20) and (4.21) were computed in [27] using vanishing flux boundary conditions that, as explained in the end of section 2.3, are exactly the AdS/CFT standard and alternative boundary conditions. However, [27] missed the existence of half of the normal mode spectrum, namely the half part that has negative frequencies.

Setup of the numerical problem
In this section we solve numerically the Dirac equation and search for linear instabilities of the Dirac solution in the AdS-RN background. Before proceeding it is important to note that: 1) the Dirac radial equation (2.24) for R 2 (r) is just the complex conjugate of the radial equation for R 1 (r) so if R 1 (r) is a solution one automatically has R 2 (r) = R 1 (r) * , and 2) the Dirac angular equations for S 1,2 are related by the symmetry θ → π − θ so if S 1 (θ) is a solution then S 2 (θ) = S 1 (π − θ). Therefore, we just need to find the solutions R 1 (r) (S 1 (θ) are just the spin-weighted s = 1/2 spherical harmonics with quantum number ). Further note that if R 1 has charge q then R 2 = R * 1 has charge −q, and complex conjugation maps quasinormal modes to quasinormal modes. It follows that if ω = ω r + iω i is a linear mode frequency of R 1 then −ω * = −ω r + iω i is a linear mode frequency of R 2 = R * 1 . Thus, if we compute the frequency spectrum of R 1 , we have the spectrum of R 2 too. It also follows that there is no loss of generality in assuming that qQ > 0 in our analysis: results for qQ < 0 are obtained simply by reversing the sign of the real part of the frequencies. Finally note that when we compute ω we have to allow both positive and negative values of ω r , i.e. if ω r = Re(ω) is a frequency of R 1 there is no symmetry in the system that requires −ω r to be also a frequency of R 1 . The only exception is if µ=0 (or e = 0) and m = 0 i.e. a massless Dirac field in Schwarzschild-AdS. In this case if For concreteness, we will set the mass of the fermion to zero, i.e. we solve the Dirac equation (2.24) for R 1 with m = 0 (and in the gauge A| ∞ = 0 where the frequency of the fermionic wave is ω; see footnote 11) subject to the physically relevant boundary conditions. The asymptotic decay of R 1 (r) is given in (2.29). For reasons discussed previously, we impose the asymptotic boundary condition (2.42) (standard quantization, ψ (4.22) Rewriting this in ingoing Eddington-Finkelstein coordinates (v, r, θ, φ), with v = t + f −1 dr, which are smooth across the future event horizon H + , we find that regularity of R 1 (r) at H + requires that we impose the boundary condition B out = 0. 21 For the numerical solution it is convenient to redefine 23) and to work with the compact radial coordinate such that the horizon is now at z = 0 and the asymptotic infinity at z = 1. This has the advantage that analytical solutions that are smooth at the horizon simply have to obey the horizon Neumann boundary condition q (z = 0) = 0. The asymptotic boundary condition for q(z) follows straightforwardly from (2.42)-(2.43) and (4.23)-(4.24).
The numerical methods that we use are very well tested [56-61, 16, 44, 29, 30, 62] and reviewed in [63]. To discretize the field equations we use a pseudospectral collocation grid on Gauss-Chebyshev-Lobbato points. The eigenfrequencies and associated eigenvectors are found using Mathematica's built-in routine Eigensystem. For a given , this method has the advantage of finding several modes (i.e., from distinct families and with distinct radial overtones) simultaneously. However, to increase the accuracy of our results at a much lower computational cost we use a powerful numerical procedure which uses the Newton-Raphson root-finding algorithm discussed in detail in section III.C of the review [63]. All our results have the exponential convergence on the number of gridpoints, as expected for a code that uses pseudospectral collocation. In particular, all the results that we present are accurate at least up to the 10th decimal digit.
The Dirac equation in AdS-RN also has the scaling symmetry that determines that the physical dimensionless quantities are those listed in (3.9).

Results
As discussed in section 4.3, for fermion mass m = 0, we can have two independent boundary conditions that yield normalizable modes: the standard (ψ (0) + = 0) and alternative (ψ (0) − = 0) boundary conditions. Moreover, for each of these boundary conditions, the eigenvector R 1 can have negative or positive real part of the frequency. It follows that, for a given harmonic , m φ (and m = 0) we have a total of two frequency spectra to discuss for each one of the two possible boundary conditions. Before proceeding to the actual physical analysis of the frequency spectrum and instabilities of the system in subsection 4.4.2, we first test our numerical code by comparing the associated numerical results with some analytical perturbative expansions that are presented in Appendix B.

Comparing numerical results with analytical expansions
To test our numerical code we first compute the quasinormal mode frequencies for a massless Dirac field in global AdS and in Schwarzschild-AdS: see Fig. 6 for the standard (2.42) and alternative (2.43) quantizations. When r + = 0 our frequencies reduce to the normal mode frequencies of global AdS (4.20) and (4.21) for the standard and alternative quantizations, respectively. Recall that is the spin-weighted harmonic number and n is the radial overtone (related to the number of nodes of the wavefunction along the radial direction). So the smallest |ω| is obtained for = 1/2 and n = 0. For the alternative boundary condition these are ωL = 3/2 and ωL = −5/2, while for the standard quantization these are ωL = 5/2 and ωL = −3/2). All the numerical results we will present describe solutions with = 1/2 and n = 0. If there is an instability it should already be present in this sector of perturbations (see the argument in section 4.1).
On the other hand, for finite r + /L our numerical curves reproduce the values first computed in [27]. As explained in the end of section 2.3, this is because the AdS/CFT standard (2.42) and alternative (2.43) quantizations have vanishing energy flux and correspond precisely to the boundary conditions imposed in [27]. To complete the spectrum of Schwarzschild-AdS (and global AdS) note that frequencies that are the negative of the complex conjugate of the frequencies (−ω * ) plotted in Fig. 6 are also eigenvalues of the system (which was missed in [27]).
Next we test our numerical code for AdS-RN. As a first test, in Appendix B.1 we use a matching asymptotic expansion method to find that for ωr + 1 and qQ 1 (and m = 0, = 1/2, n = 0) the imaginary part of the frequency is approximately given by In Fig. 7 we plot Im(ωL) as a function of r + /L, for µ = 1 2 µ ext and a fermion charge qL = 1 (blue dots). We compare these numerical results with the analytical result (4.25). Both agree for small r + /L, i.e. in the regime were the matching expansion (4.25) is valid.
As a second test, in Appendix B.2 we use a perturbative expansion in r + /L about global AdS to find that the real part of the frequency behaves as Figure 6. Real (left) and imaginary (right) part of the QNM frequencies as a function of the horizon radius for the two families of boundary conditions in Schwarzschild-AdS background (m = 0, = 1/2). The black curve that reduces to the normal mode of AdS ωL = 5/2 has the standard quantization (2.42) while the red curve (that reduces to the normal mode of AdS ωL = 3/2) has the alternative quantization (2.43). Not plotted, if we take the black curve data then −Re(ω) + i Im(ω) is also a solution with alternative quantization (2.43) (which reduces to the AdS normal mode ωL = −5/2 when the horizon shrinks). Similarly, if we take the red curve data then −Re(ω)+i Im(ω) is also a solution with standard quantization (2.42).
The plots for the real part of the frequency are displayed in figure 8. The analytical perturbative results are (4.27) and (4.28). An important feature of this result is that for a fixed chemical potential, the slope of the real part of the frequency (the term proportional  to R + ) changes sign for a certain electric charge. This coincides with our numerical findings. In the top panels we show the real part of the frequency for the boundary condition (2.43) with qL = 0.7, 0.8 and µ = µ ext (1 − 10 −3 ). In the bottom panels of figure 8 we show the real part of the frequency for the boundary condition (2.42) as a function of the horizon radius for fixed chemical potential µ = µext 2 and fermion charges qL = 1.9, 2 from left to right. The dashed red lines are the analytical results (4.28) and (4.27) respectively. Note that for these parameters the change of sign occurs for charges qL ∼ 1.95 and qL ∼ 0.71 respectively.
As a general comment, it is perhaps worthy to comment that we find that matching asymptotic and perturbative analysis like the ones provided in Appendices B.1 B.2 are valid for smaller windows of r + /L in the Dirac field case when compared with similar analysis done for the scalar field.
We have confirmed that our numerical code is generating physical data. Thus in the next section we proceed to the numerical analysis of the system.

Numerical results
Our aim is not to present the full spectrum of frequencies of a Dirac field in AdS-RN black hole. Instead, we are motivated to search for unstable modes, i.e. on eventually finding modes that, in some range of parameters, have Im(ωL) > 0. There is a wide window of parameters to explore although the instability, if it exists, should appear near extremality for fermion charges q above a critical value. Thus one needs a good strategy to hunt efficiently for unstable modes. We proceed as follows. From the near-horizon bound (4.4) arguing for the existence of an instability, we see that this bound is lower if we set m = 0 and = 1/2. So, first, we either: 1) fixed m = 0, = 1/2 and µ close to µ ext , and varied {r + /L, qL}, or 2) fixed m = 0, = 1/2 and qL and varied {r + /L, µ}. In both cases, as described in the end of section 4.3, we solved our system as an eigenvalue problem for ωL. This finds "all" the solutions of the system (as long as the hierarchies do not grow large, e.g. |ωL| 1, which makes the numerical problem hard). This allows to eventually identify unstable modes with Im(ωL) > 0 or, in the worst case, to identify modes with Im(ωL) < 0 that are closest to the marginal case for instability (Im ω = 0). Once these interesting modes are identified we then used a Newton-Raphson root-finding algorithm to follow efficiently the modes to other values of the parameter space.
In spite of our efforts, we have found no sign of an instability. Recall that for mL = 0 both the standard (2.42) and the alternative (2.43) boundary conditions yield normalizable modes. In general, we do find that the stable modes with smallest |Im(ωL)| are those that reduce to the alternative normal modes of AdS (4.21) or to the standard AdS normal modes (4.20), when the horizon radius shrinks to zero. Among these, we further find that the modes with smallest |Im(ωL)| are, for both quantizations, the ones that reduce to the positive normal mode frequencies when r + /L → 0, i.e. ωL = 3/2 (alternative boundary condition) and ωL = 5/2 (standard quantization). Therefore, to avoid distraction from the main point, in the rest of this manuscript we only discuss these two families of modes.
Probably the plots that best illustrate the main conclusions of our Dirac study are those of Fig. 9 (for alternative quantization) and of Fig. 12 (for standard quantization). Recall that in the best case scenario the expectation is that, close to extremality, modes should become unstable above a fermion charge q that should be higher than the nearhorizon bound (4.4). Thus, in these figures we fix the black hole horizon to be r + /L = 0.5 and choose a chemical potential close to extremality, µ = 0.99µ ext . Starting from qL = 0, where Im(ωL) < 0, we then increase this charge to see if there is a critical value above which Im(ωL) becomes positive. (That is to say, we adopt a similar strategy as the one followed in the scalar field case to get Fig. 2).
For the alternative quantization, the left panel of Fig. 9 shows that, starting from q = 0, as qL grows, Im(ωL) < 0 increases and approaches Im(ωL) = 0 very closely. However, no matter how large qL is we never reach a situation where Im(ωL) ≥ 0. Interestingly, there is a critical value of q, namely qL = q max L ∼ 0.9390 (vertical brown dashed line) where Figure 9. Dirac field frequency with alternative quantization (2.43) as a function of the dimensionless scalar charge qL for a AdS-RN black hole with µ = 0.99µ ext and r + /L = 0.5 (also, mL = 0 and = 1/2). Left panel: Imaginary part of the dimensionless frequency, Im(ωL) which attains a maximum for q = q max ∼ 0.9390/L (vertical brown dashed line). The inset plot zooms-in around this maximum and shows that Im(ωL) < 0 for any qL. Right panel: Real part of the dimensionless frequency, Re(ωL), measured with respect to qµL. This quantity changes sign at q = q ∼ 0.9344/L with q < q max , i.e. for a smaller qL than the one where Im(ωL) attains its maximum value (vertical brown dashed line): this is better seen in the inset plot which zooms-in the relevant region.
We have done similar experiments as those of Fig. 9 for other black hole parameter values µ and r + /L. Keeping µ fixed, black holes with distinct r + /L have plots similar to Fig. 9 with the feature that larger values of r + /L reach the maximum of Im(ωL) (but remaining negative) at smaller critical values of q = q max . On the other hand, keeping r + /L fixed, black holes with distinct µ also have similar plots to Fig. 9 with the property that larger values of µ reach the maximum of Im(ωL) (but still negative) at smaller critical values of q = q max and this maximum of Im(ωL) is increasingly closer to zero as µ approaches the extremal value µ ext . To have a complementary perspective of the system's properties, in Fig. 10 and in Fig. 11 we illustrate other attempts we have made to find an instability. In Fig. 10, we keep the alternative quantization and fix the chemical potential at µ = 0.999µ ext , and plot Im(ωL) as a function of r + /L for five different values of qL, namely, qL = 0.7, 0.8, 0.9 (from bottom to top in the left panel) and qL = 1, 1.1 (right panel). (The two plots are needed for the presentation of the results because the relevant qL = 1 case in the right panel reaches a maximum that is approximately two orders of magnitude higher than the first three cases on the left panel). The main feature in these plots is the typical presence of a local minimum and local maximum (bump). As we increase the fermion charge from zero to a value slightly above 1, the relative minimum and relative maximum of Im(ωL) raise and shift to lower values of r + /L. But the local maximum always has Im(ωL) < 0, i.e. there is no instability. However, for charges qL above a value that is in between 1 and 1.1, the local minimum and maximum are no longer present and Im(ωL) decreases monotonically with r + /L (see e.g. qL = 1.1 displayed as the yellow curve in the right panel; higher values, qL ≥ 1.1, have a similar monotonic behaviour).
As yet another illustration of experiments we made, in Fig. 11 we fix the fermion charge to be qL = 1 (which was already analysed in Fig. 10 for µ = 0.999µ ext ) and we study the effect that changing the chemical potential has by considering a total of 5 curves with 5 different values of µ. Namely, in the left plot we consider the cases µ = 0.9µ ext and µ = 0.95µ ext . These cases have no bump (no local maximum) and illustrate that it only appears close to extremality. In the right panel we show three more cases where we fix µ = 0.99µ ext , µ = 0.995µ ext and µ = 0.999µ ext (from bottom to top). The bump is now present and the local maximum increases as one approaches extremality but never becomes positive. For the case µ = 0.999µ ext this local maximum is at Im(ωL) ∼ −5.93 × 10 −6 . So far we have focused our discussion of the results for the alternative quantization case because, typically, for the same values of black hole parameters this is the case where Im(ωL) approaches Im(ωL) = 0 the most. Nevertheless, we have also tried hard to find an instability in the standard boundary condition (2.42) case. Again without success. To illustrate briefly this conclusion, in Fig. 12 we give the partner plot of Fig. 9 but this time for the standard quantization. Although the features of Fig. 12 are clearly more elaborated than those of Fig. 9 (e.g. there are several local maxima and minima), the main conclusions are still the same: i) one always has Im(ωL) < 0; ii) there is a q = q max where the solution approaches Im(ωL) = 0 the most (vertical brown dashed line); iii) Re(ωL) − qLµ changes sign at q = q < q max . It follows that we find no sign of an instability and the standard boundary condition case, much like the alternative quantization case, also gives results that are very different from the scalar field case of Fig. 2.
Altogether, all our attempts − best illustrated in Figs. 9-11− to find an instability due to Dirac field perturbations with alternative boundary condition (2.43) failed miserably. The outcome is similar when we consider the standard boundary condition (2.42). This is best summarized in Fig. 12 where we show the partner plots of the alternative quantization Fig. 9 but this time for the standard quantization. Again, and essentially, we find that Figure 12. Dirac field frequency with standard quantization (2.42) as a function of the dimensionless scalar charge qL for a AdS-RN black hole with µ = 0.99µ ext and r + /L = 0.5 (also, mL = 0 and = 1/2). Left panel: Imaginary part of the dimensionless frequency, Im(ωL) which attains a maximum of Im(ωL) ∼ −0.0037491 for q = q max ∼ 3.3555/L (vertical brown dashed line). The main inset plot zooms-in around this maximum and shows that Im(ωL) < 0 for any qL. The secondary inset plot shows the detail of the curve around qL ∼ 1.45 to show that the apparent cusp in the main plot is smooth. Right panel: Real part of the dimensionless frequency, Re(ωL), measured with respect to qµL. This quantity changes sign at q = q ∼ 3.2873/L with q < q max , i.e. for a smaller qL than the one where Im(ωL) attains its maximum value (vertical brown dashed line): this is better seen in the inset plot which zooms-in the relevant region.
Im(ωL) reaches a maximum at a critical q max but never crosses the borderline Im(ωL) = 0 that would signal an instability.

Discussion and conclusions
A scalar field in an asymptotically AdS 4 Reissner-Nordström black hole can satisfy the asymptotically AdS 4 UV Breitenlöhner-Freedman (BF) bound but violate the infrared 2dimensional BF bound associated to the AdS 2 × S 2 near-horizon geometry of the extremal black hole of the system, as reviewed in section 3. When this is the case, the AdS-RN black hole is unstable to scalar condensation and the system evolves to a new configuration in the phase diagram of solutions that preserves both the UV BF bound and the near-horizon 2-dimensional BF bound. Such a solution is a hairy black hole with a charged scalar field floating above the horizon [14,15,64,49,44,65,50]. Coulomb repulsion balances the gravitational force and the system is static. There is no doubt that the violation of the AdS 2 BF bound is the physical mechanism responsible for the near-horizon scalar condensation instability since the associated minimum bound (3.7) on the scalar field charge that triggers the instability is sharp (at extremality) as best demonstrated by Fig. 1.
Given these considerations, the study done in this manuscript for Dirac field perturba-tions in the global AdS 4 RN black hole was motivated by the following observation. Dirac fields in AdS-RN can also preserve the UV fermionic BF bound but violate the near-horizon infrared fermionic BF bound, as seen in section 4.1. From the scalar field case lessons, this suggests that the system might be unstable to fermion condensation. However, in spite of our efforts to scan the relevant parameter space near extremality, we found no sign of a linear mode instability. The sharp distinction between the scalar and Dirac field cases is best illustrated comparing the scalar Fig. 2 with the Dirac Fig. 9 (for alternative quantization) or Fig. 12 (for standard quantization). Of course our numerical study does not prove linear stability but we did such a detailed scan that we are very confident that no linear instability is present. Our stability results are also consistent with the stability study of fermions in planar AdS, where no instability was found [15,17] (see also [18][19][20][21][22][23]). Indeed, the planar AdS case is the r + /L → ∞ limit of the global AdS system. So, the planar AdS studies [15,[17][18][19][20][21][22][23] and our present study in global AdS establish that the violation of the 2-dimensional BF bound of a Dirac field in AdS-RN does not lead to a linear mode instability. But the system might still be unstable if non-linear effects play a role in the discussion. That is to say, if we perturb a AdS-RN black hole with a Dirac field in a region of parameter space where the infrared BF bound is violated, it could still be the case that the system evolves non-linearly to a new configuration that has a charged Dirac field floating above the horizon and that preserves both the UV and IR BF bounds. How difficult would it be to prove whether this scenario is correct?
One must proceed with caution. To begin with one needs to first formulate more precisely the setup of the problem. It is certainly much harder to find, if they exist, the proposed Dirac hairy black holes than it was to construct the scalar hairy black holes [14,15,64,49,44,65,50]. There is a fundamental difference between fermionic and bosonic fields. The fermion has no classical limit: Planck's constant is present in the stress tensor and associated equations of motion for a fermion. As discussed in detail in section 14.3 of Wald's textbook [66], the absence of a classical limit means that in Einstein's equation we have to promote the differential Einstein and energy-momentum operators G and T to quantum operators and the quantum version of Einstein's that gives the back-reaction of a fermion on the gravitational field is G µν = 8π( T Max µν + T µν ), where T Max µν and T µν stands for the Maxwell and Dirac stress tensor contributions and · · · stands for the expectation value of the corresponding operator.
Thus, to find the backreaction that fermions induce on the gravitoelectromagnetic background one needs to first compute the expectation value of the fermion energy momentum tensor T µν . This is a highly non-trivial task. Even worse, once we consider the quantum backreaction of fermions one also needs to consider the quantum backreaction of gravitons and photons, i.e. one also needs to compute G µν and T Max µν [66]. In a best case scenario, where we have a large number N of Dirac fields one might be able to assume that, roughly speaking, the effects of N Dirac fields are N times as relevant as that of the gravitons and photons [66]. For a 'fermion hairy black hole' , the fermionic condensate should be made of a large number of fermions. In these conditions, for large N , one might be able to neglect the quantum backreaction of gravitons and photons and work in the semi-classical limit whereby the backreaction of the Dirac field on the gravitoelectromag-netic background is simply governed by G µν = 8πT Max µν + 8πN T µν . This semi-classical system should be viewed as the leading term of a 1/N expansion of the full theory [66]. But this semi-classical computation still requires that one computes T µν . And this is still a remarkable task. An overview on the physical and technical tools required to accomplish this task can be found in [67,68,23] (and references there-in) where asymptotically planar AdS quantum electron stars are discussed as semi-classical solutions of Einstein-Maxwell theory.
Finally note that in the present manuscript we focused our attention on modes that could eventually become unstable. We have not studied in detail the full spectrum of quasinormal mode frequencies of the Dirac field in AdS-RN. Moreover, we focused on the case of a massless fermion because, as explained previously, this was enough for our purposes. However, the equations of motion and relevant boundary conditions for any fermion mass and any sector of perturbations are given in section 2. We have also computed the normal modes of massive fermions in AdS (previously only the massless spectrum was computed). It might be useful to have a more complete frequency spectra study for future studies/applications. For x 1 one is close to the outer horizon and for σ 1 the inner and outer horizon are very close, i.e. one is close to extremality. Next, we take the limit σ → 0 whilst keeping z = x σ fixed. From [24][25][26] the 'near-extremal' modes are expected to saturate the superradiant bound ω = qµ at extremality so onwards we measure the frequency difference δω with respect to this bound via the redefinition Using the condition f (r − ) = 0 for the location of the inner horizon one can find r − = r − (r + , Q, L) which is then inserted into (A.1) to express Q as a function of (r + , σ, L).
In these near-extremality conditions, we are ready to find the near-horizon solution of the Klein-Gordon equation. Concretely, introducing the above redefinitions into the Klein-Gordon equation (3.2), to leading order in σ, we obtain: where we have introduced the dimensionless quantities R + = r + /L, e = qL,m = mL, δω = L δω.
With the field redefinition 3) is rewriten as a standard hypergeometric ODE The regular (i.e. ingoing) solution at the future event horizon is given by where 1 F 2 (a, b, c; z) is the standard hypergeometric function and a − , a + are defined by: In the context of a matched asymptotic expansion, the near-region (near-horizon) solution (A.7) must now be matched with the far-region solution of (3.2) (in near-extremality conditions). As explained above we expect the 'near-extremal' modes we are looking into to have wavefunctions that die-off very quickly away from the black hole horizon (at least near-extremality). Therefore, as a first rude approximation we take the far-region to be described by a vanishing wavefunction. That is to say, in the overlapping region, we match the near-region solution (A.7) with φ = 0. In the end of the day, this approximation turns out to be quite good because the analytical approximation for the 'near-extremal' frequency that we obtain −see (A.11) − matches remarkably well the numerical solution of (3.2). This is best seen comparing the black dashed analytical curve of our expansion in Fig. 2 with the numerical blue dot results. For this reason, we do not try to improve further our matching asymptotic approximation.
Proceeding in these conditions, the leading order behaviour of the large R = r/L series expansion (z → −∞) of φ, namely φ ≈ (−z) ± √ 1+4η , needs to be matched with the farregion solution φ = 0. Before we can do it, we still need to distinguish the cases 1 + 4η ≥ 0 and 1 + 4η < 0. For our proposes (comparing with the numerical results of section 3), we want to consider the small scalar field charge case for which one finds that 1 + 4η ≥ 0 holds as long as e ≤ e c where e 2 c = 1 + 6R 2 + 8R 2 + (1 + 3R 2 + ) 6 + 4m 2 R 2 + + (1 + 2 ) 2 (A.9) (the reader also interested in the case 1 + 4η < 0 can follow the steps detailed in [26]). In these conditions, from the matching condition one finds that where p = 0, 1, 2, · · · is the radial overtone of the mode. Replacing this into (A.2) one finally finds that 'near-extremal' modes have a frequency given by: ωL eµ + σ δω + O(σ 2 ). (A.11) This is (3.10) in the main text when we set the radial overtone p = 0.

B Perturbative results
Although we solve the Dirac equation numerically in the main text (section 4.4.2), it is very good practice to testify the numerical results against analytical predictions that can be obtained within perturbation theory in some region of the parameter space. Therefore, in this appendix we find some useful analytical perturbative approximations for the Dirac frequencies. More concretely, this Appendix is divided in two parts. In Appendix B.1 we use a matching asymptotic expansion approach to find the frequency approximations (4.25)-(4.26) that, in Fig. 7 of section 4.4.1 of the main text, are compared against our numerical results. Then, in Appendix B.2 we use a systematic perturbative expansion in the dimensionless horizon radius r + /L 1 (with no further approximations) to find the analytical frequency approximations (4.27)-(4.28) which, in Fig. 8 of section 4.4.1 of the main text, are also compared with our numerical results. In both cases, there is agreement between the analytical approximation predictions and the numerical results in the regime of parameter space where the former are valid.
In this appendix, as in the main text, we solve the Dirac equation (2.24) in a AdS-RN background (2.1) (gauge choice C = 0) with regular (ingoing) boundary conditions at the future event horizon and the standard (2.42) or alternative (2.43) boundary conditions at the conformal boundary. We will work exclusively with vanishing fermion mass, m = 0.

B.1 Matched asymptotic expansion
In this section we derive an analytical expression for the imaginary part of the Dirac frequency using the method of matched asymptotic expansion introduced in [69,70] (see also e.g. [71][72][73]). We assume r + L and split our spacetime into two regions; an asymptotic globally AdS far region where the effects of the black hole can be neglected and a near region about the black hole outer horizon where the effects of the cosmological constant can be neglected. In each region the associated perturbation equation can be solved analytically, then matching the near and far region solutions in their overlapping region will fix the integration constants as well as the imaginary part of the perturbation frequency ω. More concretely, the near region is defined by r −r + 1/ω and the far region is defined by r − r + r + . It follows that the overlapping region exists for ωr + 1. A further assumption we must make is that the Coulomb interaction is weak, Qq 1, where q (Q) is the fermion (black hole) charge.

B.1.1 Far region solution
Since in the far region the effect of the black hole (BH) is assumed to be negligible, we effectively have a fermion field in the global AdS background. Thus, the general solution for the massless fermionic field R 1 is given by (4.18) that we reproduce here: F (r) = C 1 2 F 1 1 2 + , ωL + + 1, 2(1 + ), 2r r + iL +C 2 2r r + iL where C 1,2 are two arbitrary amplitudes to be determined below. Asymptotically this solution decays as (2.29), namely R 1 r→∞ ∼ α 1 + β 1 L r + O(r −2 ) with α 1 = i ωL 2 + 1 2 C 1 2 F 1 + 1 2 , + ωL + 1; 2 + 2; 2 + C 2 2 F 1 − − 1 2 , − + ωL; −2 ; 2 , This solution has to satisfy the asymptotic boundary condition (for a massless fermion). For the standard quantization this is (2.42) while for the alternative quantization the boundary condition is (2.43) which fix β 1 as a function of α 1 or, equivalently, C 1 as a function of C 2 . This yields where the upper sign refers to the standard quantisation (2.42) and the lower sign to the alternate quantisation (2.43). Note that we do not impose any boundary condition at a inner boundary since this far region solution does not extend till there.

B.1.2 Near region solution
In the near region we can approximate ∆(r) = r 2 f (r) by: where r − ≈ Q 2 2r + . This follows from the assumption that r + L and therefore in the near region we have r ∼ O(r + ) L, so we can neglect the r 2 /L 2 term in f (r). Further applying the near region assumptions to the Dirac equation (2.24) we can neglect terms of order ωr + or higher powers. Other terms appear which are dominated by a 1/∆ term in the small black hole (BH) approximation near the horizon; therefore we can evaluate the numerators of these terms at r ≈ r + .
With these approximations, and the coordinate transformation (the horizon r = r + is at z = 0) the Dirac equation is approximately given by the near region equation, Making the field redefinition Using the property of the hypergeometric function 2 F 1 (a, b, c, 0) = 1 we find that the z → 0 behaviour of the near region solution is R near 1 (z) ≈ z 1 4 (αz −iσ + βz iσ ). Requiring regularity (only ingoing modes) at the horizon implies that we must set β = 0.

B.1.3 Matching
To find the large r (z → 1) behaviour of the near region solution (B.11) with β = 0, we use the z → 1 − z transformation law of the hypergeometric function [52]: We again use that 2 F 1 (a, b, c, 0) = 1 as well as 1 − z ≈ r + −r − r (when r → ∞), to obtain: In the overlapping region, one must have R near 1 large r = R far 1 small r . That is to say, we must match independently the r 1 2 + and r − 1 2 − terms of (B.13) with those of (B.14). This matching yields: where we have used (B.9) and (B.7) forσ to restore the explicit dependence on the frequency ω. The ratio C 2 C 1 follows straightforwardly from (B.3). We want to solve the transcendental equation (B.15) to get an analytical expression for the frequency. There is no closed form solution, unless we do some educated approximations that we now discuss. Since we are working with very small and weakly charged RN-AdS black hole, one expects that the mode frequencies in such a background are close to the massless Dirac normal modes frequencies of AdS 4 , already computed in (4.20) (standard quantization) or (4.21) (alternative quantization). Denote this normal mode frequency by ω AdS 4 . However, since the background now has a horizon, the system becomes dissipative and with respect to the normal modes of AdS 4 , the frequency of the system should acquire a small imaginary contribution. Denote it by i δ. In (B.15), it is thus a good approximation to replace the frequency ω by ω = ω AdS 4 + i δ with |δ| ω AdS 4 . Our target now is to solve (B.15) at leading order for δ ω AdS 4 ∼ O(1). A posteriori, we compare the prediction of our analytical computation with the numerical result to confirm that this approximation is valid. Equation (B.15) has the additional challenges that: 1) the frequency appears in the argument of the Gamma functions and 2) Γ [−1 − 2 ] in the numerator diverges for the allowed values (2.26) of the harmonic number (recall that Γ[−p] = ∞ for non-negative p). To deal with these obstacles, we use the Gamma function property Γ[z + 1] = zΓ[z] and the assumptions of our problem, ωr + 1 and δ ω AdS 4 ∼ O(1). This allows to expand the Gamma functions whose argument depends on ω (and thus on δ) to extract δ out of the argument of the Gamma functions. In particular, this permits to find that the divergence of Γ [−1 − 2 ] in the numerator is cancelled by the Gamma function in the denominator that depends on δ.
In these conditions, one finds that the leading order solution for δ is whereK is a positive real number for each , n. For the alternative boundary condition into a near and far regions. We expand the frequency Ω = ωL and the field R 1 in each region in a power series in R + = r + /L: We now series expand the Dirac equation (2.24) in small R + . The leading, zeroth order equation is simply the Dirac equation in global AdS 4 for a massless fermion (B.1) (that we already studied in sections 4.2.2 and B.1.1). Not surprisingly, the small R = r/L expansion of this leading order solution breaks down at order R + /R. This motivates splitting our spacetime into a far region, R R + , and a near region, R + ≤ R 1. In the far region we work with the radial coordinate R but in the near region we work instead with the radial coordinate y = R/R + (since the far region small R expansion breaks down at order R + /R).
In the far region, at each order in R + , we impose the standard boundary condition (2.42) or the alternative boundary condition (2.43). In the near region we impose boundary conditions that only allow for ingoing waves at the horizon. We then perform a matching procedure at each order in R + , in the region where the far and near region overlap, to determine the frequency coefficients Ω (k) , as well as amplitudes that were not fixed by the two boundary conditions. At leading (zeroth) order, we fix Ω (0) to be the normal mode frequency for a massless fermion already obtained in (4.21) or (4.20) for the boundary conditions (2.43) or (2.42), respectively. We will do this for the mode with harmonic number = 1/2, and radial overtone to be p = 0. Our aim is then to find the first frequency correction Ω (1) due to the presence of the black hole. For concreteness, in most of our discussion below we only explicitly present details of the case where we impose the alternative boundary condition (2.43). We then present the final result also for the standard boundary condition (2.42).
In the far region the leading order R 0 + solution is (B.1) and imposing the asymptotic boundary conditions amounts to repeat mutatis mutandis the analysis done in (B.1)-(B.3). With our choice of = 1/2 and p = 0 this fixes the frequency at order zero to be Ω (0) = 3 2 ; see (4.21). To fix the normalization, we set the amplitude of the Dirac field at infinity to be 1 at all orders in R + : R f ar 1 | R→∞ = 1 + O(1/R). Introducing the near region radial coordinate y = R/R + , still at leading order R 0 + , the near region Dirac equation (for = 1/2) reads The solution which is regular at the horizon (y = 1) is We must now match the far and near regions solutions at order R 0 + in their overlapping region R + R 1. This procedure, typically fixes all other constants of the problem that were not fixed by the boundary conditions. The large R expansion of ψ near (0) is ψ near (0) | large R = α (0) + β (0) R + · · · whereas the small R expansion of ψ f ar (0) is ψ f ar (0) | small R = (−1) 3/4 R + · · · . Therefore, matching ψ near (0) | large R = ψ f ar (0) | small R requires that we set α (0) = 0 and β (0) = (−1) 3/4 . Collecting the results at order 0 for the alternative quantization (2.43) we have: ; For the standard quantisation one has a similar result with Ω (0) = 5/2. We can now proceed to the first order R 1 + contribution that enables us to find Ω (1) . We repeat the procedure outlined above for the far and near regions, imposing the alternative boundary condition in the far region and ingoing boundary condition at the horizon in the near region. The near region equation for ψ near (1) is the same as at order R 0 + (i.e. there is no source). After imposing the horizon boundary condition, we find that the large R expansion of the near region solution gives R near large R ≈ (−1) 3/4 R + β (1) R + + · · · . (B.26) The far region equation at order R 1 + can also be solved analytically for ψ f ar (1) . As usual, we find a solution with two integration constants, say C 1 and C 2 . These constants will in general depend on µ, e and Ω (1) . Firstly by imposing the alternative quantisation we find an expression for C 2 in terms of µ, e and Ω (1) . Then we again impose that our field has amplitude 1 at infinity, yielding an expression for C 1 in terms of µ, e and Ω (1) . After doing this we find the small R behaviour of the far region solution to be R f ar small R ≈ (−1) 3/4 R + R + −14 + 2eµ − 7µ 2 − 6πΩ (1) − R + R 2 − 4eµ + µ 2 + 2πΩ (1) + · · · , (B.27) In the overlapping region R + R 1, (B.27) must match (B.26). A straightforward matching of the terms R + R 0 fixes the constant β (1) in (B.26). On the other hand, there is no R + /R term in the large R series expansion of the near region solution. It follows that Ω (1) must be such that it eliminates the corresponding R + /R term in the small R expansion of the far region solution (B.27). This fixes Ω (1) to be