Critical lumpy black holes in AdS${}_p\times S^q$

In this paper we study lumpy black holes with AdS${}_p \times S^q$ asymptotics, where the isometry group coming from the sphere factor is broken down to SO($q$). Depending on the values of $p$ and $q$, these are solutions to a certain Supergravity theory with a particular gauge field. We have considered the values $(p,q) = (5,5)$ and $(p,q) = (4,7)$, corresponding to type IIB supergravity in ten dimensions and eleven-dimensional supergravity respectively. These theories presumably contain an infinite spectrum of families of lumpy black holes, labeled by a harmonic number $\ell$, whose endpoints in solution space merge with another type of black holes with different horizon topology. We have numerically constructed the first four families of lumpy solutions, corresponding to $\ell = 1, 2^+, 2^-$ and $3$. We show that the geometry of the horizon near the merger is well-described by a cone over a triple product of spheres, thus extending Kol's local model to the present asymptotics. Interestingly, the presence of non-trivial fluxes in the internal sphere implies that the cone is no longer Ricci flat. This conical manifold accounts for the geometry and the behavior of the physical quantities of the solutions sufficiently close to the critical point. Additionally, we show that the vacuum expectation values of the dual scalar operators approach their critical values with a power law whose exponents are dictated by the local cone geometry in the bulk.


Introduction
The gauge/gravity duality has become an important tool to unravel the physics of certain quantum field theories (QFTs) at strong coupling, where the traditional techniques of perturbation theory fail. With the current understanding of the duality, its power mostly lies on the 'gravity side' of the correspondence, which involves (semi-)classical theories of gravity that, in general, are more tractable than their respective field theory duals. The first example of the gauge/gravity duality was introduced by Maldacena in [1]. This was originally coined as the AdS/CFT correspondence, since it relates a conformal field theory (CFT) such as N = 4 Super Yang-Mills in 3 + 1 spacetime dimensions and gauge group SU(N) to type IIB superstring theory in AdS 5 × S 5 . Although the equivalence was conjectured to be true for generic values of N and g YM (the gauge theory coupling constant), the duality is better understood in the large-N and large 't Hooft coupling λ = g 2 YM N limit, where the dual superstring theory reduces to classical supergravity in AdS 5 × S 5 [2,3]. By now, the duality has been extended to other dimensions and even to non-conformal settings. The literature on the subject is vast and we will not review it here (see the textbook [4] for a review and extensive references).
In this paper we are interested in equilibrium static black holes in asymptotically global AdS p × S q spacetimes. By the gauge/gravity duality, these black holes are dual to thermal phases of certain gauge theories. The field theory is defined at the boundary, B p−1 , of the pdimensional Anti-de Sitter (AdS) space, which enjoys a timelike conformal structure at infinity. 1 For instance, if B p−1 is taken to be Minkowski space, then the bulk spacetime is asymptotically AdS p in Poincaré coordinates. In the Poincaré patch there exist two homogeneous solutions of the Einstein's equations, namely pure AdS p and the planar Schwarzschild-AdS p black hole. The latter always dominates the canonical ensemble and hence there are no phase transitions. This is consistent with the fact that a CFT cannot have phase transitions on a scale invariant background such as Minkowski space. On the other hand, considering B p−1 to be a space with non-zero curvature [5] can lead to new and interesting physics. In the simplest case, one can take B p−1 to be the Einstein static universe (ESU), Ê × S p−2 , which corresponds to considering asymptotically global AdS p spacetimes in the bulk. In this case, there exists the well-known Hawking-Page phase transition [6] between global AdS p and (large) Schwarzschild black holes in that background. This phase transition was later re-interpreted by Witten [7] (see also Ref. [8]) as a confinement/deconfinement phase transition in the dual CFT.
So far the discussion has ignored the structure of the internal manifold, namely the S q in our case. There are situations where the internal geometry leads to new and interesting dynamics. This is famously the case for the Gregory-Laflamme (GL) instability of black branes [9], whose endpoint has led to counterexamples of the weak cosmic censorship conjecture in higher dimensions [10]. In fact, soon after the discovery of the GL instability, it was suggested that black holes in global AdS would undergo a similar dynamical instability that should lead to localization in the internal space [11,12]. The instability of black holes in global AdS was addressed in [13], who found the threshold zero (i.e., time-independent) mode(s) for the Schwarzschild-AdS 5 × S 5 black hole. In that paper it was shown that small (compared to the AdS radius) Schwarzschild-AdS 5 × S 5 black holes have zero modes (preserving a SO(5) isometry) when the negative modes of the Schwarzschild-AdS black hole [14] coincided with the momentum along the internal space.
For the latter, being a S 5 , modes are expanded in terms of spherical harmonics, which are labeled by the harmonic number ℓ. Then there exists an infinite number of zero modes, one for each (discrete) value of ℓ and the corresponding specific horizon radius.
The existence of zero modes suggests two things: First, the presence of a dynamical instability and, second, the existence of new branch of static black holes. Indeed, in the case of the GL instability of black strings, [15,16] first showed that there exists a new branch of non-uniform black strings that emanates from the uniform string branch precisely at the onset of the GL instability. Similarly, [17] constructed asymptotically global AdS 5 × S 5 black holes with nonuniform horizons along the S 5 that emanate from the Schwarzschild-AdS 5 × S 5 black holes at the onset of GL instability for the ℓ = 1 and ℓ = 2 modes respectively. These so-called lumpy black holes have the same horizon topology as their Schwarzschild-AdS 5 × S 5 counterparts, but no longer preserve the full isometry group of the internal space.
In the examples of [17], only a round S 4 inside the internal S 5 is preserved. The role of the lumpy black holes in the space of black hole solutions in AdS 5 × S 5 is analogous to that of the non-uniform black strings in the context of black holes in standard Kaluza-Klein (KK) theory. Similar to the fact that there exist localized black holes in KK theory [18], there also exist black holes in AdS 5 × S 5 that are localized on the S 5 . The latter were numerically constructed in [19]. The localized black holes of [19] should presumably merge with the ℓ = 1 lumpy black holes, but the branches of solutions constructed in [17] and [19] are still quite far from the merger point. One of the goals of this paper is to get close to this point from the lumpy black holes side.
Mergers in the space of static black hole solutions in general relativity have been thoroughly studied in the context of KK theory. Kol [20] proposed that the horizon topology change between the two branches is locally described by conifold-type-of transition controlled by a (singular) self-similar double-cone geometry. This double-cone model has been confirmed for the localized black hole/non-uniform black string transition [21] (see also the recent studies [22][23][24][25][26]) and in asymptotically flat space [27,28]. One non-trivial prediction of the double-cone model is that the behavior of any physical quantity Q near the critical point is given by where Q c is the value of Q at the critical point, ρ 0 is any length that measures the deviation from the cone, and s + and s − are the critical exponents that govern the perturbations away from the double-cone geometry. These critical exponents are dimension dependent, and they are complex for D < 10 and purely real for D > 10; D = 10 is a degenerate case that leads to a power law with a logarithmic correction. This behavior of the physical quantities near the critical point has been beautifully confirmed for D < 10 [22] and D = 10 [23,24]. Whilst physical quantities may be defined on the whole horizon, eq. (1.1) indicates that sufficiently close to the critical point in a topology-changing transition, they are actually governed by the regions of the horizon that pinch-off. So far no topology change transitions analogous to the localized black hole/non-uniform black string have been studied in AdS. In this paper we fill this gap by numerically constructing critical lumpy black holes. We review, correct and extend the work of [17] by constructing lumpy black holes in AdS 5 × S 5 and AdS 4 × S 7 with non-uniformity parameters (see Section 3 for the definition) λ ∼ 25 and λ ∼ 104 (for ℓ = 1) respectively. 2 Ref. [17] only considered the AdS 5 × S 5 case, and from their plots we have estimated that their most critical solutions have λ ∼ 2. For the ℓ = 1 lumpy family, the phase diagrams that we obtain exhibit some differences compared to those reported in [17], with our branches being shorter (see Section 4.1). This may appear to be a little surprising since we have managed to construct solutions with significantly larger values of λ, and hence closer to the critical point; these discrepancies can be attributed to the different quality of our numerical solutions and theirs. To get close to the critical point, we employ the sophisticated methods that we successfully used in [24] (see also [22,23] and Section 3.4). As we will show in Section 4.4, in the context of asymptotically AdS p × S q black holes, the critical geometry is controlled by a triple cone which is not Ricci flat anymore. We have checked that, to leading order, the behavior near the critical point of physical quantities of the different families of lumpy black holes that we have constructed (more specifically ℓ = 1, 2, 3) is given by (1.1). Furthermore, the critical exponents obtained from our numerical solutions (see Tables 3 and 4) match the prediction of the triple-cone model. Therefore, we are confident of the correctness and accuracy of our numerical solutions.
Another novelty of our work is that we use the holographic dictionary to extract how the physics of the topology change is imprinted in the dual field theory observables. Using the tools of KK holography [29], we extract the vacuum expectation values (vev's) of the boundary stress-energy tensor and certain scalar operators (see Section 4). 3 The latter parametrize the deformations of the internal space and, unsurprisingly, near the critical point they exhibit the behavior predicted by the triple-cone model (1.1). This result can be interpreted as a prediction from holography of how the dual field theory may detect a topology change in the bulk. We should point out that our results from KK holography differ from those of [17]. In particular, the expressions for the holographic stress-energy tensor and the scalar vev's are different, and we also find that O T 4 = 0. To correct possible errors and typos in the literature, we felt the need to provide an extensive review of KK holography and the details of our calculations in Appendices A and B.
So far we have focused on asymptotically AdS 5 × S 5 lumpy black hole solutions since these are the relevant ones for type IIB supergravity in ten spacetime dimensions. In this paper we also consider lumpy black holes in asymptotically AdS 4 × S 7 spacetimes, hence they are solutions to 11-dimensional supergravity/M-theory. 4 In this case, the only bosonic fields are the metric and the 4-form field strength. The equations of motion for the bosonic fields in this theory are simpler than in the IIB case (because there is no self-duality condition for the field strength), which allows us to get much closer to the merger point and test the triple-cone model with exquisite detail. Lumpy black holes in this theory should be dual to certain thermal phases of ABJM theory [30].
in the field. This should not be confused with the 't Hooft coupling.

The gravity problem
Our problem consists of finding static black holes that are solutions to certain supergravities and are asymptotically AdS p × S q (so D = p + q is the total number of spacetime dimensions) such that the S p−2 inside AdS and an S q−1 inside the S q are round. Asymptotically, the solutions must tend to the metric: ds 2 AdS×S = ds 2 (p) (AdS) + R 2 dΩ 2 (q) , (1.2) where R denotes the radius of the S q . R is fixed by supersymmetry in terms of the radius L of the AdS factor. For (p, q) = (5,5) and (4,7), the background metric (1.2) is a solution to the bosonic sector of type IIB SUGRA with only a RR self-dual 5-form field strength or 11D SUGRA respectively. The actions for the theories considered in this paper can be derived requiring local supersymmetry and we truncate them to the desired field content [31]. We have considered the following normalizations for the field strengths: where F (n) = dA (n−1) . In type IIB case, one has to additionally impose the self-duality condition F (5) = ⋆F (5) after deriving the equations of motion. A few comments are in order. First, notice that in D = 10, this is the minimal field content consistent with the symmetries. Since the S 3 within the AdS 5 factor is assumed to be round, a 3-form field strength would also be compatible with the symmetry, but it is not necessary. For simplicity, in this paper we have turned it off. We leave for future work the problem of exploring for phase diagram of asymptotically AdS 5 ×S 5 black hole solutions with an additional 3-form field strength. Second, in D = 11 the action contains a Chern-Simons term, ∼ F (4) ∧ F (4) ∧ A (3) , whose coefficient is fixed by supersymmetry. However, for static configurations and with the symmetry assumptions of this paper, the Chern-Simons term does not contribute and we omit it.
The equations of motion that one obtains from (1.3) can be written as follows: The indices run over the total number of spacetime dimensions. None of the equations for the metric has a 'bare' cosmological constant. Instead, this emerges from the flux of the gauge field which allows the background (1.2) to be a solution of (1.4) with: (1.5)

Numerical methods
We use the Einstein-DeTurck method to numerically integrate the equations of motion. This was originally introduced in [32] to solve static problems in vacuum Einstein's gravity, but it can be straightforwardly generalized to include matter. In fact, matter does not modify the principle part of Einstein's equations and therefore the non-vacuum Einstein-DeTurck equation: is manifestly elliptic on static metrics. ξ M is the so-called DeTurck vector and contains the usual Levi-Civita connection Γ compatible with the spacetime metric g, and a Levi-Civita connection Γ compatible with some reference metricḡ that we are free to prescribe. Additionally, there are also the matter field equations which will not be elliptic due to the underlying local gauge symmetry where Λ (n−2) is an arbitrary (n − 2)-form. The matter equations can also be modified by adding a DeTurck-like term, but in our case it is easier to fix the gauge algebraically with a suitable ansatz for the gauge field. For static spacetimes that are either asymptotically flat, KK or AdS, and whose boundary conditions are compatible with ξ M vanishing at the boundaries of the manifold (if any exist), [33] showed that all solutions to R H M N = 0 are necessarily Einstein. Even though the proof of [33] is for the pure gravity case, it should be straightforward to generalize to the present case. In fact, as we shall see below, for our particular problem, the boundary conditions that we impose fall within the class considered in [33], and hence solving (1.6) should be equivalent to solving the Einstein equations. In practice, we use the square of the DeTurck vector, ξ 2 ≡ ξ M ξ M , to quantify the validity of our numerical solutions and to perform convergence tests. More details about the DeTurck method can be found in [32] and in the reviews [34,35].
From the technical point of view, the equations are always discretized using pseudo-spectral methods on a Chebyshev grid, where it has been proven that exponential convergence is achieved for sufficiently smooth solutions. Then we solve them by a standard iterative Newton-Raphson method; the resulting linear system of equations is solved using LU decomposition (implemented by subroutine LinearSolve in Mathematica.)

Organization of the paper
The rest of the paper is organized as follows. In Section 2 we review the simplest black hole solution with the desired asymptotics, namely the Schwarzschild-AdS p × S q black hole. We do so for generic values of p and q, and we present the expressions for the physical quantities. In Section 3 we present the metric and gauge field ansatz, again for generic values of p and q, that we use for the numerical construction of the lumpy black holes. This section also includes a discussion about the boundary conditions that yield a well-posed boundary value problem. The one-parameter family that we use as ansatz captures, by construction, the threshold GL zero mode for any ℓ. We use this mode as the initial seed to find the actual non-linear solutions. The results are presented in Section 4. This section constitutes the core of the paper, where the phase diagrams are shown and discussed. Special emphasis is put on the isometric embeddings, the topology changes and the critical behavior of the solutions close enough to the merger point. In particular, we compute the critical exponents and find the critical values of various thermodynamic quantities. We close in Section 5 with a discussion and outlook of our results. Finally, we include a few appendices with technical details. KK holography in AdS 5 × S 5 is reviewed thoroughly in Appendix A, and detailed application of this procedure to our solutions can be found in Appendix B. This requires a basis of spherical harmonics that preserve an SO(5) subgroup of the full SO(6) rotational group of the S 5 ; this is presented for a generic value of q in Appendix C. In Appendix D we show the embedding plots in D = 10, which are not included in the main body of the paper. In Appendix E we include further details of the numerical construction, especially of the grids employed, and we perform some convergence tests.
2 Schwarzschild-AdS p ×S q black hole The simplest black hole solution with the desired asymptotics is the Schwarzschild-AdS p × S q black hole, which is known explicitly. In this section we review those aspects of this solution that are relevant for our analysis for generic values of p and q. Most of the material in this section is well-known and can be skipped by the expert reader.
The horizon of Schwarszchild-AdS p × S q black hole is uniformly smeared along the internal S q and it thus respects the full SO(q + 1) symmetry; the horizon topology is S p−2 × S q . The metric reads: where dΩ 2 (p−2) and dΩ 2 (q) = dθ 2 +sin 2 θ dΩ 2 (q−1) are the line elements of the round spheres appearing in the AdS p and S q factors respectively, and The field strength and radius of the S q coincide with those for the background AdS p × S q and are given in (1.5).
It is straightforward to compute the temperature and entropy of this black hole. The energy may be found by integrating the first law of black hole mechanics and then the free energy is given by F = E − T S. We find: is the vacuum energy of Anti-de Sitter space, which is non-vanishing for p = 3, 5 (it is zero for p = 4), and whose precise expression cannot be found from the first law, as it appears as an integration constant. One can find those 'Casimir' energies, for instance, with the prescription of [36] for defining the stress tensor of asymptotically AdS spacetimes. Here we instead use the results that follow from the KK holography analysis (see Appendix B and in particular equation (B.16)). Letting L be the AdS radius, Schwarzschild-AdS p × S q black holes are called small if r + L or large if r + L. We may invert the relation between the temperature T and the horizon radius r + in (2.3) to find the two branches of black holes at the same temperature T : The positive branch describes large AdS black holes, while the negative branch describes small AdS black holes. The temperature has a minimum and diverges in the limits of either large or small r + . This minimum temperature separates the two branches of black holes and there are no homogeneous static black holes in AdS p × S q (for fixed values of p and q) for T < T min . This is summarized in Fig. 2.1.
With the expressions for the physical quantities in (2.3), we can now review the phase diagram of the homogenous asymptotically AdS p × S q solutions. These include the background AdS p × S q itself and the two branches of Schwarzschild-AdS p × S q black holes. Since these spacetimes are static, we can trivially consider their Euclidean sections. There are two relevant thermodynamic ensembles: the microcanonical and the canonical ensembles. In the former, we compare solutions with the same energy and those with the largest entropy dominate; in the latter, the temperature is fixed and the solutions with the lowest free energy are preferred.
In the microcanonical ensemble the phase diagram is trivial: the background has no horizon and hence vanishing entropy; consequently, there is only one solution at a given energy, namely the corresponding Schwarzschild-AdS p × S q black hole. On the other hand, the phase diagram in the canonical ensemble is much richer. First notice that global AdS p × S q exists at any temperature: the period β of the Euclidean time circle can be chosen to have any value. The geometry has no horizon and therefore vanishing entropy. Consequently, F AdSp = E AdSp . From now on we subtract this contribution in (2.3), and consider the difference ∆F = F − F AdSp . For the Schwarzschild-AdS p × S q , regularity of the Euclidean section requires that the period of the Euclidean time circle is given by β = 1/T with T given in (2.3). Furthermore, from (2.3) we see that ∆F < 0 for r + > L, and hence large Schwarzschild-AdS p × S q black holes above this threshold dominate the canonical ensemble. On the other hand, ∆F > 0 for r + < L, and thus thermal AdS p × S q is preferred. At r + = L, corresponding to Hawking-Page (HP) temperature, there is a phase transition between two distinct equilibrium states: large Schwarzschild-AdS p ×S q black holes and thermal AdS p × S q . This is summarized in Fig. 2 Table 1: Relevant thermodynamic values in the phase diagram of the solutions preserving the background's isometries. In D = 10, the G 10 has been expressed in terms of the dual gauge theory parameter N using (4.1). In D = 11, G 11 is a fundamental constant and can be taken to be one.
yielding a non-vanishing latent heat ∆E, and the phase transition is of first order. For temperatures T < T min , there only exists thermal AdS. For T min < T < T HP , homogenous black hole solutions do exist but do not dominate the canonical ensemble; for T > T HP , large Schwarzschild-AdS p × S q black holes are the preferred phase, with both phases co-existing at T HP . Small Schwarzschild-AdS p × S q black holes are always subdominant. 5 So far, the physics only presents dependence on p. For the values (p, q) of interest, these quantities are shown in Table 1. Schwarzschild-AdS p × S q describes the universal sector to which Ref. [5] refers, only involving gravity in AdS and with the internal space playing no role. This is indeed the case if we only consider large black holes and the background. On the other hand, small enough Schwarzschild-AdS p × S q black holes are dynamically unstable under perturbations that break the isometries of the internal space [13]. Whilst the thermodynamic stability of the small black hole branch changes at becoming unstable for r + < r min + , this does not necessarily signal the onset of a dynamical instability. Indeed, the onset of the dynamic instability occurs when the Euclidean eigenmode λ E of the Schwarzschild-AdS p black hole coincides exactly with the 'momentum' in the transverse space. In the present case, where we take the internal space to be a S q , the threshold unstable mode will occur when λ E = −ℓ(ℓ + q − 1), where ℓ is labels the spherical harmonics on the S q . This is how the zero mode that signals the instability was first obtained. In turn, it suggests the existence of new branches of black holes that emerge from this point and that have deformed horizons along S q . These new branches of black holes are completely analogous to the nonuniform strings in the context of the GL instability of black strings [15,16]. From the Ddimensional point of view, for each ℓ there is a GL zero mode that breaks the SO(q + 1) isometry of S q down to SO(q). These zero modes can be uplifted to fully non-linear solutions, leading to the so-called lumpy black holes, first constructed in [17] in AdS 5 × S 5 . The aim of this paper is to study in detail several branches of lumpy black holes in AdS 5 × S 5 and AdS 4 × S 7 .
Linear perturbations about the Schwarzchild-AdS p × S q spacetime can be decomposed into a sum of scalar, vector and tensor perturbations, depending on how they transform under coordinate transformations on the S q . The GL zero modes are in the scalar sector and are classified in terms of ℓ. Here, instead of proceeding as in [13] and expand the AdS sector of the metric perturbations in spherical harmonics, we will linearize the full D-dimensional general ansatz for the lumpy black holes (see the next section) around Schwarzchild-AdS p × S q . The problem of numerically finding the lumpy black holes requires to construct the linearized Einstein-DeTurck operator, 6 and hence we can easily adapt our code to find the spectrum of time independent perturbations of Schwarzschild-AdS p × S q black holes. Since the latter turn out to be thermodynamically unstable for r + < r min + , we started our search for zero modes at r + = r min + and gradually decreased this value until a zero mode is found. The first three zero modes, corresponding to ℓ = 1, 2, 3, occur at where y 0 = r + /L . In practice we use these zero modes as initial seeds for finding the corresponding non-linear families of lumpy black holes. We note that in the D = 10 case, our values coincide with those found in [17] for ℓ = 1, 2. In the D = 11 case, as far as we know, our results are new. The spectrum of negative modes of the Schwarzschild-AdS 4 black hole was computed in [14]. The threshold radii for the zero modes that we have found agree, within numerical error, with the values of the radii where these negative modes coincide with the eigenvalues −ℓ(ℓ + 6) of the scalar harmonics on the S 7 .
3 Lumpy black holes in AdS p ×S q In this section we present the numerical construction of lumpy black holes in AdS p × S q . The different cases that we explore differ by the equations of motion and by the n-form gauge field, but the ansatz for the metric is general and can be written down for any values of p and q. Lumpy black holes in AdS p × S q only retain a SO(q) subgroup of the full isometry group of the internal S q . In the space of static solutions, they emanate from the small Schwarzschild-AdS p ×S q black hole branch at the GL threshold point, and therefore they have the same horizon topology: S p−2 × S q . They are classified by ℓ, which is the harmonic number corresponding to the scalar harmonic on the S q that labels the zero mode. Each lumpy solution of the ℓth family depends on two coordinates: the radial coordinante r in the AdS p factor, and the polar angle θ on the S q ; using standard spherical coordinates, the polar angle ranges from 0 to π, with the 'north' pole of S q at θ = 0, the equator at θ = π/2 and the 'south' pole at θ = π.
Uniform horizons on S q trivially possess a reflection symmetry θ ↔ π − θ with respect to the equator of this S q . However, under this parity transformation scalar harmonics pick up a factor of (−1) ℓ . Recall that the various families of lumpy black holes are non-linear solutions of the field equations that are continuously connected (in the space of solutions) to the linear deformations of the uniform black holes parametrized by the scalar harmonics. Therefore, the lumpy families labelled by the harmonic ℓ transform in the same way as the scalar harmonics under reflections about the equator of the internal S q . Consequently, lumpy black holes labelled by an odd ℓ are not symmetric with respect to reflections about the equator, and one has to consider the whole range of the polar angle θ on the S q to construct them. On the other hand, for even ℓ lumpy black holes enjoy a reflection symmetry about the equator and it suffices to consider θ ∈ [0, π/2].

Metric ansatz
Given the previous considerations, a suitable ansatz for the metric of the lumpy black holes consists of a general enough deformation of the Schwarzschild-AdS p × S q black hole that breaks the SO(q + 1) symmetry down to SO(q). To solve the equations of motion numerically, it turns out to be convenient to compactify the radial coordinate and to redefine the polar angle as follows: In these coordinates, the horizon is located at y = 0 and the conformal boundary of AdS p is at y = 1, while the north and south poles of the S q are at a = ±1 respectively, and the equator is at a = 0. Then for odd/even ℓ the integration domain is simply We will be parametrize the ansatz in terms of the dimensionless parameter y 0 = r + /L. If we further rescale the time coordinate t → Lt, then the AdS radius L drops out of the equations of motion, leaving y 0 as the only control parameter. This parameter effectively sets the temperature of the black hole. Therefore, for any given harmonic ℓ, the lumpy black holes form a one-parameter family of solutions specified by the temperature. Applying these changes of coordinates and redefinitions to (2.1), the metric reads: This is the metric that we will use as the reference metric in the Einstein-DeTurck equations (1.6). We can now present the ansatz for the lumpy black holes. For any value of p and q, we consider: where Q(y, a) ≡ {Q 1 , . . . , Q 6 }(y, a) are the unknowns for the metric field. For Q i (y, a) = 1, ∀i = 3 and Q 3 (y, a) = 0, we recover (3.2). These are actually the Dirichlet boundary conditions imposed at y = 1, as required by the asymptotics. The ansatz (3.4) preserves the full SO(p − 1) symmetry of the S p−2 within the AdS p factor and the unbroken SO(q) symmetry of the S q−1 inside the S q .

Gauge field
The different theories that admit AdS p × S q asymptotic solutions also contain a non-trivial nform gauge field strength; n = 5 in D = 10 and n = 4 in D = 11. In order to solve the equations of motion numerically, we need to provide a suitable ansatz for this field. Instead of considering the form field strength F (n) , we will work with the potential form field A (n−1) , such that F (n) = dA (n−1) . The reason is that the equations of motion for A (n−1) are of 2nd order and elliptic. To find a suitable ansatz, we consider a deformation of the gauge field A (n−1) corresponding to the Schwarzschild-AdS p × S q black hole (1.5) written in terms of the coordinates (y, a); see Eq. (3.1). In particular, we consider a gauge such that the (n − 1)-form potential vanishes at the horizon, which has the benefit of simplifying the boundary conditions there. Taking these considerations into account, we have the following ansatz for A (n−1) : where dσ (r) is the volume r-form of a round unit r-sphere, and Q 7 = Q 7 (y, a) is an unknown function that encodes the deformations away from Schwarzschild-AdS p × S q . For Q 7 (y, a) = 1 (and a particular function W for (p, q) = (5, 5)), these expressions yield the field strengths given in (1.5) in the (y, a) coordinates. Having fixed the gauge a priori, the equation of motion for the potential is elliptic; this equation together with the Einstein-DeTurck equations for the metric (1.6) and for suitable boundary conditions (see Section 3.3) form a well-posed boundary value problem.
There is a small caveat in the case (p, q) = (5,5). Notice that we have not specified the particular form of the function W (y, a) appearing in the ansatz. The reason is that the selfduality condition F (5) = ⋆F (5) allows us to eliminate the first derivatives of W in terms of the first derivatives of Q 7 . Then, the equation dF (5) = 0 can be written solely in terms of Q 7 and its first and second derivatives. This is not surprising, since the self-duality condition by itself is non-dynamical.

Boundary conditions
The boundary value problem has to be supplemented with suitable boundary conditions to form a well-posed problem. These are obtained by requiring regularity at the horizon and at the poles (or reflection symmetry at the equator for even values of ℓ), and by requiring that the solution is asymptotically AdS p × S q . In full, the boundary conditions that we impose are: • Horizon at y = 0: regularity at the horizon implies that all functions Q must be even in y and therefore we impose a Neumann boundary condition on all of them, except Q 3 which is Dirichlet. The condition Q 1 (0, a) = Q 2 (0, a) ensures that the geometry is free of conical singularities and fixes the surface gravity of the lumpy black holes to be the same as that of the reference metric.
• North and south poles of S q (for ℓ odd) at a = ±1 respectively: all Q's must be even in (1 ∓ a) respectively, and thus we impose Neumann boundary conditions on all of them except the crossed term Q 3 , which must satisfy a Dirichlet boundary condition, Q 3 = 0. In addition, to avoid conical singularities at the poles we impose: Q 4 (y, ±1) = Q 6 (y, ±1).
• Equator and north pole of S q (for ℓ even) at a = 0 and a = 1 respectively: all Q's must be even there and thus we impose Neumann boundary conditions on all of them except the crossed term Q 3 which must satisfy a Dirichlet boundary condition, Q 3 = 0. In this case, we only need to further impose the absence a conical singularity at the pole (a = 1) by requiring Q 4 (y, 1) = Q 6 (y, 1).

Further considerations
Now we are in position to construct the numerical solutions for the lumpy black holes. We start with the Schwarzschild-AdS p × S q black hole, corresponding to Q i (y, a) = 1, ∀i = 3, and Q 3 (y, a) = 0 in (3.4) and (3.5), and we add a bit of the ℓth zero mode, depending on the lumpy family that we aim to find, to construct a good enough initial guess of the Newton-Raphson loop. Whereas for ℓ odd, adding or subtracting the zero mode connects to the same lumpy family, this is not the case for ℓ is even. This is better understood in terms of the perturbations as written in (2.10): under a parity transformation, the perturbations transform as h ℓ µν → (−1) ℓ h ℓ µν . Therefore, we see that for odd ℓ, the positive sign in the perturbation is equivalent to the negative sign by transforming θ → π − θ, which is nothing but a redefinition of what we call 'north' and 'south' poles of the S q . However, this redundancy does not apply when ℓ is even, and the choice of the sign yields physically different perturbations that connect to different families of lumpy black holes. In this paper we have focused on finding the first four families of lumpy black holes, i.e., ℓ = 1, 2 + , 2 − and 3. The upper sign for ℓ = 2 refers to the family that branches off from the Schwarzschild-AdS p × S q black hole by adding (+) or subtracting (−) the zero mode for seeding the first solution.
Given the boundary conditions above, the temperature of our lumpy black holes is the same as the Schwarzschild-AdS p × S q black hole used as a reference metric and is controlled by the parameter y 0 in (3.4). We move along a given branch of lumpy black holes by varying y 0 . In analogy with non-uniform black strings in KK theory [16], we define a 'lumpiness' parameter λ (see Eq. (4.9)), that measures the size of deformations away from the Schwarzschild-AdS p × S q black hole. Slightly non-uniform (on the S q ) lumpy black holes have λ ≪ 1, while λ → ∞ as we approach the critical regime, where lumpy black holes merge with another family of black holes with a different horizon topology. For non-critical lumpy black holes (λ 1), one single patch of Chebyshev grid points is enough to obtain accurate numerical results. However, critical solutions need a more careful treatment. In order to access the critical regime we employed three techniques that have been used successfully in the past [23,24] . For critical solutions but not the most critical ones, we supplement the single domain with mesh-refinement in the region where it is necessary. As we approach the merger point, some functions develop very steep gradients near one (or both) pole(s) and mesh-refinement alone is insufficient to accurately construct such solutions; to resolve them, we simply redefine the unknowns as Q new = 1/Q old , for some functions (see below), and we split the domain into several (smaller) patches and apply mesh-refinement to increase the resolution where it is needed. The boundary conditions for the new unknowns are the same as the ones listed before, but the new functions remain bounded as we approach the critical regime.
In general, we have tracked the value of ξ 2 to determine when the resolution needs to be increased. We have always required ξ 2 < 10 −10 , despite the most critical solutions we reached do not satisfy this bound. Still, the numerical error is very localized and under control. In addition, we have also monitored some physical observables (see Section 4) along a given family; when numerical errors are unacceptably large they are reflected in the physical quantities. With these considerations, we have managed to confidently construct families of lumpy black holes up to following maximum values of λ: To construct solutions in the critical regime for the various families, we generically redefine the unknowns Q 1,2,7 as described above. For the families ℓ = 1, 2 + and 3, the singular behavior appears near a pole and near the horizon. To resolve it, we increase the resolution in that corner by splitting the integration domain into several patches, with higher resolution and applying mesh-refinement in the regions where it is needed. We have considered up to 4 patches for the most critical solutions. See Appendix E for more details and some examples of the types of domains that we have used to construct critical solutions.
It turns out to be much harder to construct critical solutions for the ℓ = 2 − branch. The reason is that the singular behavior of the functions is not localized at a corner but somewhere near the horizon and in between the equator and the poles. Additionally, for the two cases that we have considered, namely (p, q) = (5, 5), (4,7), there is a turning point along the branch in the non-critical regime. Indeed, we find that starting at y 0 = y ℓ=2 0 (see equation (2.11)) and decreasing y 0 , at some point we can no longer find solutions. Examining the behavior of the physical quantities near this point, it is clear that we are not at the end of the branch; in fact, we observe that the derivative of any thermodynamical variable as a function of the temperature becomes infinite. It turns out that this point corresponds to a local maximum of the temperature along the branch. To go past it and to larger values of λ, we use the trick explained in Ref. [35].

Thermodynamics
In this section we display and discuss the phase diagrams of asymptotically global AdS p × S q spacetimes, including the lumpy black holes solutions that we have numerically constructed. To make the microcanonical and canonical phase diagrams easier to visualize, we plot the dimensionless differences of the entropy or the free energy with respect to the small Schwarzschild-AdS p ×S q black hole with the same energy or temperature, depending on the ensemble.
In the D = 10 case we can use the holographic relation, to express the physical quantities in terms of gauge theory units. In this way, we shall display S/N 2 versus EL/N 2 in the microcanonical ensemble and F L/N 2 versus T L for the canonical ensemble. In D = 11, the Newton's constant is a fundamental quantity and it is only related to Planck's length. We can safely take units where G 11 = 1, and then S/L 9 , E/L 8 and F/L 8 are dimensionless. Before proceeding to display the phase diagrams, we collect the expressions for the various physical quantities of interest. Given our ansatz for the metric (3.4) and the boundary conditions that we impose, the temperature of lumpy black holes is the same as the reference Schwarzschild-AdS p × S q black hole and it is given by (2.3). The entropy of a black hole is proportional to its horizon area, S = A H /(4G D ). For generic values of p and q, the area of the horizon for the lumpy black holes is given by for odd ℓ. For even ℓ the integral ranges from 0 to 1 and we have to multiply by a factor of 2 due to the reflection symmetry about the equator of the internal S q . As a consistency check, note that setting Q 4 = Q 5 = Q 6 = 1, yields the horizon area of the Schwarzschild-AdS p × S q solution given in (2.3). The most straightforward way to compute the energy of the solutions is by integrating the 1st law of thermodynamics. We can easily find the curve S(T ), since it only involves the computation of horizon area, and the families of lumpy black holes are parametrized by T . Then the energy is given by: where E GL and T GL are the energy and temperature of the Schw-AdS p × S q black hole at the ℓth GL threshold point. This procedure, however, fails if there is a turning point along a given family, since at this point the derivative dS(T ′ )/dT ′ blows up. This is the case for the ℓ = 2 − branches. In this case, we can still integrate over the rest of the branch, where the first derivative is finite, and leave the contribution coming from the solutions around the turning point region as an undetermined constant. Then we can try to 'guess' such a constant using the phase diagrams: in the microcanonical ensemble the phase diagram must be smooth, whereas in the canonical ensemble the turning point must be a cusp. This trick works well as long as we have enough solutions sufficiently close to the turning from both sides. In both D = 10 and D = 11, the phase diagrams for the ℓ = 2 − families have been completed using this trick. An alternative way to compute the energy E of the solutions is by performing an asymptotic expansion of the fields near the AdS boundary and apply KK holography [29] and holographic renormalization [37]. In this paper we only have worked out this alternative way of finding E for D = 10 case. The calculation is long and tedious and the details can be found in Appendices A and B. The final expression is given by: Note that for Q 4 = Q 5 = 1, this expression reduces to the energy in (2.3) for (p, q) = (5, 5), after using (4.1) and identifying the well-known dimensionless vacuum energy of global AdS 5 [36]. The drawback of using (4.4) to compute E is that it implies computing four derivatives of the numerical solutions, which inevitably leads to a loss of accuracy.
(p,q) = (5,5) In Fig. 4.1 we display the phase diagrams for both the microcanonical and canonical ensembles for the ℓ = 1 lumpy family. The branch emanates from the GL threshold point, and it moves towards higher energies or lower temperatures. As this figure shows, the ℓ = 1 lumpy solutions are always subdominant in both the microcanonical and canonical ensembles. Furthermore, the phase diagrams do not present any turning points for the range of temperatures that we explored. As we shall see in Section 4.4, the local cone model for the merger predicts that in D = 10 there should not be any turning points sufficiently close to the merger point, i.e., in the limit λ → ∞. For this particular family it turns out that there are no turning points at all. The phase diagram for the ℓ = 1 lumpy AdS 5 × S 5 black holes was previously discussed in [17]. We now compare their results with ours. In the region near the GL threshold (i.e., λ 1), our results and theirs are in perfect agreement. However, we find some disagreement for λ ∼ O(1). According to Fig. 9 in their paper, the solution with the lowest temperature (i.e., the most critical) that they managed to find has T L = 0.49444 and energy EL/N 2 ∼ 0.388; using the data in their figure, we see that Min[R S 3 ] ∼ 0.09 and Max[R S 3 ] ∼ 0.56, which implies λ ∼ 2.6. On the other hand, the most critical solution that we have found for this branch has λ ≃ 25, which corresponds to a temperature T L = 0.49575 and energy EL/N 2 = 0.38540. Clearly our results and theirs are not compatible. The extent of the ℓ = 1 family in the phase diagram that we find is significantly shorter than previously reported [17]. The most reasonable explanation for this discrepancy is that the resolution used in [17] is insufficient to accurately construct the lumpy black holes with large λ. Indeed, we had to use various tricks (see Section 3.4) to be able to numerically construct critical lumpy black holes, whilst [17] used a single patch with moderate resolution. According to our estimates, the grid setup used in [17] rapidly becomes insufficient to accurately resolve the lumpy black holes with λ ∼ O(1). Therefore, it seems reasonable to conclude that part of the phase diagram reported in [17] is a numerical artefact and hence unphysical. In the canonical ensemble, the ℓ = 2 − lumpy black holes are also dominant with respect to the small Schwarzschild-AdS 5 × S 5 black hole for temperatures near the GL point, but thermal AdS 5 × S 5 is always the preferred phase in this range of temperatures.
The phase diagrams for the branches ℓ = 2 + , 3 are qualitatively the same as those for the ℓ = 1 family and we do not discuss them any further. We simply stress that these other families Phase diagram in the canonical ensemble containing all asymptotically AdS 5 × S 5 spacetimes that we consider in this paper. The GL threshold point for the various modes is indicated with a solid black disk. The zero of the free energy is taken to be the background thermal AdS 5 × S 5 , which has free energy F AdS = 3N 2 /(16L). All lumpy black holes are subdominant, even though the ℓ = 2 − family has lower free energy than the small Schwarzschild-AdS 5 × S 5 black hole for some temperatures.
of lumpy black holes are also always subdominant in both the microcanonical and canonical ensembles. In order to 'locate' in the phase diagram the various families of lumpy black holes that we have studied, in Fig. 4.3 we display the phase diagram in the canonical ensemble with all solutions that we have found. In this case we compare the free energy with respect to that of global AdS 5 × S 5 . As this figure illustrates, lumpy black holes occupy a very small portion of the entire phase diagram. Before we move on to discuss the phase diagrams for the ℓ = 2 − branch, we briefly comment on the previous results by [17] on the ℓ = 2 + and ℓ = 2 − branches. In particular, their most critical solution has temperature T L = 0.57278. We infer from their Fig. 11 (right panel) that this solution has λ ∼ 3.1. On the other hand, the left panel of their Fig. 11 indicates that their most critical solution for the ℓ = 2 − branch has temperature T L = 0.60849 and λ ∼ 0.5. In both cases, this agrees with our results and up to these solutions, our phase diagrams and theirs match very well. With our more accurate numerics, we have been able extend both the ℓ = 2 + and ℓ = 2 − families to larger values of λ, closer to the critical point.
It turns out that for the ℓ = 2 − branch, the solution with T L = 0.60849 is very close to a turning point, and more specifically to a local maximum of the temperature. This solution has λ ∼ 0.54, and hence it is still very far (in solution space) from the critical regime. We have managed to continue beyond this turning point, and the full phase diagram for the ℓ = 2 − family is displayed in Fig. 4.2. Note that whilst the ℓ = 2 + branch never dominates any ensemble, the ℓ = 2 − lumpy black holes can be favored with respect to the small Schwarzschild-AdS 5 × S 5 black hole. In the microcanonical ensemble, the ℓ = 2 − phase has larger entropy than the Schwarzschild-AdS 5 × S 5 black hole (for the same energy) in the region near the GL point. Likewise, in the canonical ensemble the ℓ = 2 − lumpy black holes dominate over the small Schwarzschild-AdS 5 × S 5 black hole near the GL threshold point and until slightly beyond the turning point. However, in this range of temperatures, thermal AdS is always the preferred phase in the canonical ensemble. It is interesting to note that for this particular branch of lumpy black holes, the microcanonical and the canonical ensembles seem to give slightly different results. This minor observation adds to the puzzles raised in [38] related to the apparent differences between the two ensembles in the thermodynamic limit.
(p,q) = (4,7) The phase diagrams in the microcanonical and canonical ensembles for ℓ = 1 lumpy black holes in AdS 4 × S 7 are depicted in Fig. 4.4. According to this figure, the ℓ = 1 lumpy phase exists beyond the corresponding GL threshold point towards larger energies or lower temperatures, depending on the ensemble under consideration. Just as in the (p, q) = (5, 5) case, this phase is always subdominant with respect to the small Schwarzschild-AdS 4 × S 7 black hole (and of course with respect to thermal AdS 4 × S 7 , which is the dominant phase in the canonical ensemble in this range of temperatures). It is interesting to note that even though we have managed to follow this branch up until λ = 104, and hence very close to the merger point, we did not find any turning points along the family. In the critical regime (i.e., large values of λ) the local cone model predicts the absence of turning points and our numerical construction confirms this model beautifully (see Section 4.4 for more details). The phase diagrams for the ℓ = 2 + , 3 lumpy AdS 4 × S 7 black hole families are qualitatively similar to those of the ℓ = 1 black holes and we do not present them here.
On the other hand, the phase diagrams for the ℓ = 2 − family exhibit some notable differences. In Fig. 4.5 we display the microcanonical and canonical phase diagrams for the ℓ = 2 − lumpy AdS 4 × S 7 black holes. Using the trick explained in the discussion after equation (4.3), we can find the complete phase diagram beyond the turning point and up to the most critical solution that we have been able to construct, even if we are integrating the 1st law. Fig. 4.5 includes the whole ℓ = 2 − family, up to λ = 6.00374 corresponding to T L = 0.37907.
In the microcanonical ensemble (left panel in Fig. 4.5), the ℓ = 2 − black holes dominate for energies near the GL point. This suggests that this family of black holes could be dynamically stable in some region of parameter space, and a certain (perhaps fine-tuned) class of perturbations of the Schwarzschild-AdS 4 × S 7 black holes could evolve into them at the non-linear level. There is a minimum energy along this branch at E min /L 8 = 519.30461; presumably, at this point the stability properties of these black holes change. The branch continues to higher energies again while the entropy decreases until it reaches E * /L 8 = 521.80046. At this energy, the lumpy black holes and the Schwarzschild-AdS 4 × S 7 black holes have the same entropy. Continuing along the ℓ = 2 − branch, for energies higher than E * , the Schwarzschild-AdS 4 × S 7 black hole is the phase with the highest entropy. The last solution we have found has energy 710.39314 × L 8 , corresponding to λ = 6, but the phase diagram extends way beyond this value.
In the canonical ensemble, right panel in Fig. 4.5, the ℓ = 2 − black holes dominate over the small Schwarzschild-AdS 4 ×S 7 near the GL threshold. There is a turning point along this branch at T turn L = 0.400824 (λ turn = 1.37424) corresponding to a local maximum of the temperature. The difference of free energies increases at lower temperatures beyond the turning point and it becomes zero at T L = 0.399124. At even lower temperatures, the small Schwarzschild-AdS 4 × S 7 dominates over the lumpy black hole phase. However, in the range of temperatures that the ℓ = 2 − black holes exist, thermal AdS 4 × S 7 is always the dominant phase in the canonical ensemble.
The phase diagram in the canonical ensemble showing all the black hole phases that we have considered in this paper is displayed in Fig. 4.6. As this figure shows, all families of lumpy black holes that we have studied are subdominant with respect to thermal AdS 4 × S 7 in the range of temperatures that they exist. However, the ℓ = 2 − branch locally dominates over the small Schwarzschild-AdS 4 × S 7 black holes for a small range of temperatures. It is clear from Fig. 4.6 that the families of lumpy black holes occupy a very small portion of the phase diagram.

Kaluza-Klein holography
An important class of observables that characterize the lumpy black holes are the vev's of the dual scalar operators of different conformal dimensions. These essentially parametrize the deformations of the internal space and they can be obtained using KK holography [29].   subsection we will focus on the (p, q) = (5, 5), i.e., type IIB supergravity in AdS 5 × S 5 , since in this case the details of this procedure have been fully worked out in [29]. We provide the details of our calculation in the Appendices A and B to correct some typos in the literature and because we find some disagreements with the calculation presented in [17]. The generalization to the (p, q) = (4, 7) is beyond the scope of this paper, even though there are some partial results in the literature [39,40]. The lumpy black holes in AdS 5 × S 5 have horizons that are non-uniformly smeared over the S 5 , breaking the SO(6) symmetry of the latter down to SO (5). Holographically, this corresponds to the spontaneous breaking of the R-symmetry in the dual N = 4 SYM theory and is reflected in the condensation of an infinite tower of scalar operators. The vev's of these operators encode the deformations of the internal S 5 and they can be computed, using KK holography, from the asymptotic expansion of the unknown functions Q 1 , . . . , Q 7 near the boundary of AdS 5 (see Appendix B). We find that the operators with the lowest conformal dimension that develop a non-trivial vev are S 2 , S 3 and S 4 . From (B.18), the vev's of these operators are given by: where β 2 , γ 3 and δ 4 are constants that can be obtained from our numerical solutions. On the other hand, the other family of scalar operators of low dimension that could have contributed, There are various ways to extract the coefficients β 2 , γ 3 and δ 4 (and δ 0 ) appearing in (4.5) from our numerical solutions. One possibility consists of taking derivatives of the Q's at the boundary and then project onto the harmonic basis. For instance, from (B.5) the constant β 2 is given by: where Y 2 (x) is the ℓ = 2 SO(5) scalar harmonic (see Appendix C) and the precise coefficient in front of the integral depends on the function Q that one uses to do the extraction. Another possibility to find the undetermined constants is to fit our numerical data according to the asymptotic expansion (B.5). This latter option has been shown to be more accurate and stable towards the end of the branches, specially for the coefficients that appear at 3rd and 4th order. The reason is that too close to the merger, the accuracy of the numerical solution deteriorates and we could not reliably compute the numerical derivatives needed to calculate the coefficients. In practice, we compute β 2 using both: (4.6) for all Q i 's, i = 1, 2, 4, 5, 6, 7 and performing data fits to check the consistency of the results and estimate the error in this quantity. The other coefficients, γ 3 and δ 4 (and δ 0 ), have generically been found fitting our numerical data. The temperature dependence of O i ≡ O S i , i = 2, 3, 4, along the ℓ = 1 lumpy black hole family is shown in Fig. 4.7 (right). The vev's vanishes at the GL threshold point, which lies on the Schwarzschild-AdS 5 × S 5 solution, and attains a finite value at the merger. The behavior of O 1 is consistent with the plot shown in Ref. [17], but O 3 differs in sign and value. We find that the behavior of the vev's along the ℓ = 3 branch is qualitatively the same as that for the ℓ = 1 family. As we show in Section 4.4, the approach of all three vev's to the merger is controlled by the local (singular) cone geometry in the bulk; in particular, it follows the same power law (with a logarithmic correction) and the same critical exponents as the other physical observables.
For the families with even ℓ we have only been able to accurately compute the vev O 2 and the result is shown in Fig. 4.8 (right). Again, it vanishes at the GL threshold point and becomes 7 From the expressions in equation (A.42) of [17], it would seem that O T 4 = 0. We believe that this is incorrect. In fact, we have checked from the near boundary expansions in equation (A.10) of [17] that the equations of motion for some of the gauge invariant scalar fields at linear and quadratic order are not satisfied. Therefore, it is possible that some of the expressions for the vev's given in [17] are not correct. positive along the extend of our branch. Unfortunately, in this case we could not find enough critical solutions to recover the critical exponents predicted by the local cone model.
We have also included an energy comparison when it is computed via integration of the 1st law and using the holographic expression (4.4). In Fig. 4.7 (left) and 4.8 (left) we present the energy as a function of the temperature for ℓ = 1, 2 + and 2 − lumpy black holes. For the ℓ = 1 and 2 + families we check the consistency of the energy values almost for the whole branch. For the ℓ = 2 − family we could only do it partially, but far enough to cover the turning point, which is the main interest in this case. Clearly, the turning point is well-resolved by the energy computed using the 1st law.

Geometry
In this subsection we study in detail the geometry of the horizon of the lumpy AdS p × S q black holes. We construct embedding diagrams of the horizon geometry into Euclidean space to visualize how the shapes of the various spheres change as one moves along a given family of solutions. This allows us to get a better understanding of the possible topology changes at the critical points, where the lumpy black holes merge with new families of black holes.
The induced metric on the spatial cross-sections of the horizon of the lumpy AdS p × S q black holes is: At the horizon, the radii of the S p−2 within the Anti-de Sitter factor and the radii of the S q−1 within the internal S q are given by: respectively. These are gauge invariant quantities within the symmetry class of spacetimes that we consider. The behavior of these two quantities (and in particular their zeroes) as functions of the polar angle a dictate the topology of the horizon. Indeed, the geometry (4.7) can be understood as a fibration of an S p−2 and an S q−1 over an interval parametrized by the coordinate a with a ∈ [−1, 1]. For all lumpy black holes, we see from (4.8) that R q−1 (a) smoothly shrinks to zero size at both ends of the interval whilst R p−2 is always finite there. This results in the familiar horizon topology S p−2 × S q , as the Schwarzschild-AdS p × S q black hole. In the latter case, both radii R p−2 (a) and R q−1 (a) are constant everywhere. As we shall see momentarily, for the lumpy black holes these radii develop a non-trivial profile as we move along a given family, away from the GL threshold. In particular, R p−2 (a) approaches zero at a point (or several points, depending on the harmonic ℓ labelling the family) on the horizon, becoming precisely zero at the merger. At the zero(s) of R p−2 (a), the horizon becomes singular and the branch of lumpy black holes presumably merges with a new branch of topologically distinct black holes in a conifold-type-of transition. In analogy with the non-uniform black strings [16], we may define the 'lumpiness' parameter in terms of R p−2 (a) as a measure of the (inverse) 'distance' to the merger: For solutions close to the threshold of the GL ℓth zero mode, R p−2 (a) is approximately constant on the horizon and consequently λ ∼ 0, with λ = 0 precisely at the GL threshold point.
Conversely, as we move along the family of lumpy black holes and towards the merger point, max[R p−2 (a)] remains finite whilst min[R p−2 (a)] → 0 and hence λ → ∞, just as in the standard KK setting. The vanishing of R p−2 at a certain point a = a * at the merger happens in a singular way. We should expect that at the other side of the transition, the vanishing of R p−2 occurs in a smooth manner, as it should if there exists a smooth family of black holes with a different horizon topology. We point out that even though the radius R q−1 (a) of the horizon S q−1 vanishes smoothly at the endpoints of the interval parametrized by a for all lumpy black holes, for some families (but not all) the horizon S q−1 can also become singular at the merger. However, for all cases that we have considered, the horizon S p−2 shrinks to zero size in a singular manner somewhere on the horizon precisely at the merger, and hence the lumpiness parameter defined in (4.9) is a good measure of the degree of deformation of the horizon geometry.
In Fig. 4.9 we plot λ and the minimum of R p−2 on the horizon, as functions of the temperature. This figure shows that λ defined in (4.9) behaves in the same way as the analogous quantity for the non-uniform black strings: λ increases monotonically along the family and diverges at the merger point. The behavior of λ for the other families of lumpy black holes labelled by different harmonics ℓ is qualitatively similar and we will not display them here. Note that for the ℓ = 2 − family, λ (and min[R p−2 (a)]) are not uni-valued functions of the temperature because there is a turning point at some temperature. However, the corresponding λ still increases monotonically along the family and it diverges in a qualitatively similar manner as in the other families as we approach the merger point.
One may also define the horizon length along the S q : (4.10) Again, for ℓ = 2 solutions the integral is from 0 to 1 and there is an extra factor of 2. The Schwarzschild-AdS p × S q black hole has L hor = πR, which is nothing but the geodesic distance between the poles of the S q of radius R. The behavior of this physical quantity along the family of lumpy black holes is quite generic. For a given ℓ, the curve starts at the GL point where L hor /L = π or L hor /L = 2π in D = 10 or 11 respectively, and it is a monotonically increasing function until it attains a finite value at the merger point. In Table 2 we list the different values at the merger for the families for which we have constructed critical enough solutions. As this table shows, the value of L hor /L at the merger point is unique to each family. Critical localizedtype black holes that merge with the various family of lumpy solutions studied here should give the same values of L hor /L. Ref. [17] considered the possible topologies of the localized black hole solutions that would merge with a given family ℓ of lumpy black holes in D = 10. This reference based their analysis on the behavior of the induced Ricci scalar at the poles of the S q , and concluded that the divergence of the Ricci scalar there indicated localization. In this paper, we follow a different route and instead we consider embedding diagrams into Euclidean space as they are more intuitive and the topology changes become apparent.
One may consider two possible embeddings: at fixed coordinates of the S p−2 or fixed coordinates of the S q−1 . For fixed S q−1 coordinates, the horizon geometry may be embedded into E p , whereas for fixed S p−2 coordinates the horizon geometry may be visualized as an embedded   Table 2: Dimensionless horizon's length values of the most critical solutions that we have found for the families ℓ = 1, 3 in D = 10 and ℓ = 1, 2 + , 3 in D = 11. We did not construct critical enough solutions of the ℓ = 2 + family in D = 10.
surface in E q+1 . Generically, the metric of E n reads The embedded surface has the form: X = X(a) and Y = R p−2 (a) or Y = R q−1 (a) for fixed S q−1 or S p−2 coordinates respectively. Then, matching (4.11) with (4.7), we obtain X(a): (4.12) Recall that for ℓ = 2 lumpy black holes the range of the integral is from 0 to 1, and to represent the full embedding we just add the mirror image with respect to the equator (a = 0). In our plots we present X/L versus Y /L for the three families ℓ = 1, 2 + and 2 − . For each of these families, we show together the embedding plots for fixed S p−2 and fixed S q−1 coordinates respectively, to make the topology of the horizon apparent. The qualitative behavior of the ℓ = 3 family is essentially the same as that of ℓ = 1; we comment on it below but we do not display any embedding plots for this case. In Figs. 4.10, 4.11 and 4.12 we show the embeddings for different families of the AdS 4 × S 7 lumpy black holes. The analogous embeddings for the AdS 5 × S 5 lumpy black holes are qualitatively similar and we display them in Appendix D. The only exception is the ℓ = 2 − family, which we comment on below.
In Fig. 4.10 we show representative embeddings for the ℓ = 1 AdS 4 × S 7 lumpy black holes for fixed S 2 coordinates (top) and fixed S 6 coordinates (bottom). This figure shows that whilst the S 6 smoothly shrinks to zero size at both poles of the S 7 , as λ increases it develops a conical singularity at the south pole. Similarly, the horizon S 2 has finite size everywhere on the S 7 for all lumpy black holes, but as λ increases it also develops a conical singularity at the south pole, suggesting that it pinches-off there. At the other side of the transition, both of these singularities will be resolved, resulting on a S 6 that smoothly shrinks to zero size at the north pole of the S 7 , whilst it has finite size on the south pole; conversely, the S 2 will continue to have finite size at with the poles located at a = 0 and a = π/2, and non-trivial smooth and finite functions f i (a), i = 1, 2, 3, such that f 3 (0) = f 1 (0) and f 2 (π/2) = f 1 (π/2). Indeed, this geometry has the topology of S D−2 and hence it corresponds to a single localized black hole, as expected. In Fig. 4.11 we consider the embeddings for different values of λ for the ℓ = 2 + family of AdS 4 × S 7 lumpy black holes. From the top panel in this figure, (i.e., fixed S 2 coordinates,) we see that the S 6 will develop conical singularities at the merger point at both the north and the south pole of the S 7 . Likewise, the bottom panel (i.e., fixed S 6 coordinates) shows that the S 2 will also develop conical singularities at both poles for λ → ∞. Then, we expect that the resolved geometry at the other side of the merger will have an S 6 that is finite in the whole a-interval, whilst the S 2 will smoothly shrink to zero size at both endpoints, giving an S 3 . Therefore, we conclude that the topology of the horizon of the family of black holes sitting at the other side of the merger is S 3 × S 6 , corresponding to a black belt.
The sequence of embeddings for the ℓ = 2 − family are shown in Fig. 4.12. For fixed S 2 coordinates (top panel), we see that see that the S 6 becomes more deformed, but it continues to shrink to zero size smoothly at the endpoints of the a-interval, even in the λ → ∞ limit. On the other hand, for fixed S 6 coordinates (bottom panel), we see that for λ → ∞ the S 2 will pinch off at two different points, a N * and a S * , near the north and the south poles of the S 7 respectively. Therefore, in the resolved geometry at the other side of the merger, we expect that the a-interval will be divided in three subintervals, −1 < a S * < a N * < 1, such that: For −1 ≤ a ≤ a S * , the S 6 shrinks to zero size smoothly at a = −1, whilst it remains finite at a = a S * ; on the other hand, the S 2 is finite at a = −1 and smoothly shrinks to zero at a = a S * . Therefore, the horizon topology in this region is S 9 . Likewise, a similar reasoning shows that the horizon topology in the a N * ≤ a ≤ 1 interval is also that of an S 9 . On the other hand, in the middle interval, a S * ≤ a ≤ a N * , we have that the S 6 is finite everywhere, whilst the S 2 smoothly shrinks at the endpoints, giving a horizon topology of S 3 × S 6 . Putting everything together, we conclude that the ℓ = 2 − family of lumpy black holes merges with a family of black holes whose horizon topology is a connected sum S 9 #(S 3 × S 6 )#S 9 , that is, two localized black holes with a black belt joining them. Interestingly, the ℓ = 2 − family in the AdS 5 × S 5 seems to merge with a different family of black holes. Indeed, as [17] pointed out and our embedding diagrams confirm (see Fig. D.1), in this case the lumpy black holes merge with a double localized black hole, with horizon topology S 8 #S 8 . It would seem that it should also be possible to connect the two localized black holes with a black belt, giving a topology S 8 #(S 4 × S 4 )#S 8 . However, if such a family of black holes exists, it does not merge with the ℓ = 2 − lumpy AdS 5 × S 5 black holes.
The ℓ = 3 solutions do not differ much with respect to the ℓ = 1 ones in terms of the topology change across the transition. Although the profile of the embeddings presents more oscillations (since it belongs to a higher harmonic family), our data seems to indicate that the branch merges with another localized black hole with S 9 horizon topology, but potentially having lumpiness along its horizon and hence distinct from the localized black hole family that merges with the ℓ = 1 lumpy branch. Presumably this new family of deformed localized black holes would merge with new families of multiple localized black holes and black belts. This is quite analogous to the situation in higher dimensional asymptotically flat singly spinning black holes, where different families of topologically spherical bumpy black holes, corresponding to different harmonics, merge with different families of black rings and bumpy black rings [28].
Summarizing, we have the following topology-change phase transitions for the different families of lumpy black holes considered in this paper: (4.14) Given the richness of the various topology changes that we have uncovered, it is not inconceivable that other topology changes are possible for values of (p, q) not considered in this paper. We see that except for the ℓ = 2 − family, the topology-changing transitions depend in a trivial way on the number of spacetime dimensions. The details of the underlying theory in each of the cases that we have considered, either type IIB supergravity or 11D supergravity, do seem not play a significant role. However, the particular way in which physical quantities approach the critical point does depend on the spacetime dimension, and is predicted by the local cone model governing the change of topology. Embedding plots for fixed S 2 /S 6 coordinates of lumpy black holes in AdS 4 × S 7 with ℓ = 2 − . This family merges with a black hole whose horizon topology is the connected sum S 9 #(S 3 × S 6 )#S 9 .

Critical behavior
Topology-changing transitions in the space of equilibrium solutions have been previously studied in higher-dimensional general relativity. The first and most studied system is the transition between non-uniform black strings and localized black holes in standard KK theory. Kol [20] proposed that, in this case, the transition is governed by a Ricci-flat double cone over S 2 × S D−3 . The S D−3 appears explicitly in both the non-uniform black strings and in the localized black holes; the S 2 on the other hand comes from the fact that the Euclidean time is fibered over an interval whose endpoints are on the horizon. The S D−3 is non-contractible in the black string phase, while it is contractible in the localized black hole phase; the opposite happens with the S 2 , which is contractible in the black string phase and non-contractible in the black hole phase.
There is now detailed numerical evidence that this model is correct in D = 5, 6 [22,23] and in D = 10 [24,25] (see [26] for an analytical approach in the large-D limit of general relativity). Beyond KK theory, Ref. [27,28] considered mergers of singly spinning black hole solutions in asymptotically flat space and described local Lorentzian double cone models for the critical geometries controlling them. These papers considered mergers between rotating black holes, black rings, di-rings and black saturns in D ≥ 6, providing strong evidence that such doublecone models do indeed describe the topology-changing phase transitions for asymptotically flat stationary black holes as well.
The topology-changing transitions that have been considered so far involve black hole spacetimes in pure Einstein gravity, and hence it is natural to consider local Ricci-flat cones to model the local singularity at the merger point. However, in this paper we consider Einstein gravity in asymptotically AdS p × S q coupled to a certain form field strength so it is not obvious that local Ricci-flat cones should still control the topology-change transitions in this case. In this section we will show that even though the topology change is still controlled by a conical geometry, the latter is no longer Ricci flat due to the presence of fluxes in the internal S q .
We will focus on the local merger between the ℓ = 1 lumpy black holes with horizon topology S p−2 ×S q and a localized black hole with horizon topology S p+q−2 , but the discussion also applies to the other ℓ's. Based on previous studies, the double-cone inherits a huge amount of spherical symmetry from the full solution: the S p−2 and S q−1 that appear explicitly in the metric, plus an additional S r that comes from the fibration, in the Wick rotated geometry, of the Euclidean time circle over some base manifold that depends on the particular merger. In the canonical nonuniform black string/localized black hole merger, the base is an interval with the two endpoints on the horizon; this, together with the fibered Euclidean time circle, gives r = 2; in the black ring/black hole merger the base is a 2-disk giving r = 3. In the case under consideration the base is also an interval. This consists of any path that connects two points on the horizon along the polar direction and thus r = 2, just as in non-uniform black string/localized black hole transition. Similarly, this S 2 is contractible in the lumpy black hole phase, but it is not contractible in the localized phase.
By symmetry, the merger should be described by a cone over a triple direct product of spheres, S 2 × S p−2 × S q−1 : With the presence of a non-trivial fluxes, it is not obvious that the cone that controls the merger should be Ricci flat and, as we shall see shortly, our data indicates that it is not. We will leave to future work the first principles understanding of the types of cones that control the topology changing transitions in the presence of fluxes. The numerical factors in (4.15) have been chosen so that the embeddings agree with our numerical solutions in the critical regime. The factors in front of the S 2 and the S p−2 in (4.15) ensure that the components of the Ricci tensor of the cone along these directions vanish; on the other hand, our data suggests that which implies that the components of the Ricci tensor along the S q do not vanish.
To test the validity of the local model (4.15), we compare the embeddings of the cone: Fixed Ω (q−1) : The plots on the left correspond to fixed S 6 coordinates while the plots on the right correspond to fixed S 2 coordinates. The black lines correspond to the embeddings of the local cone that mediates the topology change transition. For fixed S 6 coordinates, the agreement between the cone model and the local horizon geometry is very good. For fixed S 2 coordinates, demanding that the cone model agrees with the local geometry suggests α = q − 2 + D 4 = 31 4 . In these two plots the embeddings dictated by the Ricci-flat cone are indicated with dashed lines. Clearly, the latter do not agree with our data.
AdS 5 × S 5 one. For fixed S 6 coordinates, Fig. 4.13 shows that there is a very good agreement between the local cone model and the actual geometry of the horizon. On the other hand, for fixed S 2 coordinates along the AdS factor, we can demand that the local geometry of the horizon near the pinch off region is well described by (4.15) and thereby fix the constant α. By doing this, we find that (4.16) is preferred. The dashed line in this figure shows the embedding of a Ricci flat cone. Clearly, our numerical data is clearly incompatible with a Ricci flat cone. The local cone model not only captures the shape of the horizon near the pinch off, but it also predicts the behavior of the physical quantities near the merger [41]. The latter follows from considering perturbations of (4.15) that preserve the various spheres while they change their relative sizes: where ξ and η are free parameters. Imposing that the ρρ, S 2 and S p−2 components of the Ricci tensor of the perturbed cone (4.18) vanish implies that or Either choice leads to the same linear equation of motion for ǫ(ρ) and, unsurprisingly, it coincides with the equation of motion for the perturbation mode of the double-cone of [20]. Note that in our case, the first solution is the relevant one since it implies that the metric on the internal S q is not modified at linear order; this is consistent with the vanishing of the components of the Ricci tensor along the S 2 and the S p−2 directions. On the other hand, the components along the S q are affected by the leftover fluxes on this internal sphere, and hence they should be sourced at second order in perturbation theory. The general solution is given by: with critical exponents For D < 10, the solutions have a non-zero imaginary part, which implies that they will inspiral infinitely many times as they approach the origin ρ = 0; D = 10 is a marginal dimension for this model, where the imaginary part vanishes and the solution has a quartic power-law with an additional logarithmic term. For D > 10, the critical exponents become purely real and the approach to the origin is govern by two powers. In D = 11, we have s + = 6 and s − = 3. The prediction of [20] is that the deformations of the cone metric (4.21) should be reflected in the behavior of the physical quantities of the black holes sufficiently close to criticality, with the zero mode ǫ(ρ) measuring the deviation from the cone. According to this prediction, any physical quantity Q near the critical solution should behave as (p, q) = (5, 5) : , 7) : where δQ ≡ Q − Q c and ρ 0 is any typical length scale that measures the deviation from the actual cone. In principle, sufficiently close to the merger, the topology change develops in the same manner independently of the boundary data that one considers. To compute the critical exponents from our numerical data and compare them with the predictions of the cone model (4.23), we have fitted the temperature, energy, entropy and free energy of critical enough solutions according to (p, q) = (5, 5) : (p, q) = (4, 7) :  measuring the deviation from the cone. Notice that in D = 11 we have already assumed that one of the exponents is A + 3. This exponent corresponds precisely to the smallest contribution to f (x). In this case one would ideally want to fit the near critical data to a function of the form f (x) = f c + Ax B + Cx D , so that both exponents can be independently extracted and tested. However, in practice this turns out to be quite hard. The reason is that to carry out the fits, we only consider the solutions that are close to the merger and for those we have x ∼ 10 −3 . Therefore, the behavior of f (x) near the critical point is completely dominated by the smallest exponent, 3 in this case, and the corrections introduced by the other exponent (6 in D = 11) are too small to be reliably detected. In Figure 4.14 we present the fits for the free energy and the temperature (normalized with respect to the values at the marginal GL point) for the ℓ = 1 lumpy family in D = 10 and D = 11. The other physical quantities behave in a qualitatively similar way and we do not present the fits here. We have critical enough solutions for the families ℓ = 1, 3 in D = 10 and   Table 4: Critical exponents and other parameters obtained from the fits of the temperature, energy and entropy of the different lumpy solutions in AdS 4 × S 7 . Note that for a given ℓ, the critical exponents for the energy and the entropy are exactly the same. This is not surprising, since they are not independent: recall that in D = 11 the energy is found by integrating the 1st law of black hole mechanics. ℓ = 1, 2 + , 3 in D = 11. The fit parameters that we have extracted from the various physical quantities are shown in Tables 3 and 4 respectively. As these tables show, the values of the exponents obtained from fitting different physical quantities are perfectly compatible with each other and match the predictions of the cone model. With our fitting procedure we also extract the values of the different thermodynamic quantities (normalized with respect to the values of the ℓth family at the marginal GL point) at the merger point. In the different cases studied in this paper, the agreement with the cone model is very good, with deviations of 0.1% in the worst case. Working with a denser grid and/or including more data points, i.e., more solutions near the merger region, may improve the accuracy of the critical exponents. It should be clear from our studies that (4.15) does indeed capture both the local geometry of the horizon and the behavior of the physical quantities in the critical regime. In D = 10 we have also explicitly checked that the vev's of the dual scalar field S 2 , S 3 and S 4 , follow the same power law behavior near the merger as the other physical quantities. From the point of view of the cone model, this should not be surprising since any observable should exhibit the same behavior near the merger. However, from the point of view of holography, this is quite interesting. Firstly, the fact that the field theory observables, such as the vev's of the scalar operators, exhibit this scaling behavior indicates how the field theory knows about the topology change in the bulk: none of the field theory observables seem to become singular at the merger point; they simply develop a scaling behavior with exponents determined by the local cone model in the bulk. Secondly, one might naively expect that the vev's, being associated to scalar fields that are coupled to 5-dimensional gravity in AdS, develop a spiraling behavior since in D = 5 the critical exponents have a non-zero imaginary part. However, our result shows that this is not the case; the behavior near criticality is inherited from the full higher-dimensional solution.

Discussion and outlook
In this paper we have numerically constructed various families of asymptotically AdS p × S q black holes with a non-trivial field strength; we have considered the specific cases of (p, q) = (5, 5) corresponding to IIB supergravity in AdS 5 × S 5 , and (p, q) = (4, 7) corresponding to 11D supergravity in AdS 4 × S 7 . These black holes have horizon topology S p−2 × S q and were called 'lumpy' black holes in [17]. Each family is specified by an integer ℓ that labels the harmonic on the internal S q in the decomposition of the linear perturbations around the Schwarzschild-AdS p × S q black hole. For each harmonic, at a specific value of the horizon radius of the Schwarzschild-AdS p × S q black hole given in (2.11), there exists a time-independent linear solution (i.e., a zero mode) of the Einstein's equations corresponding to a black hole with a slightly deformed horizon on the S q . These linear modes can be uplifted to full non-linear solutions of the Einstein's equations, corresponding to the lumpy black holes. These zero modes signal the onset of the corresponding GL instability of the small Schwarzschild-AdS p × S q black hole.
Ref. [17] was the first one to construct lumpy black holes in AdS 5 × S 5 ; on the other hand, our solutions in the AdS 4 × S 7 case are new. Arguably, the solutions of [17] are only slightly non-uniform on the S 5 , with their non-uniformity parameter (4.9) being λ ∼ O(1) at most. In this article we have used similar numerical techniques as in our previous paper [24] to extend the various families of lumpy black holes to very large values of λ, all the way to near the critical regime where they are expected to merge with other families of black holes with different horizon topology. While our results and those of [17] agree very well near the GL threshold point (λ ≪ 1), we found some inconsistencies further along the branch towards more non-uniform lumpy black holes. We attribute this disagreement to the fact that the resolutions used in [17] were insufficient to accurately describe their most non-uniform black holes.
We have constructed the phase diagrams in both the microcanonical and the canonical ensembles all the way to the critical regime. In the microcanonical ensemble, we find that the lumpy black holes are always subdominant with respect to the Schwarzschild-AdS p × S q black hole with the same energy, except for the ℓ = 2 − family, which dominates in the region near the GL threshold point. In the canonical ensemble, for the range of temperatures that the lumpy black holes exist, thermal AdS p × S q is always the dominant phase whilst the lumpy black holes are always subdominant, even though the ℓ = 2 − family dominates over the small Schwarzschild-AdS p × S q black hole. The apparent disagreement between these two thermodynamic ensembles in AdS 5 × S 5 was highlighted in [38].
The main goal of this paper was to explore the topology-changing transitions that presumably connect, in the space of static solutions, AdS p × S q lumpy black holes with other families of black holes with different horizon topologies. These transitions have been well-studied in the past in the context of black holes in KK theory. Kol [20] predicted that such transitions were of the conifold-type and were governed by a local Ricci flat cone. According to this model, the local cone not only controls the change of topology, but it also determines the near-merger behavior of the physical quantities of the black holes on both sides of the transition. These predictions have been recently numerically verified [22][23][24]. In this paper we have extended these local cone models to black holes in AdS 5 × S 5 and AdS 4 × S 7 coupled to a (self-dual) five-form and a fourform field strength respectively. By constructing lumpy black holes with very large 'lumpiness' parameters (see (3.6)), we have verified, for the first time, that such local cone models also accurately describe the topology-changing transitions in AdS p × S q with the presence of nontrivial fluxes. Interestingly, we have found that the local cones (4.15) are no longer Ricci flat because of the presence of fluxes in the internal S q−1 in this limit. Work towards understanding the new class of conical geometries that arise when there are fluxes in the internal spheres is in progress.
For the asymptotically AdS 5 × S 5 lumpy black holes, we have extracted the dual field theory observables using the procedure of KK holography [29]. Our computation fixes some typos and inconsistencies in [17]. We give the expressions for the vev's of the scalar operators S 2 , S 3 and S 4 in (B.18); on the other hand, the vev's of T ℓ , ℓ = 0, 1, 2, 3, 4, vanish identically in these backgrounds (and in the localized AdS 5 × S 5 black holes). The vev's of the dual scalar operators exhibit the same scaling behavior near the merger as the other physical quantities, thus showing how these topology changes are 'detected' by the dual CFT. It is important to note that the behavior of the scalar operators near the merger is dictated by the full 10D geometry. The 11D solutions can also be dimensionally reduced on the S 7 . The KK map in this case, to leading and quadratic order, was worked out in [39], and [40] respectively. We will leave the computation of the holographic data for this case for future work.
The most natural extension of this work is to consider localized black holes in AdS 5 × S 5 and in AdS 4 × S 7 . Localized black holes in AdS 5 × S 5 have been previously constructed in [19], but their solutions are still quite far from the merger with the lumpy black holes considered in this paper. Using the same numerical techniques as in [24], we should be able to construct localized black holes in the critical regime and verify the predictions of the local cone model from the other 'side' of the transition. Work in this direction is in progress.
(Grant No. UF140319 and URF\R\201026). We would like to thank José Juan F. Melgarejo for useful comments on earlier drafts of this paper. B.C. would like to thank Universidad de Murcia for warm hospitality during the final stages of this work.

A Kaluza-Klein holography
In this appendix we review the procedure of holographic renormalization for solutions of IIB supergravity that are asymptotically AdS 5 × S 5 . This procedure was worked out in [29] (see also [42,43] for related work) and named, for obvious reasons, Kaluza-Klein holography, or KK holography for short. We review it here to set up our conventions and fix the numerical factors. In Appendix B we apply it to the specific class of solutions considered in this paper. These calculations were previously reported in [17]. However, since we found some disagreements with their results, we decided to write down our calculations in full detail so that they can be straightforwardly and independently checked, should anyone wish to do so. That said, we have verified that our calculations pass all consistency checks and hence we are quite certain of their correctness.
Let us begin recalling that in the context of gauge/gravity duality one can compute the vev's of the dual gauge invariant fields in the bulk by differentiating the (renormalized) onshell supergravity action once with respect to the sources, a procedure known as holographic renormalization [37] (see also [44,45]). In this way, one obtains the vev's of gauge invariant operators of the dual field theory in terms of certain coefficients in the asymptotic expansion of the bulk fields. Originally, this procedure was worked out for solutions that are asymptotically AdS. The extension to generic solutions that asymptote to the direct product AdS × X, where X is some internal (compact) manifold, requires to first carry out the dimensional reduction on X. In general, this dimensional reduction will involve an infinite tower of massive KK fields that come from the deformations of X, and which are dual to an infinite tower of gauge invariant operators in the CFT. Here we will concentrate on the AdS 5 × S 5 case, which is the relevant one for this paper.
Given that the information needed is encoded in the asymptotic behavior of the solution near the AdS boundary, one may work out the dimensional reduction analysis as a linear perturbation about the AdS 5 × S 5 background and then apply holographic renormalization. This is what the authors of [29] did in the most general way, using gauge invariant variables. This paper finds the effective 5-dimensional field equations and an explicit non-linear map between the higher-dimensional fields (which is the numerical data we have access to) and the solutions to these 5-dimensional equations. Integrating the latter equations one obtains a 5-dimensional (super)gravity action, which involves certain scalar fields; from here one can compute the 1-point correlation functions of the dual field theory applying the standard holographic renormalization prescription.
At this point the idea is clear and [29] provides an algorithm to implement it in practice. However, the expressions that arise are quite lengthly and involved. This is essentially because the dimensional reduction cannot be truncated to the massless sector when the isometries of the internal space are broken. In other words, the mapping between higher-and lower-dimensional fields, ψ k and Ψ k respectively, in general is going to be highly non-linear: where A i jk and B i jk are numerical coefficients. Luckily, it turns out that to find the vev of an operator of a given dimension, one only needs the map truncated at certain order in the number of fields [29]. Let us elaborate a bit more on this. If Ψ i is dual to an operator of dimension i, one would need to expand the field up to order Z i in the standard Fefferman-Graham (FG) coordinate to extract the operator's vev: If only a finite number of terms in the near boundary expansion are needed to compute the vev of (A.1), then the r.h.s. of this expression can only involve a finite number of fields. Recall that after all, the outcome of this mapping is directly related to a physical observable which must be finite. As we will see in the following sections and in Appendix B, the scalar operators whose vev's we shall compute have at most conformal dimension 4. Hence we will need to compute the near boundary expansion of the corresponding bulk fields up to order Z 4 in the FG radial coordinate. This truncates the mapping (A.1) to quadratic order in the number of fields.
In the rest of this appendix we will provide the details of this procedure in the context of this paper. This has been done before in [29] and in [17]. Here we repeat the analysis in these references because we use slightly different conventions and we take this opportunity to fix typos.

A.1 Kaluza-Klein dimensional reduction on S 5
The first step consists of dimensionally reducing the higher-dimensional solutions on the S 5 . To this end, we need to expand a general perturbation of AdS 5 × S 5 into a complete set of scalar, vector and tensor harmonics. Due to the symmetries of the spacetimes that we consider, we only have to consider scalar harmonics on the S 5 ; from these we can obtain the corresponding scalar derived vector and tensor harmonics. The details of the basis of harmonics that we use in this paper are summarized in the Appendix C.
We consider solutions that are fluctuations about the AdS 5 × S 5 background (denoted by the superscript "(0)"), The indices are split up into M = (µ, a), with µ ∈ AdS 5 and a ∈ S 5 . From the point of view of transformations on the S 5 , the perturbations of the metric components with two legs along the AdS factor behave as scalars, and the ones with one/two leg/s along the S 5 , behave as vectors/tensors. Then we may expand the metric as follows: Notice that the internal metric perturbation is split into the trace, which behaves as a scalar, and a traceless symmetric part which behaves as a tensor. Here we use the notation T (ab) to denote the symmetric and traceless part of a given 2-tensor T ab : T (ab) = 1 2 (T ab + T ba ) − 1 5 g ab T a a . The 5-form field strength is expanded in the following way: 8 where ω abcde denotes the volume 5-form of the S 5 and Λ ℓ = −ℓ(ℓ + 4) is the harmonic eigenvalue with respect to the scalar Laplacian on the S 5 . Given the self-duality condition of the 5-form, it turns out that the perturbation functions b νρστ (z) and b(z) are algebraically related, so we may just focus on one of them; b(z) for instance. A few comments are in order. First, in the expansions (A.5) we have only included fields that couple to scalar-derived vector and tensor harmonics. Although the analysis can be carried out in full generality, including transverse vector and tensor harmonics, the additional terms would break the SO(5) isometry that the class of spacetimes that we consider have. In addition, not all the expansions in (A.4) and (A.5) start with the lowest harmonic, ℓ = 0; there are a few exceptions. Given that Λ 0 = D a Y 0 (x) = D (a D b) Y 0,1 (x) = 0, then A 0 µ (z), B 0 (z), B 1 (z) and b 0 (z) are not defined, and the expansion for the perturbation fields h µa (z, x), h (ab) (z, x) and f abcde (z, x) starts at ℓ = 1, ℓ = 2 and ℓ = 1, respectively. The next step consists of inserting the expansions (A.4) and (A.5) into the linearized field equations, work out the mapping between the higher and the lower-dimensional fields, and obtain the 5-dimensional equations of motion for the dimensionally reduced fields. Before that, let us count the number of fields that we have. From the lower-dimensional perspective we have, for each value of ℓ, a graviton, h ℓ µν , a vector, A ℓ µ , and three scalars: B ℓ , C ℓ and b ℓ . However, some of these fields are pure gauge modes. By finding the gauge transformation rules, we can determine which ones are pure gauge and which ones correspond to the physical modes; a straightforward way to do this is to explicitly fix the gauge. For instance, Refs. [46][47][48][49] performed the dimensional reduction in the de Donder gauge, which can be shown to be equivalent to setting to zero the fields A ℓ µ and B ℓ . In the analysis we do here we follow [29], and these are not set to zero but instead they are used to define gauge invariant physical fields. This is much more convenient in practice since, in general, solutions to 10-dimensional supergravity equations will not be in the de Donder gauge.
Consider a gauge transformation: δx M = −ξ M . The shift up to linear terms in the perturbation fields is given by 9 Here the gauge parameter ξ M (z, x) is also expanded in harmonics, , we obtain the gauge transformations of the perturbations at linear order, and from these we can identify the gauge invariant combinations. At linear order these transformations can be directly read off by looking at the coefficients of each harmonic. At quadratic order, one needs to project onto the basis of spherical harmonics to identify the transformation rules.

A.1.1 Gauge invariance and equations at leading order
At leading order one obtains the following gauge transformations of the perturbations: Given that the different harmonics ℓ are decoupled at linear order, we analyze the cases ℓ = 0, 1 and ℓ ≥ 2, separately. From (A.9) one can easily identify the gauge invariant combinations, which we summarize in the following table: Fields Gauge invariant fields Table 5: Fields and gauge invariant combinations at linear order. For ℓ = 0 there is only one gauge invariant field and h 0 µν is a deformation of the AdS part of the background metric (since it carries trivial harmonic index). To construct gauge invariant combinations for ℓ > 0, one needs to introduce the auxiliary fieldÃ 1 Notice that the definition of the fieldb ℓ in Table 5 includes a factor 1/4 which is not present in [29]. This simply corresponds to a different normalization of the field strength: in our equations of motion (1.4) the coupling in front of the field strength is 1/96, whereas [29,[46][47][48][49] considered 1/6. To be consistent with the mapping, we need to rescale our field strength which amounts to the 1/4 factor in the dimensionally reduced b ℓ field.
To obtain the equations of motion for the gauge invariant fields, we plug the expansions (A.4) and (A.5) into the linearized type IIB supergravity equations of motion and project onto a given spherical harmonic. The equations for the fieldsh ℓ µν ,C ℓ norb ℓ are of the Klein-Gordon type, but these fields are not mass eigenstates. Diagonalizing, one finds the following gauge invariant combinations [46,47]: satisfying:

A.1.2 Gauge invariance and equations at quadratic order
At quadratic order the expressions are much more involved; for this reason we first discuss the scalar sector (B ℓ , C ℓ and b ℓ ) and then the tensor sector (h ℓ µν and A ℓ µ ). As mentioned earlier, at quadratic order the gauge transformations are obtained by projecting onto the basis of spherical harmonics. For instance, the transformation rule for the scalar field B ℓ is given by: with δh (ab) = 2D (a ξ P h b)P + ξ P D P h (ab) ; z(ℓ) and q(ℓ) are defined in Appendix C. Let us discuss again different values of ℓ separately. For ℓ = 1 we have two scalars, C 1 and b 1 , whose gauge transformations are the same as those at linear order, equation (A.9), and hence the gauge invariant quantities at quadratic order are the same as those in Table 5. For ℓ = 0 we only have one scalar, i.e., C 0 ; its gauge transformation is: 13) and the corresponding gauge invariant combination is given by: whereC ℓ andÃ ℓ µ are defined in Table 5. Notice that the sum above starts at ℓ = 1, but B 1 is not defined. We have written it this way just for notational compactness, and from this section onwards whenever we have B 1 we actually mean C 1 /Λ 1 . Here Ω 5 denotes the surface area of the unit round S 5 .
For ℓ ≥ 2 we have the whole scalar sector, B ℓ , C ℓ and b ℓ . At quadratic order the gauge transformations are given by with the notation for normalization factors and triple integrals defined in Appendix C: z ℓ 1 ≡ z(ℓ 1 ), a ℓ 1 ℓ 2 ℓ 3 ≡ a(ℓ 1 , ℓ 2 , ℓ 3 ), etc. From these transformations one infers the following gauge invariant fields: 10C (A. 16) In the tensor sector we are only interested in the massless KK graviton, i.e., h 0 µν (see [29] for the analysis of massive gravitons). At quadratic order, this field transforms according to: Then the combination 10 A few typos have been corrected with respect to the expressions given in [17,29].
where g (0) µν denotes the AdS 5 background metric, is not gauge invariant but provides a correction to the spacetime metric that transforms nicely. Notice that what plays the role of 'gauge invariant' field at linear order here is the combination:h µν , although it is not gauge invariant either (recall that δC 0 = 0 but δh 0 µν = D µ ξ 0 ν +D ν ξ 0 µ ). The shift inC 0 can be understood as a Weyl transformation that brings the 5-dimensional action to the Einstein frame. Thenh (0) µν satisfies the proper linearized Einstein equation [29].
From the gauge invariant scalar fields, one can define the same combinations as in (A.10), now replacing the 'bar' fields by the double 'bar' ones, which satisfy the modified equations of motion [47]: where the double 'hat' notation is to emphasize that these are defined in terms of gauge invariant fields at quadratic order (i.e., they are built up from double 'bar' fields). Note that the masses are the same as those of the linear fieldsŝ andt, equation (A.11). The box operator is again the D'Alembertian in AdS 5 and the Q's that appear in the r.h.s. of the equations are given by: and

A.1.3 Non-linear mapping: 5D action
If we pack together fields and masses:ψ I ℓ = {ŝ ℓ ,t ℓ }, m 2 I,ℓ = {ℓ(ℓ − 4), (ℓ + 4)(ℓ + 8)} (I = 1, 2), the scalar field equations at quadratic order (A.19) may be written as: for some coefficients D I ℓ 1 ℓ 2 ℓ 3 , E I ℓ 1 ℓ 2 ℓ 3 and F I ℓ 1 ℓ 2 ℓ 3 that are easily identified using the expressions in (A.20) (see Table 6). Notice that the scalar fields in the l.h.s. are the gauge invariant combinations at second order whilst the fields in the r.h.s. are the gauge invariant fields at linear order (given that the r.h.s. side is quadratic in the fields). In the case whereŝ 2 is the only non-zero field at linear order (as it is for our solutions), the gauge invariant combinations can be written in the form:ψ The four coefficients A I ℓXX are also given in Table 6. The higher derivative terms on the r.h.s. of equation (A.23) can be removed via a field redefinitionψ I ℓ → Ψ I ℓ = {S ℓ , T ℓ } of the form: Then the field equations for the fields Ψ I ℓ read All quadratic and cubic scalar couplings have been determined and can be found in the literature [46][47][48][49], but in our case it is sufficient to retain only the quadratic coupling to the fieldŝ 2 (or S 2 ). The numerical factors up to this quadratic coupling can be found in Table 6.
A very similar approach can be carried out for the massless KK graviton; the equation of motion forh 0 µν at quadratic order is sourced by higher derivative interactions [29] which, just as in the case of the scalars, can be removed through the field redefinition: The equations of motion for the scalar fields and the massless graviton after the field redef-   28). Notice that our numerical coefficients differ with respect to those given in [29]. More precisely, all the non-zero coefficients for ℓ = 2 or ℓ = 4 differ by the same factor, 4 2/5 and 7/5 respectively, compared to the corresponding coefficients in that paper. This is due to the fact that we use a different harmonic representation, which also includes odd values of ℓ. The coefficients are the same for ℓ = 0. For ℓ odd all coefficients vanish.
theory corresponds to the bosonic sector of D = 5, N = 8 gauged supergravity, whose action is given by: This yields the equation of motion for G µν : , (A.32) and the non-homogeneous Klein-Gordon equations (A.27) for the scalar fields.
Let us remark again that if one considered the problem in full generality, the resulting action would involve an infinite tower of KK fields. However, to find 1-point correlation functions we only need to expand the 5-dimensional fields near the AdS boundary up to a certain order. This truncates the infinite KK tower to a finite number of fields. When applying this procedure to our solutions we will see that the only field that contributes to the dual scalar vev's and the dual stress tensor is the field S 2 . This amounts to having the following non-zero potential in the action (A.30):

A.2 Holographic renormalization
Having a five-dimensional theory of gravity coupled to scalars in AdS allows us to apply the standard holographic renormalization prescription [37] to calculate the dual field theory observables. It is convenient to introduce the Fefferman-Graham gauge for the 5D solution, where the 4-dimensional AdS boundary B 4 is at Z = 0, and the additional coordinates X k run along B 4 . Solving the bulk equations of motion order by order in a near boundary expansion yields to a schematic form of the metric G ij (Z, X) and scalar field Ψ k (Z, X), where k is the dimension of the dual operator, ij (X) + . . . , The leading terms G ij (X), Ψ 2 (0) and Ψ k (0) parametrize the Dirichlet boundary conditions and correspond to the dual sources for the CFT stress-energy tensor, and the operators of dimension 2 and k respectively. This data, together with the normalizable modes G (4) ij ,Ψ 2 (0) and Ψ k (2k−4) , is sufficient to determine the series solution completely. Both pieces of data are related to the vev's of the 1-point correlators [29,44,45]: ij Tr(G (2) ) In (A.36), 'lower' indicates terms with index smaller than 2k − 4, and they only appear when the weight k of the operator can be written as a sum of weights of other lower operators. Given that in our case the dual operator with smallest conformal dimension is 2, this non-linear term vanishes for k = 3 and it is proportional to O 2 2 for k = 4. In the latter case, the proportionality factor can be computed for our harmonic representation using the expressions in [29]; see their equation (B.18).

B Harmonic expansion, stress-energy tensor and scalar vev's
In this section we perform the near boundary dimensional reduction for the class of spacetimes that we consider. We first expand our solutions near the AdS boundary in spherical harmonics on the S 5 according to (A.4) and (A.5). We work with the compactified radial coordinate y in our ansatz (3.4), and then transform to the FG gauge.
Our ansatz for the metric (3.4) involves products of the variables Q 2 , Q 3 and Q 4 , which are not convenient to perform the decomposition in spherical harmonics. For instance, the (aa)component is not just the polar component of the S 5 times the unknown function Q 4 , but it also includes a shift proportional to the product Q 2 Q 2 3 , coming from the crossed y − a term. The Q 2 3 factor in this term would introduce scalar-derived vector harmonic modes into the component h aa (y, a), which is not compatible with the decompositions in (A.4). We overcome this issue through a simple redefinition of variables: so that the new variables Q ′ 2 , Q ′ 3 and Q ′ 4 have a transparent decomposition in spherical harmonics. Then the spacetime metric reads 2) Given the asymptotic boundary conditions that our solutions satisfy, the unknown functions appearing in the metric (B.2) have the following near boundary expansion: (k = 1, 2, 5, 7), where the various q(i) ℓ 's are constants that can be determined by solving the equations of motion order by order near the AdS boundary (y = 1). The factors in the expansions for Q ′ 4 and Q 6 have been chosen such that the internal part of the metric can be written as where h (ab) and h are the symmetric-traceless and trace parts respectively. They have an expansion in harmonics as in (A.4), where now the coefficients of the harmonics are functions of the compact radial coordinate y. γ 1 , γ 3 , δ 0 , δ 4 and α 0 . Gauge freedom of both the metric and gauge field allows us to set γ 1 and α 0 to any value and we choose to set them to zero. Note that even if they are not set to zero, observables are gauge invariant and therefore, at the end, they cannot depend on these two constants. The rest of coefficients, {β 2 , γ 3 , δ 0 , δ 4 }, correspond to data that can only be determined from the full bulk solution.
To bring the metric into the FG gauge (z, θ) in 10D we perform the coordinate transformation y = y(z, θ), a = a(z, θ) in a series expansion near the AdS boundary z = 0. At each order in the FG coordinate z i , we have two functions to solve for and the only requirement is that the AdS part of the metric takes the FG form, i.e., g zz = 1 z 2 (recall that L = 1) and g zi = 0. This completely determines the form of the coordinate transformation. The first terms are given by: a(z, θ) = θ + 1 20 To obtain the metric up to order z 4 we need this change of coordinates up to order z 7 for y(z, θ) and up to order z 5 for a(z, θ). Now we can identify the different fields appearing in (A.4) and (A.5). 12 To this end, it is useful to write the background AdS 5 × S 5 in FG coordinates since we have to subtract it: and F (B.9) Finally, the 5-dimensional scalar fields arising from the 5-form are: (B.10) From these expressions we can find the gauge invariant combinations at linear order from Table 5, and then the corresponding gauge invariant mass eigenfields defined in (A.10). To leading order in the number of fields, the equations of motion are satisfied up to order z 2 , since higher order terms in z can (and will!) receive contributions from the non-linear interactions. The leading order fieldŝ 2 is given bŷ The rest of the fieldsŝ ℓ , ℓ = 3, 4 are at least O(z 3 ) andt ℓ , ℓ = 0, 1, 2, 3, 4, are O(z 4 ). We do not present the explicit expressions for them in this appendix, but they can be non-trivial at this order and they are needed to find the dimensionally reduced scalar fields S ℓ and T ℓ . The relevant fields to the second order areŝ 2 andh 0 µν . Using (A.24) and (A.18) we find: where δ ij = diag(1, σî), and {h ′ 0 , h ′ 1 } = {1/3, −61/45} for indices along the S 3 and one otherwise. We are interested in the scalar fields S ℓ , T ℓ given in (A.25) and the massless KK graviton G µν (A.29). These are 5D fields involving only quadratic contributions inŝ 2 and its derivatives. After adding the background in the case of the metric, we find: with η ij = diag(−1, σî), (G 0 , G 1 , G 2 ) = (−1, 15/31, −1/3) for the components along the S 3 and one otherwise. The fields T ℓ , ℓ = 0, 1, 2, 3, 4 are O(z 5 ). Therefore, up to order z 4 , there are only two other scalars which can acquire non-trivial vev's: (B.14) Having the 5D fields, we are now in a position to apply the standard holographic renormalization prescription (see §A.2). Clearly, the metric in (B.13) is not in FG gauge; introducing the 5D FG coordinate Z through z(Z) = Z + 13 86400 we can bring it to the canonical FG form. Finally, identifying the different terms that appear in the FG expansion, we obtain the stress-energy tensor from (A.37): The boundary geometry B 4 is conformally flat and therefore there is no conformal anomaly: the stress-energy tensor is traceless T i i = 0. It is also trivially conserved, D i (0) T ij = 0, and in the limit β 2 , δ 0 → 0, it reduces to the stress-energy tensor of the Schwarzschild-AdS black hole, as expected. The energy of the spacetime is simply given by E = Ω 3 T 00 , which includes the well-known Casimir energy of the dual N = 4 SYM on the ESU 4 . Using the expressions in (B.5) one has the relation β 2 which is very useful in order to check the numerics.
Recall that the conformal dimension of operators dual to S ℓ is ∆ + = ℓ (∆ − = 4 − ℓ); thus their vev's in our case are given by:

E Numerics and convergence tests
In this appendix we provide further details of our numerical construction of the lumpy black hole solutions and present some convergence tests.
Our numerical approach to construct critical solutions relies on the redefinition of certain unknown functions to ensure that the new functions have bounded values at the merger. However, steep gradients are unavoidable and therefore, very close to critical limit, these redefinitions are useless if they are not accompanied by a numerical grid that can properly resolve these gradients. To this end, it is convenient to divide the original grid into other small subregions in order to be able to increase the grid resolution just where it is necessary. Of course, it is not desirable to consider many subdomains as it may result in a loss of computational efficiency. A moderate number of patches together with mesh-refinement overcomes the efficiency issue while redistributing the density of points in a useful manner.
In this paper we have used up to 4 subdomains and a mesh-refinement function as in [24] mesh(X; A, B, in the patch where the singular behavior appears, which is typically at a pole in the near horizon region. Therefore, where needed, we redefine the coordinates asx = mesh(x; 0, x * , χ x ) and y = mesh(y; 0, y * , χ y ). In total, this introduces four new parameters: x * and y * , corresponding to the coordinate values where the different patches meet, and χ x and χ y , which control the 'strength' of the density of the new grid points along each coordinate direction. Typical values that we used are χ x , χ y = 4. At each shared boundary one must impose continuity of the functions and their first normal derivatives. An example of such grids is shown in Fig. E.2. As for the resolution used, for all families we began with a single patch with 60 × 60 grid points, and then we increased the number of patches up to four as we moved along the branch of solutions towards the critical regime. The smallest patch where the peaks of the functions typically appear has 50 × 50 resolution, whereas for the other patches lower resolutions suffice, e.g., ∼ 30 × 50 or 30 × 30.
To check whether the solution on a given grid is good enough, we monitor the numerical error. If this is greater than a few percent, then the resolution needs to be increased. To this end, we use the squared norm of the DeTurck vector ξ 2 and its non-zero components ξ x and ξ y since they should all vanish in the continuum limit and hence they are good measures of the error. In the case of using a pseudo-spectral numerical approach, as in this article, the error should be exponentially suppressed in the continuum limit, i.e., with increasing the grid size. Insufficient resolution can also reflect itself in some unphysical behavior of the physical parameters in the phase diagram of the family of black holes under consideration. For this reason, we also monitor the physical quantities along the branch.
In Fig. E.1 we compare the maximum values of ξ 2 and ξ y for a reference solution (with λ ∼ 0.2) at different resolutions in each of the branch of solutions we have constructed. To produce this figure we picked up a reference solution of each branch, we interpolated it at different resolutions and then we filtered through the Newton-Raphson loop. Solutions with other values of λ behave in a qualitatively similar manner. Of course, the larger the λ is the larger the error is with the same grid structure and resolution. Indeed, in all cases shown we observe that those quantities vanish in the continuum limit. Note that for the component ξ y , the values themselves, being coordinate dependent, are not physical. However, the fact that they become zero for large N confirms that we are finding an Einstein metric rather than a Ricci soliton. Small deviations from exponential convergence to the continuum in some cases can be explained by the interpolation errors induced at the boundaries Right: meshrefined grid pointsx = mesh(x; 0, x * , χ x ),ỹ = mesh(y; 0, y * , χ y ), with χ x,y = 4, in the blue patch. In the green and red patches only x-and y-directions need to be refined respectively.
of the patches and potential non-smoothness due to the deTurck gauge fixing. 13 Notice also that beyond a certain resolution the different curves flatten out due to the finite machine precision.