Critical behavior of the black hole / black string transition

TL;DR

This paper numerically studies the critical transition between localized black holes and black strings in higher dimensions, revealing a spiral phase diagram and complex critical exponents consistent with theoretical predictions.

Contribution

It provides high-precision numerical results on the phase diagram and critical behavior of black hole/black string transitions in five and six dimensions, extending previous analyses.

Findings

Phase diagram exhibits a spiral structure near the critical point.

Physical quantities follow a discrete scaling symmetry with complex critical exponents.

Numerical critical exponents match those from the double-cone metric theory.

Abstract

We numerically construct static localized black holes in five and six spacetime dimensions which are solutions to Einstein's vacuum field equations with one compact periodic dimension. In particular, we investigate the critical regime in which the poles of the localized black hole are about to merge. A well adapted multi-domain pseudo-spectral scheme provides us with accurate results and enables us to investigate the phase diagram of those localized solutions within the critical regime, which goes far beyond previous results. We find that in this regime the phase diagram possesses a spiral structure adapting to the one recently found for non-uniform black strings. When approaching the common endpoint of both phases, the behavior of physical quantities is described by complex critical exponents giving rise to a discrete scaling symmetry. The numerically obtained values of the critical exponents agree remarkably well with those derived from the double-cone metric.

Tables

Table 1. Table 1: Parameter values for mass, relative tension, entropy and temperature when fitted with ( 27 ).

$D=5$
localized black holes
	$f_{\text{c}}$	$a$	$b$	$c$	$d$
$M/M_{\text{GL}}$	1.6771933	2.4700	1.4997	1.9362	2.0766
$n/n_{\text{GL}}$	0.7795283	0.5762	1.4986	1.9359	4.2842
$S/S_{\text{GL}}$	2.6718298	7.3502	1.5001	1.9359	2.0752
$T/T_{\text{GL}}$	0.6738645	0.7869	1.4990	1.9367	5.3444
non-uniform black strings
	$f_{\text{c}}$	$a$	$b$	$c$	$d$
$M/M_{\text{GL}}$	1.6771932	0.7161	1.4995	1.9364	3.6215
$n/n_{\text{GL}}$	0.7795282	0.1691	1.5010	1.9375	5.8367
$S/S_{\text{GL}}$	2.6718297	2.1232	1.4994	1.9369	3.6237
$T/T_{\text{GL}}$	0.6738646	0.2295	1.4998	1.9358	0.6010

Table 2. Table 2: Parameter values for the horizon length when fitted with ( 29 ).

$D=5$
	$L_{\text{c}}$	$a_1$	$b_1$	$a_2$	$b_2$	$c_2$	$d_2$
$L_{\text{polar}}/L$	1.428268	0.5548	1.0024	0.7976	1.4976	1.9267	1.9116
$L_{\text{hor}}/L$	1.428265	0.2441	1.0041	0.2319	1.5046	1.9500	3.5614

Table 3. Table 3: Parameter values that we chose to construct the majority of the numerical localized black hole solutions (in particular those in the critical regime) with the choice of reference metric ( A.1 ). We used these values both for D = 5 𝐷 5 D=5 and D = 6 𝐷 6 D=6 .

$L$	${\varrho }_0$	${\varrho }_1$	${\varrho }_i$	${\phi }_i$	${\xi }_1$	${\xi }_2$
8	0.5	1.5	1	0.1	0	0.8

Equations

$$R_{\mu \nu }=0$$

$$\mathrm{d}{s}_{\mathrm{flat}}^2=-\mathrm{d}{t}^2+\mathrm{d}{x}^2+\mathrm{d}{y}^2+x^2\mathrm{d}{\Omega }_{D-3}^2.$$

$$\mathrm{d}{s}^2=-T_{\mathrm{a}}\mathrm{d}{t}^2+A_{\mathrm{a}}\mathrm{d}{x}^2+B_{\mathrm{a}}\mathrm{d}{y}^2+2F_{\mathrm{a}}\mathrm{d}x\mathrm{d}y+x^2{S}_{\mathrm{a}}\mathrm{d}{\Omega }_{D-3}^2.$$

$$T_{\mathrm{a}}=A_{\mathrm{a}}=B_{\mathrm{a}}=S_{\mathrm{a}}=1\mathrm{and}{F}_{\mathrm{a}}=0.$$

$$A_{\mathrm{a}}=S_{\mathrm{a}},F_{\mathrm{a}}=0\mathrm{and}\frac{\partial T_{\mathrm{a}}}{\partial x}=\frac{\partial A_{\mathrm{a}}}{\partial x}=\frac{\partial B_{\mathrm{a}}}{\partial x}=\frac{\partial S_{\mathrm{a}}}{\partial x}=0.$$

$$\frac{\partial T_{\mathrm{a}}}{\partial y}=\frac{\partial A_{\mathrm{a}}}{\partial y}=\frac{\partial B_{\mathrm{a}}}{\partial y}=\frac{\partial S_{\mathrm{a}}}{\partial y}=0\mathrm{and}{F}_{\mathrm{a}}=0.$$

$$x=\varrho \sin \phi ,y=\varrho \cos \phi .$$

$$\mathrm{d}{s}^2=-T_{\mathrm{h}}\mathrm{d}{t}^2+A_{\mathrm{h}}\mathrm{d}{\varrho }^2+{\varrho }^2{B}_{\mathrm{h}}\mathrm{d}{\phi }^2+2\varrho F_{\mathrm{h}}\mathrm{d}\varrho \mathrm{d}\phi +{\varrho }^2{\sin }^2\phi S_{\mathrm{h}}\mathrm{d}{\Omega }_{D-3}^2,$$

$$T_{\mathrm{h}}={\kappa }^2(\varrho -{\varrho }_0{)}^2{\tilde{T}}_{\mathrm{h}}$$

$${\tilde{T}}_{\mathrm{h}}=A_{\mathrm{h}},\frac{\partial {\tilde{T}}_{\mathrm{h}}}{\partial \varrho }=\frac{\partial A_{\mathrm{h}}}{\partial \varrho }=\frac{\partial }{\partial \varrho }({\varrho }^2{B}_{\mathrm{h}})=\frac{\partial }{\partial \varrho }({\varrho }^2{S}_{\mathrm{h}})=0\text{and}{F}_{\mathrm{h}}=0.$$

$$B_{\mathrm{h}}=S_{\mathrm{h}},\frac{\partial {\tilde{T}}_{\mathrm{h}}}{\partial \phi }=\frac{\partial A_{\mathrm{h}}}{\partial \phi }=\frac{\partial B_{\mathrm{h}}}{\partial \phi }=\frac{\partial S_{\mathrm{h}}}{\partial \phi }=0\mathrm{and}{F}_{\mathrm{h}}=0.$$

$$R_{\mu \nu }-{\nabla }_{(\mu }{\xi }_{\nu )}=0,$$

$${\xi }^{\mu }:=g^{\alpha \beta }({\Gamma }_{\alpha \beta }^{\mu }-{\bar{\Gamma }}_{\alpha \beta }^{\mu }).$$

$$f(x)\approx \sum _{k=0}^{N-1}{c}_k{T}_k\left(\frac{2x-b-a}{b-a}\right).$$

$$x_k=(b-a){\sin }^2\left[\frac{\pi k}{2(N-1)}\right]+a\text{with}k=0,1,\dots ,N-1.$$

$$T_{\mathrm{a}}\simeq 1-\frac{c_t}{x^{D-4}},B_{\mathrm{a}}\simeq 1+\frac{c_y}{x^{D-4}}.$$

$$M=\frac{L{\Omega }_{D-3}}{16\pi G_D}\left[(D-3)c_t-c_y)\right],$$

$$n=\frac{c_t-(D-3)c_y}{(D-3)c_t-c_y}.$$

$$A=2{\varrho }_0^{D-2}{\Omega }_{D-3}{\int }_0^{\pi /2}\sqrt{B_{\mathrm{h}}{S}_{\mathrm{h}}^{D-3}}(\sin \phi {)}^{D-3}\mathrm{d}\phi ,$$

$$R_{\text{eq}}={\varrho }_0\sqrt{S_{\mathrm{h}}}\text{at}(\varrho ,\phi )=({\varrho }_0,\pi /2),$$

$$L_{\text{polar}}=2{\varrho }_0{\int }_0^{\pi /2}\sqrt{B_{\mathrm{h}}}\mathrm{d}\phi \text{at}\varrho ={\varrho }_0.$$

$$L_{\text{axis}}=2{\int }_{{\varrho }_0}^{L/2}\sqrt{A_{\mathrm{h}}}\mathrm{d}\varrho \text{at}\phi =0$$

$$\mathrm{d}{s}^2=\mathrm{d}{X}^2+\mathrm{d}{Y}^2+X^2\mathrm{d}{\Omega }_{D-3}^2.$$

$$X(\phi )$$

$$Y(\phi )$$

$$\mathrm{d}M=T\mathrm{d}S+\frac{nM}{L}\mathrm{d}L.$$

$$(D-2)TS=(D-3-n)M.$$

$$f(Q)=f_{\text{c}}+a{Q}^b\cos \left(c\log Q+d\right),$$

$$Q_{\text{lbh}}$$

$$Q_{\text{nbs}}$$

$$L_{\text{H}}(Q)=L_{\text{c}}-a_1{Q}^{b_1}+a_2{Q}^{b_2}\cos \left(c_2\log Q+d_2\right),$$

Full text

1. 1 Introduction

2. 2 Physical and numerical setup

3. 2.1 Metric ansätze and boundary conditions

4. 2.2 Numerical implementation

5. 2.3 Physical quantities

6. 3 Results

7. 3.1 Thermodynamics

8. 3.2 Geometry

9. 3.3 Critical behavior

10. 3.4 Accuracy

11. 4 Conclusions

12. A Details on the numerical implementation

13. A.1 Choice of reference metric functions

14. A.2 Decomposition of the domain of integration

15. A.3 Parameter values

16. B Short review of the double-cone metric

1 Introduction

The study of black holes in higher dimensions $D>4$ has become a continuously growing topic over the last decades [1]. Considering $D$-dimensional asymptotically flat spacetimes with one periodic spatial dimension, e.g. topologically an $\mathbb{S}^{1}$, there are at least two different types of solutions: black holes and black strings. These two types can be distinguished by their horizon topology. Black strings wrap the compact dimension, and hence have horizon topology $\mathbb{S}^{D-3}\times\mathbb{S}^{1}$, while black holes are localized on the $\mathbb{S}^{1}$ and have horizon topology $\mathbb{S}^{D-2}$.

Uniform black strings in $D$ dimensions may be obtained by adding a compact, periodic dimension (an $\mathbb{S}^{1}$) to the solution of a ($D-1$)-dimensional Schwarzschild black hole.444The generalization of a Schwarzschild solution to higher dimensions is also known as Schwarzschild-Tangherlini solution. Gregory and Laflamme showed in the seminal papers [2, 3] that if the size of the compact dimension $L$ is large enough compared to the radius of the Schwarzschild black hole, small perturbations will grow exponentially, hence breaking translational invariance along the $\mathbb{S}^{1}$. In other words, assuming $L$ to be fixed, a uniform black string is stable for large masses $M>M_{\text{GL}}$ and unstable for small masses $M<M_{\text{GL}}$. Here, $M_{\text{GL}}$ denotes the mass at the Gregory-Laflamme point, where the uniform black string is marginally stable.

A few years later it was shown that a new branch of static solutions emanates from the Gregory-Laflamme point. Since these solutions naturally break the translation invariance along the compact periodic dimension, they were called non-uniform black strings. First, they were constructed only in $D=5$ by considering small perturbations around the Gregory-Laflamme point [4], but later this procedure was adapted to higher dimensions [5, 6]. Beyond the perturbative regime solutions could only be obtained numerically, which was done in a series of papers [5, 7, 8, 9, 10, 11, 12, 13], covering the dimensions from $D=5$ up to $D=15$.

Besides black strings there exist black hole solutions with a spherical horizon topology, thus not wrapping the $\mathbb{S}^{1}$. These solutions were first discussed in [14]. If the size of the event horizon of such a localized black hole555Sometimes these solutions are also referred to as caged black holes. is small compared to the size of the compact dimension $L$, the spacetime in the vicinity of the black hole is well approximated by a $D$-dimensional Schwarzschild solution. This fact was used to construct such solutions perturbatively in the regime of small masses [15, 16, 17]. Again, for solutions beyond perturbation theory we have to use numerics. This was done in $D=5,6$ [18, 19, 20, 21, 9] and very recently in $D=10$ [13].

Already in 2002, when numerical data for both the non-uniform black strings and the localized black holes was rare, Kol conjectured that both branches merge into a topology-changing singular transit solution [22]. Moving along the localized black hole branch this conjecture implies that the black hole will spread over the $\mathbb{S}^{1}$ until its poles merge and the compact dimension is completely wrapped. From the non-uniform black string side the transition is approached when the radius of the string at a certain point on the $\mathbb{S}^{1}$ shrinks to zero and the horizon pinches off. A central role in Kol’s conjecture plays the so called double-cone metric, which is supposed to be a local model of the transit solution at the critical point, where the merger or the pinch-off happens (depending on from which side the transition is approached).

All numerical results mentioned above are in accordance with Kol’s conjecture, though all implementations break down before the transition is reached. Furthermore, references [23, 8] provided numerical evidence that the spacetime of non-uniform black strings in $D=6$ dimensions may locally approach the double-cone metric in vicinity of the critical point.

Recently, non-uniform solutions were constructed far beyond previous results [12] and it was possible to investigate the regime of non-uniform black strings very close to the transition. Indeed the authors provided further evidence that the horizon of highly non-uniform black strings in $D=5,6$ converges towards an envelope close to the critical point, exactly described by the double-cone metric. Moreover, in this regime the associated thermodynamic quantities begin to oscillate with a rapidly decreasing amplitude. This led to typical spiraling behavior. In fact, evidence was presented in favor of a distorted logarithmic spiral, implying that not only its amplitude shrinks exponentially with each turn, but also its frequency. Consequently, there may be infinitely many oscillations before the pinch-off is reached.

In the present work we complement the results of [12] by constructing localized black holes in $D=5$ and $D=6$ dimensions. We are particularly interested in the regime, where the poles of the horizon are about to merge on the $\mathbb{S}^{1}$. We provide evidence that the associated thermodynamic quantities display similar behavior as in the non-uniform black string case. In particular, we find that these quantities are oscillating around their final values close to the critical solution. The corresponding spiral curves adapt extremely well to the ones of the non-uniform black string branch.

With the data of both branches at hand, for the first time, we are able to extract critical exponents of the thermodynamic quantities close to the localized black hole/non-uniform black string transition.666Reference [8] also tried to find critical behavior of non-uniform black strings but with hindsight the data was not close enough to the transition to observe the scaling. Moreover, the corresponding critical exponents are in remarkable agreement with the values that Kol derived from the double-cone metric [22, 24].

The paper is structured as follows: We introduce the physical setup for localized black holes in section 2. In addition, we outline the numerical implementation which is appropriate to obtain highly accurate solutions even close to the transition. In this critical regime the numerics are highly demanding, since the solution approaches a curvature singularity. The heart of our numerical implementation is a multi-domain pseudo-spectral method. In section 3 we present our main results and compare them to the non-uniform black string data from [12]. We conclude in section 4. Moreover, we provide supplementary material in appendix A concerning the numerical implementation. Finally, we briefly review the double-cone metric in appendix B.

2 Physical and numerical setup

In this paper we study solutions of Einstein’s vacuum field equations

[TABLE]

in $D$ dimensions approaching asymptotically $\mathbb{R}^{D-3,1}\times\mathbb{S}^{1}$ at spatial infinity. In other words, one of the dimensions is compactified to a circle of length $L$. One solution with the correct asymptotics is flat spacetime, given by the following line element

[TABLE]

The coordinate $y$ represents the compact dimension. We choose $y\in[-L/2,L/2]$ with period $L$. The $D-2$ infinitely extended spatial dimensions are expressed in spherical coordinates with $x\in[0,\infty]$ denoting the radial direction and $\mathrm{d}\Omega^{2}_{D-3}$ denoting the line element of a unit $(D-3)$-sphere.

In the following we construct localized black hole solutions approaching the flat metric (2) in the asymptotic limit $x\to\infty$. We simplify the problem as follows: First, we assume spherical symmetry of the spacetime with respect to the spatially extended dimensions, hence giving rise to the $(D-3)$-sphere discussed above. Second, we restrict ourselves to static solutions. Due to these two assumptions, our problem is effectively two-dimensional with non-trivial behavior in $x$- and $y$-direction. Third, we assume mirror symmetry of the spacetime with respect to the $y=0$ axis. Consequently, due to the periodicity in $y$-direction, we also have mirror symmetry with respect to $y=L/2$.

Localized black hole solutions have a spherical horizon topology, that is $\mathbb{S}^{D-2}$. Therefore, the horizon of such objects does not wrap the compact dimension. Let us place the center of the black hole at the origin of the $x$-$y$ plane. It is possible to choose coordinates such that the horizon has a spherical shape with radius $\varrho_{0}<L/2$. The corresponding domain of integration is depicted in figure 1 and contains five boundaries:777As a consequence of the mirror symmetry with respect to the $y=0$ axis we neglect the lower half of the domain of integration (where $-L/2\leq y\leq 0$) for convenience.

• •

the horizon $\mathcal{H}=\{(x,y)\colon x\geq 0\,,~{}y\geq 0\,,~{}x^{2}+y^{2}=\varrho_{0}^{2}\}$,

• •

the exposed axis of symmetry $\mathcal{A}=\{(x,y)\colon x=0\,,~{}\varrho_{0}\leq y\leq L/2\}$,

• •

the lower mirror boundary $\mathcal{M}_{0}=\{(x,y)\colon x\geq\varrho_{0}\,,~{}y=0\}$,

• •

the upper mirror boundary $\mathcal{M}_{1}=\{(x,y)\colon x\geq 0\,,~{}y=L/2\}$,

• •

the asymptotic boundary $\mathcal{I}=\{(x,y)\colon x\to\infty\,,~{}0\leq y\leq L/2\}$.

In addition we denote the $y=0$ plane as the equatorial plane and the point $(x,y)=(0,\varrho_{0})$ as the north pole of the horizon. Likewise, we refer to $(x,y)=(0,-\varrho_{0})$ as the south pole.

We proceed as follows in the remaining part of the section: in section 2.1 we discuss two different coordinate charts, which are adapted to different regions of the domain of integration. The overall numerical scheme is outlined in section 2.2, with more details provided in appendix A. Finally, we define the relevant physical quantities in section 2.3.

2.1 Metric ansätze and boundary conditions

The integration domain consists of five boundaries as depicted in figure 1. It is possible to cover the domain of integration with just a single coordinate chart, for example by constructing appropriate coordinate transformations, as done in [20]. However, note that such transformations have singular points and due to their complexity the resulting equations of motion will be lengthy.

Instead, we follow the approach of [9] (see also [19]) and introduce two different coordinate charts: one chart which is adapted to the asymptotic region and another chart which is adapted to the near horizon region.

The asymptotic chart

In order to describe the asymptotic behavior of the metric it is convenient to use the $x$-$y$ coordinates as evident from the line element (2). The most general ansatz which incorporates the required symmetries reads

[TABLE]

The five metric functions $T_{\mathrm{a}}$, $A_{\mathrm{a}}$, $B_{\mathrm{a}}$, $F_{\mathrm{a}}$ and $S_{\mathrm{a}}$ depend on $x$ and $y$. We impose the following boundary conditions for these functions on $\mathcal{I}$, $\mathcal{A}$, $\mathcal{M}_{0}$ and $\mathcal{M}_{1}$.

• •

The asymptotic boundary $\mathcal{I}$ ($x\to\infty$):

In order to approach the flat spacetime (2) asymptotically the metric functions have to satisfy

[TABLE]

• •

The exposed axis $\mathcal{A}$ ($x=0$):

At this boundary the metric degenerates. Nevertheless, the geometry has to be regular, i.e. there is no deficit angle at $x=0$. Therefore we impose

[TABLE]

• •

The mirror boundaries $\mathcal{M}_{0}$ ($y=0$) and $\mathcal{M}_{1}$ ($y=L/2$):

Since the metric should be symmetric with respect to those two boundaries we demand

[TABLE]

Note that the flat metric given by the line element (2) already satisfies the conditions above, a fact which will become important later. However, the flat metric obviously does not contain a horizon, where we have to impose additional conditions. Hence, we introduce another coordinate chart which is suitable for the horizon geometry.

The near horizon chart

The horizon is given as a circular contour in the $x$-$y$ plane as depicted in figure 1. Therefore, it is convenient to use polar coordinates in the near horizon geometry

[TABLE]

The metric reads in these coordinates

[TABLE]

where $T_{\mathrm{h}}$, $A_{\mathrm{h}}$, $B_{\mathrm{h}}$, $F_{\mathrm{h}}$ and $S_{\mathrm{h}}$ are functions of $\varrho$ and $\varphi$. Comparing equations (3) and (8) we conclude that $T_{\mathrm{h}}=T_{\mathrm{a}}$ and $S_{\mathrm{h}}=S_{\mathrm{a}}$, and that the functions $A_{\mathrm{h}}$, $B_{\mathrm{h}}$ and $F_{\mathrm{h}}$ are linearly connected to the functions $A_{\mathrm{a}}$, $B_{\mathrm{a}}$ and $F_{\mathrm{a}}$. We rewrite the function $T_{\mathrm{h}}$ as

[TABLE]

where $\tilde{T}_{\mathrm{h}}$ is regular at the horizon. This ensures that the black hole horizon is located at $\varrho=\varrho_{0}$. The relevant boundary conditions in the horizon chart are listed below.

• •

The horizon boundary $\mathcal{H}$ ($\varrho=\varrho_{0}$):

The metric has to be regular at the horizon and its surface gravity is given by $\kappa$. This leads to the following conditions

[TABLE]

• •

The exposed axis $\mathcal{A}$ ($\varphi=0$):

Regularity of the metric requires that

[TABLE]

• •

The mirror boundaries $\mathcal{M}_{0}$ ($\varphi=\pi/2$) and $\mathcal{M}_{1}$ ($\varrho\cos\varphi=L/2$):

It is straightforward to translate the conditions (6) into conditions for the functions defined in the near horizon chart.

2.2 Numerical implementation

Early numerical works approached the problem by solving Einstein’s vacuum field equations (1) directly [19, 20, 21]. Instead, we employ the well-established DeTurck method [9, 25] (see [26, 27] for reviews). In particular, we solve the Einstein-DeTurck equations (sometimes referred to as generalized harmonic equations)

[TABLE]

where $R_{\mu\nu}$ and $\nabla_{\mu}$ are the Ricci tensor and the covariant derivative associated with the target metric $g_{\mu\nu}$. The DeTurck vector field is defined as

[TABLE]

Here, $\Gamma^{\mu}_{\alpha\beta}$ is the Christoffel connection associated with $g_{\mu\nu}$, whereas $\bar{\Gamma}^{\mu}_{\alpha\beta}$ is the Christoffel connection associated with a prescribed reference metric $\bar{g}_{\mu\nu}$. One of the main advantages of the DeTurck method is that the Einstein-DeTurck equations (12) are strictly elliptic, in contrast to the original Einstein’s field equations (1).

Suppose that the metric $g_{\mu\nu}$ is a solution to the Einstein-DeTurck equations (12) and in addition the DeTurck vector field vanishes. Then $g_{\mu\nu}$ solves Einstein’s vacuum field equations (1). In particular, we have to ensure that the reference metric satisfies the same boundary conditions as the target metric. This is a necessary condition for a vanishing DeTurck vector field. This condition does not completely rule out the occurrence of so-called Ricci solitons, i.e. solutions to (12) which do not solve (1). However, as shown in [25], Ricci solitons do not exist in the static case considered here. Nonetheless, in practice we solve (12) numerically and we explicitly check afterwards if the DeTurck vector field decreases to reasonably small values.

In order to implement the DeTurck method, we have to find an appropriate reference metric $\bar{g}_{\mu\nu}$ which is consistent with the boundary conditions. For this purpose we follow the strategy outlined in [9] with some minor modifications. We use the flat metric (2) as the reference metric outside the ball of radius $\varrho_{1}$ centered in the origin of the $x$-$y$ plane, where $\varrho_{1}$ is chosen such that $\varrho_{0}<\varrho_{1}<L/2$. As stated beforehand, this metric satisfies the conditions on the outer boundaries. Inside this ball we have to use a reference metric that contains a horizon at $\varrho=\varrho_{0}$ with surface gravity $\kappa$. A natural candidate for such a spacetime is the $D$-dimensional Schwarzschild metric. However, this metric itself is not appropriate, since it does not match with the flat metric at $\varrho=\varrho_{1}$. Instead we modify the Schwarzschild metric appropriately. We refer the reader to appendix A.1 for a detailed description of the modification.

Note that the surface gravity of the reference metric has to be identical to the surface gravity of the desired target metric, when considering the ansatz (8) and the redefinition (9). Consequently, by varying $\kappa$ we construct physically inequivalent localized black hole solutions. In contrast, the parameters $\varrho_{0}$ and $\varrho_{1}$ can only change the gauge but they have no influence on physical properties of the solution, provided that they are chosen reasonably.

Having a reference metric at hand, we solve the Einstein-deTurck equations in the two different charts, that are discussed in the previous section. For this purpose, the domain of integration has to be divided in an asymptotic region, where the asymptotic chart is appropriate, and a near horizon region for which we use the near horizon chart. For convenience, we choose the boundary between these regions to be the contour $x=L/2$. In order to cover the whole domain of integration on a numerical grid, a coordinate transformation compactifying infinity $\mathcal{I}$ has to be performed in the asymptotic region. Moreover, we decompose the near horizon region, which has five boundaries, in such a way that each subdomain only has four boundaries. This is beneficial since there exists a coordinate transformation in each subdomain which is non-singular and covers the corresponding subdomain without an overlap to other domains. We further subdivide the integration domain using the special contour $\varrho=\varrho_{1}$ where the two different reference metrics have to be matched. The resulting basic arrangement of domains is depicted in figure 2. In appendix A.2 we give the corresponding coordinate transformations and discuss some further adaptions.

In the numerical implementation we omit one of the boundary conditions at the horizon $\mathcal{H}$ and at the exposed axis $\mathcal{A}$, respectively. This is due to the fact that there are six conditions for only five unknowns. In theory we are free to omit any of these conditions [27]. At the end, the numerical solution satisfies the missing condition, at least up to numerical precision. In practice we omit $\partial\tilde{T}_{\mathrm{h}}/\partial\varrho=0$ at the horizon $\mathcal{H}$ and $\partial S_{\mathrm{h}}/\partial\varphi=0$ at the exposed axis $\mathcal{A}$.

Finally, we discretize the resulting Einstein-DeTurck equations and boundary conditions by means of a pseudo-spectral method. The basis of our scheme is the expansion of any function $f\colon[a,b]\to\mathbb{R}$ in terms of Chebyshev polynomials $T_{k}(y)=\cos[k\arccos(y)]$, where $y\in[-1,1]$ and

[TABLE]

In the process, the function values are considered on Lobatto grid points which include the boundaries $a$ and $b$,

[TABLE]

We determine the spectral coefficients $c_{k}$ by utilizing a slightly modified version of the FFTW algorithm [28]. We apply the Newton-Raphson method for solving the discretized system describing the collection of non-linear partial differential equations (12). In the several iterative steps of this scheme a linear system involving a Jacobian matrix has to be solved. This is done by means of the so-called BiCGSTAB method [29] in combination with a preconditioner that utilizes a finite difference representation of the Jacobian. The sparse linear system arising within the preconditional step is solved efficiently (in terms of memory and time consumption) with the help of the SuperLU library [30, 31].

2.3 Physical quantities

In this section we determine the physical quantities of interest for us. In particular, we are interested in two asymptotic charges and several geometric quantities measured on the horizon $\mathcal{H}$ or on the axis of symmetry $\mathcal{A}$.

Asymptotic charges

We already stressed that in the asymptotic limit $x\to\infty$ the metric approaches flat spacetime (2). In the presence of a black hole the relevant leading order corrections to the metric at infinity are [32, 33]

[TABLE]

From the two coefficients $c_{t}$ and $c_{y}$ we can determine the mass of the black hole

[TABLE]

as well as the relative tension

[TABLE]

Here, $\Omega_{D-3}$ represents the surface area of a ($D-3$)-sphere and $G_{D}$ is the gravitational constant in $D$ dimensions.

Geometric quantities

According to the laws of black hole thermodynamics, the temperature $T$ of the black hole is proportional to the surface gravity, $T=\kappa/(2\pi)$, and the associated entropy $S$ is related to the horizon area $A$ by $S=A/(4G_{D})$. Note that we do not have to extract the surface gravity $\kappa$ from the numerical solution since we prescribe it by choosing the appropriate near horizon reference metric. We determine the horizon area $A$ within the near horizon chart (8) by computing

[TABLE]

where the metric functions are evaluated at $\varrho=\varrho_{0}$. Furthermore, we characterize the horizon with the following quantities: The horizon areal radius at the equator

[TABLE]

and the proper distance from north to south pole along the horizon

[TABLE]

The proper distance between the poles along the axis of symmetry

[TABLE]

is useful to distinguish between physically inequivalent solutions. For $L_{\text{axis}}\to L$ the radius of the black hole vanishes and hence we obtain flat spacetime (2). Moreover, we approach the transit solution in the limit $L_{\text{axis}}\to 0$.

For the purpose of illustration we embed the ($D-2$)-dimensional surface of the horizon into ($D-1$)-dimensional flat space

[TABLE]

Comparing this expression to (8) we deduce

[TABLE]

where both $X$ and $Y$ are evaluated at the horizon $\varrho=\varrho_{0}$. We have chosen an arbitrary constant in the integration of (24b) such that $Y=0$ at the equator in analogy to the coordinate $y$.

First law

For static black holes in spacetimes with one compact dimension the first law of black hole thermodynamic reads

[TABLE]

Since we keep the asymptotic length of the compact dimension $L$ fixed in our case, the first law reduces to $\mathrm{d}M=T\mathrm{d}S$. There is also an integrated form of the first law [32, 33], known as Smarr’s relation, which reads

[TABLE]

Note that equation (26) may serve as a non-trivial consistency check for a given numerical solution, since horizon values are related to the asymptotic value $c_{t}$.888Note that $c_{y}$ drops out of the right hand side of Smarr’s relation (26).

3 Results

In this work we use the previously described scheme to construct localized black hole solutions in five and six spacetime dimensions. The rather complicated domain setup allows us to increase the numerical resolutions where it is needed – in particular near the axis $\mathcal{A}$, the horizon $\mathcal{H}$ and the asymptotic boundary $\mathcal{I}$, see appendix A.2. Consequently, we are able to extend the branch of solutions far beyond previous results, while coming extremely close to the transit solution, where a transition to the non-uniform black string branch is supposed to happen. We have to increase the resolution considerably when approaching the transit solution to maintain accuracy. The computing time for a single numerical solution still remains small enough to be carried out on an ordinary PC within a couple of minutes.

First, we show the qualitative behavior of the thermodynamic and geometric quantities in sections 3.1 and 3.2. In particular, section 3.2 illustrates how the shape of the horizon changes when approaching the transit solution. In section 3.3 we present the main results of this work: We investigate the critical regime in more detail and find that the behavior of different physical quantities is governed by the same complex critical exponents, whose values were already conjectured [22, 24]. Due to the complex nature of the components we find evidence for a discrete scaling symmetry.

Moreover, in this section we compare the localized black hole solutions with the non-uniform black string solutions of [12]. In particular, we provide numerical evidence that indeed both branches are merging in the same transit solution. Finally, we discuss the accuracy of the localized black hole solutions obtained in this work in section 3.4.

3.1 Thermodynamics

In figure 3 we use the relative tensions $n$ to parametrize the solutions along both the localized black hole branch and the non-uniform black string branch. We normalize the relative tension and the other thermodynamic quantities by their corresponding values of a marginally stable uniform black string at the Gregory-Laflamme point, indicated by the subscript GL. The parametrization with respect to $n$ is convenient since the $n/n_{\text{GL}}=1$ line corresponds to the uniform black strings, from which the non-uniform black string branch emanates with values of $n/n_{\text{GL}}$ smaller than one.999All uniform black strings have the same tension, i.e. $n/n_{\text{GL}}=1$. The localized black hole branch starts where $n/n_{\text{GL}}$ is close to zero and the black hole is very small compared to the size of the compact dimension. By increasing the value of $n/n_{\text{GL}}$, the size of the localized black hole increases. However, as previous studies [21, 9] already pointed out, at some point the horizon area (or entropy) and the mass reach their maxima and then they decrease while heading towards the non-uniform black string branch.

The present study closes the gap between the localized black holes and black strings, as it is illustrated in figure 3. Following the localized black hole branch, we see that the first clearly pronounced turning point of the displayed thermodynamic variables is followed by at least three more. This leads to the spiraling behavior of the curves in figure 3. Moreover, these spirals adapt perfectly to those of the non-uniform black string branch. At each twist of the spirals the two phases approach each other. This strongly supports the expectation that both branches eventually end at the same point in the phase diagram.

In [12] it was argued for the non-uniform black string phase that the extend of the spirals shrinks exponentially with each turn (like for a logarithmic spiral), which led to the conjecture that there will be infinitely many turns before the endpoint is reached. This will be confirmed in section 3.3.

For completeness we determine the phase diagram in the microcanonical ensemble. In figure 4 we display the entropy $S$ as a function of the mass $M$. The thermodynamically preferred solutions for a given mass $M$ are the ones with largest entropy, i.e. for small masses the localized black hole solutions and for large masses the uniform black strings. In the inset of figure 4 we magnify the region around the transit solution and show that the localized black holes turn into non-uniform black strings.

3.2 Geometry

In this section we discuss geometrical aspects of the transit solution between non-uniform black strings and localized black holes and we provide additional evidence in favor of the conjectured double-cone metric. Figure 5 shows the behavior of some geometric quantities as a function of the relative tension $n$. In particular, we display the horizon areal radius at the equator $R_{\text{eq}}$ (20), the inter-polar distance along the axis of spherical symmetry $L_{\text{axis}}$ (22) and the polar distance along the horizon $L_{\text{polar}}$ (21) for the localized black hole solutions. In addition, we also show the corresponding quantities of black strings, i.e. the maximal and minimal horizon areal radius on the $\mathbb{S}^{1}$, $R_{\mathrm{max}}$ and $R_{\mathrm{min}}$, and the proper length of the horizon along the $\mathbb{S}^{1}$, $L_{\text{hor}}$.

Figure 6 depicts the spatial embeddings of the horizon (see (23) and (23)) for several localized black holes and non-uniform black strings. The horizons of the localized solutions spread more and more along the compact dimension until the poles are about to touch the boundaries of the compact dimension, and hence each other due to the periodic nature of the compact dimension. From the point of view of the black string, its horizon shrinks at the boundaries of the compact dimension until it is about to pinch-off. Note that the solutions which are closest to the phase transition, corresponding to embeddings four and five in figure 6, can hardly be distinguished. We refer to the region where the poles of the localized black hole are about to merge or the black string is about to pinch-off as the critical region.

In fact, we can only observe a difference between black strings and black holes in the critical region if we magnify it as in figure 7. Both, the localized black hole horizons and the non-uniform black string horizons, locally converge to straight lines, when approaching the transit solution.101010For non-uniform black strings this was already shown in [12]. This was part of the conjecture of Kol, namely that the two branches meet at a singular transit solution, which is locally given by the double-cone metric [22], see also appendix B for a short review. The predicted $D$-dependent opening angle of the horizon of the transit solution, see equation (44), is nicely approached from both types of solutions. However, there seems so be a qualitative difference in how the double-cone geometry is approached. In $D=6$ the double-cone geometry appears to be an envelope for both the localized black hole as the non-uniform black string horizons, whereas for $D=5$ this is only true for the black strings. The horizons of the localized solutions in $D=5$ are slightly crossing the double-cone horizon.

3.3 Critical behavior

Another interesting prediction, which follows from the analysis of the double-cone metric, is a critical behavior when approaching the transit solution [22, 24], see appendix B for a short review. In particular, different physical quantities may be expressed by the same critical exponents, which only depend on the dimension $D$. Originally, the critical behavior and their corresponding exponents appear in quantum and statistical field theories close to phase transitions. However, such a behavior was observed in some gravitational systems as well. The most famous example appears at the threshold of black hole formation in the context of a collapsing scalar field [34]. Surprisingly, there seems to be a relation between this system in $D-1$ dimensions and the black hole/black string system in $D$ dimensions [24].111111See [35] for an investigation of the spherical collapse of a scalar field in higher dimensions.

In order to show such a critical behavior in the black hole/black string transition, we fit the data of different physical quantities close to the transition with an appropriate ansatz, that is

[TABLE]

with the free parameters $f_{\text{c}}$, $a$, $b$, $c$ and $d$, cf. equation (49) in appendix B. The function $f$ denotes any physical quantity of interest, such as mass, entropy, temperature or relative tension. Furthermore, $f$ depends on $Q$, which uniquely parametrizes the localized black holes or non-uniform black strings. We choose $Q$ such that the critical transit solution is given by $Q=0$. For localized black holes, $Q$ is the inter-polar distance along the exposed axis, $L_{\text{axis}}$, and for non-uniform black strings we identify $Q$ with the minimal horizon areal radius $R_{\mathrm{min}}$. To be more precise, we define

[TABLE]

where $R_{\text{GL}}$ is the horizon areal radius of a uniform black string at the Gregory-Laflamme point. Due to the chosen normalization, $Q$ approaches one for the starting point of the corresponding branch, i.e. for an infinitesimal localized black hole or a marginally stable uniform black string.

In equation (27) $f_{\text{c}}$ denotes the value of $f$ at the critical transit solution, i.e. for $Q\to 0$, while $b$ and $c$ are the real critical exponent and log-periodicity, respectively. Following [22, 24] the values of $b$ and $c$ can be derived from the double-cone metric and hence should be the same for both branches and for different physical quantities $f$, see appendix B, and in particular equations (47) and (49) therein. The predicted values are $b=3/2$ and $c=\sqrt{15}/2\approx 1.9365$ for $D=5$ and $b=c=2$ for $D=6$, respectively.

We used Mathematica’s fit routine to obtain the values of the free parameters $f_{\text{c}}$, $a$, $b$, $c$ and $d$. As an example, in figure 8 we show data points and the corresponding fit for the mass $M$ (normalized by its value of a marginal stable uniform black string, $M_{\text{GL}}$). The left column displays the actual function $M(Q)$ and we see that the fit indeed approaches the data points for small values of $Q$. Since the amplitude and the period of the oscillations in (27) are decreasing exponentially for decreasing $Q$, a rescaling of the axes is performed in the right column of figure 8. This shows more explicitly the oscillating behavior and the good agreement of data points and fit for small $Q$.

We list the numerical values for the fit parameters of mass, relative tension, entropy and temperature in table 1. The predicted values of $b$ and $c$ are excellently reproduced, with deviations of less than $0.5\%$, and, moreover, for a given dimension they are the same for different quantities and for both branches. Furthermore, the respective values $f_{\text{c}}$ of the critical transit solution coincide for both branches up to 7 digits after the decimal point. This is therefore the best estimation of their actual values so far.

Let us stress that not all of the data points were used for the fit. Of course, for $Q\approx 1$ the quantities do not behave like the fitting function (27) and so only data points with reasonably small $Q$ are considered for the fitting procedure. To obtain the values listed in table 1 we only take the data points of about the last cycle into account. Extending the fitting range by including more data points results in slightly bigger deviations of the parameters $b$ and $c$ from the predicted ones, as expected since the ansatz (27) is only valid for small $Q$. The standard error of the fit routine for each parameter is of the order of the last digit (or even smaller) we gave in table 1, respectively.

Using the values of table 1 we are also able to perform some consistency checks. For example, we checked that the critical values of the thermodynamic quantities $f_{\text{c}}$ indeed satisfy Smarr’s relation (26) to the order of $10^{-7}$. Moreover, due to the first law of black hole thermodynamics (25), the extreme points of mass and entropy have to coincide. Given the ansatz (27) this implies that the phase shift $d$ should be the same for mass and entropy which is indeed the case with deviations of less than $1\%$, see table 1. Considering only the lowest order approximation (27) for the thermodynamic functions there are three further conditions on the parameters from Smarr’s relation and the first law, and they are also satisfied to a similar accuracy.

In addition, we are able to definitely answer the question raised in [9] whether there are localized black hole solutions with positive specific heat. Reference [9] provided evidence in favor of positive specific heat close to the first maximum of mass in figure 3. Note that the turning point of the mass does not coincide with the first minimum of the temperature as a function of relative tension, hence giving rise to positive $\partial M/\partial T$ between these points, in agreement with [9].

Furthermore, we find evidence for infinitely many regions with positive specific heat. According to table 1 the values of the phase shift $d$ of mass and temperature differ from each other (modulo $\pi$). Note that this is not a numerical discrepancy, e.g. compare with the good agreement of phase shifts of the observables mass and entropy. Thus, we conclude that in each of the tiny intervals between the corresponding turning points of mass and temperature the specific heat is positive. Of course, this argument holds also for the non-uniform black string branch.

From figure 5 it becomes clear that there is at least one quantity at both sides of the transition for which the ansatz (27) is not suitable even in the critical regime. It is $L_{\text{polar}}$ for localized black holes and $L_{\text{hor}}$ for non-uniform black strings. For convenience let us denote both as the horizon length. Obviously, near the transition there are no (leading) oscillations of the horizon length. Nevertheless, we tried to fit the horizon length as a function of the parameter $Q$ by adding a non-oscillating leading term to the ansatz (27), that reads

[TABLE]

where $L_{\text{H}}$ represents either $L_{\text{polar}}$ or $L_{\text{hor}}$ (normalized by $L$). This ansatz has seven parameters $L_{\text{c}}$, $a_{1}$, $b_{1}$, $a_{2}$, $b_{2}$, $c_{2}$ and $d_{2}$. Apparently, $L_{\text{c}}$ is the horizon length of the critical transit solution and, of course, it should take the same value no matter from which branch the critical solution is approached. Note that $b_{1}<b_{2}$ since by definition $Q^{b_{1}}$ is the leading term in the ansatz (29) for small $Q$.

The resulting fit functions are again in excellent agreement with the actual behavior of the horizon length in the critical regime as evident from figure 9. In the left column the leading behavior is illustrated by a double logarithmic rescaling. We see that in these diagrams the functions approach straight lines, with small wiggles on top of them. These are caused by the subleading term of (29). To show the oscillating nature of the subleading term more explicitly we rescaled the plots appropriately in the right column of figure 9.

We list the numerical values of the obtained fit parameters in table 2. Interestingly, the value of $b_{1}$ is approximately one for both $D=5$ and $D=6$. In other words, the horizon length is proportional to $Q$ to first order in the critical regime. Moreover, the values of $b_{2}$ and $c_{2}$ are in agreement with the critical exponents $b$ and $c$ of the thermodynamic quantities.

3.4 Accuracy

In the previous section we have seen that the predictions for the critical exponents are confirmed with deviations of less than $0.5\%$, which provides evidence that the solutions are quite accurate even in the critical regime. However, this cannot be regarded as a measure of accuracy of the solutions itself, since only the leading behavior in the limit $Q\rightarrow 0$ is taken into account to determine the critical exponents.

A commonly used method to measure the accuracy of a spectral algorithm is to compare a reference solution with high resolution with solutions of lower spectral resolution, obtained by the same procedure as the reference solution. In particular, we determine all the solutions on a fine grid, using spectral interpolation techniques, and then calculate the differences of a reference solution to the solutions of lower resolution at each of these grid points. We refer to the largest magnitude of these differences as the residue $\mathcal{R}_{N}$, where $N$ indicates the resolution of the less resolved solution. By displaying $\mathcal{R}_{N}$ as a function of $N$ we observe that the residue $\mathcal{R}_{N}$ is rapidly decaying for increasing $N$, see figure 10 for the residue of localized black hole solutions close to the critical transition. The residue $\mathcal{R}_{N}$ saturates at values of order $10^{-9}$ (for $D=5$) or $10^{-10}$ (for $D=6$) which is due to numerical limitations caused by finite machine precision and rounding errors.121212The accuracy of the non-uniform black string solutions whose data we have used in the previous sections is discussed in [12].,131313The numerical calculations were performed by using long double precision of the programming language C.

Furthermore, figure 10 shows the difference $\Delta_{\text{Smarr}}$ between right and left hand side of Smarr’s formula (26) for the corresponding solutions. Again, as the resolution increases, $\Delta_{\text{Smarr}}$ rapidly decreases until a saturation value is reached. Apparently, in $D=5$ this value is about two orders of magnitude smaller than in $D=6$. This is not surprising since in $D=6$ we have to perform two numerical derivatives at the asymptotic boundary to get the value of the coefficient $c_{t}$, which enters into the right hand side of Smarr’s relation (26), while in $D=5$ only the first derivative is of interest, cf. section 2.3.

Finally, we note that the non-trivial components of the DeTurck vector field $\xi$ as given by (13) are always smaller than $10^{-10}$ in magnitude on all grid points even for the numerical solutions closest to the critical transition. This number is negligible compared to the values of the metric functions.

4 Conclusions

We constructed localized black hole solutions in five- and six-dimensional asymptotically flat Kaluza-Klein gravity with one compact periodic spatial dimension using pseudo-spectral methods. Due to the high-precision numerics and a clever choice of coordinates and domain decomposition we are able to investigate those localized black hole solutions which are about to merge into non-uniform black string solutions.

Near the transit solution the thermodynamics associated with localized black hole solutions display a spiraling behavior that adapts to the thermodynamics of highly non-uniform black strings [12].141414For four-dimensional black hole systems that show a spiraling behavior see for instance references [36, 37]. We are able to resolve four turning points. Moreover, we fit the physical observables of localized black holes close to the transit solution to scaling ansätze with only discrete scaling symmetry. The critical exponents agree remarkably well with the ones of the conjectured double-cone metric [22], hence giving compelling evidence in favor of it. The same exponents characterize also the critical behavior of the non-uniform black string solutions. This suggests an infinite inspiral behavior of the two types of solutions. Moreover, due to this spiraling behavior and the associated phase shifts we identify infinitely many regions with positive specific heat.

Critical behavior appears in other higher dimensional black hole configurations as well. Reference [38] showed the existence of critical exponents in the context of hairy black holes in $AdS_{5}\times S^{5}$. Several other works showed the beginning of a spiral curve in the phase diagram of hairy black holes in global $AdS_{5}$ [39] or lumpy black holes in asymptotically flat spacetime [40, 41]. Note that the superradiant instability of Reissner-Nordström black holes gives rise to the hairy black holes mentioned above, while the ultraspinning instability of Myers-Perry black holes leads to the lumpy black hole branch. Recall that in the black hole / black string system discussed here the non-uniform black string branch emanates from the Gregory-Laflamme instability of the uniform black strings. In all of these situations the zero-mode of the corresponding instability leads to the formation of a spiral curve in the phase diagram, which seems to be a generic feature.

Gregory-Laflamme instabilities towards non-uniform black strings and the competition between localized black holes and non-uniform black strings are generic features of Kaluza-Klein gravities irrespective of the boundary conditions such as asymptotically flat or AdS spacetimes. Hence we expect that our results and in particular our high-precision numerical methods may be useful also for other systems such as localized black holes in Anti-deSitter spacetime [42, 13]. In addition, the application of these methods to higher-dimensional versions of the system at hand should be straight-forward.151515However, when going to higher dimensions $D>6$ at least one more technical issue will arise: the accurate extraction of the mass and the relative tension, see the discussion in appendix section A.2. So far, solutions for $D>6$ are rare despite the recent result for $D=10$ [13]. This dimension is of particular interest for two reasons. On the one hand, using solution generating techniques, $D=10$ dimensional asymptotically Kaluza-Klein solutions correspond to type IIa supergravity solutions which are dual to certain thermal states of super-Yang-Mills theory on a circle, by the virtue of the AdS/CFT correspondence. On the other hand, according to the analysis of the double-cone metric, the critical behavior near the transit solution will change for $D\geq 10$. Hence it raises the question, whether there are oscillations of physical quantities in this regime and, consequently, whether the thermodynamics will show a spiraling behavior.

Finally, we note that there are several other systems where the horizon topology changes and the local geometry is expected to be controlled by the double-cone metric. These are for example the transition from pinched rotating black holes to black rings in $D\geq 6$ and the transition from circular pinched rotating black holes to black saturns also in $D\geq 6$, see [43] for more information. According to the results of the work at hand we expect a critical behavior near the transition in all of these examples. In particular, the perturbative analysis of the double-cone geometry predicts the values of the critical exponents depending only on the total number of spacetime dimensions.

Acknowledgments

We dedicate this paper to Marcus Ansorg who pushed forward the application of pseudo-spectral methods to different fields of theoretical physics including relativistic astrophysics, numerical relativity and quantum field theory. Furthermore, Marcus Ansorg utilized these methods in higher dimensional gravity which led to the work at hand.

We thank Burkhard Kleihaus and Jutta Kunz for drawing our attention to this problem. We are grateful to Barak Kol and Julian Leiber for valuable discussions. Furthermore, we thank Toby Wiseman for providing us the data of [9] for comparison. MK and SM acknowledge support by the Deutsche Forschungsgemeinschaft (DFG) graduate school GRK 1523/2.

Appendix A Details on the numerical implementation

The overall numerical scheme is described in section 2.2. In this appendix we present more details. This involves a discussion about the reference metric functions in A.1, and a description of the domain structure and the corresponding coordinate transformations in A.2. Finally, in section A.3, we discuss the parameter values that enter the numerical algorithm.

A.1 Choice of reference metric functions

We have to specify a reference metric in order to implement the DeTurck method. As a first requirement the reference metric has to satisfy the desired asymptotic behavior and boundary conditions, see section 2.1. For $\varrho>\varrho_{1}$ we use the flat metric (2) as a reference, since it already satisfies the conditions on all boundaries except the horizon $\mathcal{H}$, see figure 1. In polar coordinates (7) the flat metric takes the form

[TABLE]

A good starting point for the reference metric near the horizon $\varrho=\varrho_{0}$ is the Schwarzschild metric, or to keep things simple

[TABLE]

where we have used $\mathrm{d}\Omega^{2}_{D-2}=\mathrm{d}\varphi^{2}+\sin^{2}\varphi\,\mathrm{d}\Omega^{2}_{D-3}$. The metric (31) approximates a Schwarzschild black hole with surface gravity $\kappa$ in the vicinity of the horizon.

We have to match these different metrics at $\varrho=\varrho_{1}$. For this purpose we write the reference metric as

[TABLE]

with

[TABLE]

For $\varrho\geq\varrho_{1}$ we obviously recover the flat metric (30). For $\varrho<\varrho_{1}$ the ansatz automatically satisfies the boundary conditions on the exposed axis $\mathcal{A}$ and the lower mirror boundary $\mathcal{M}_{0}$, since the functions $H_{\mathrm{hor}}$ and $G_{\mathrm{hor}}$ do not depend on $\varphi$. The expansion of these two functions close to the horizon has to take the form

[TABLE]

in order to satisfy the conditions on the horizon $\mathcal{H}$, see equation (10), and to approximate the line element (31). Moreover, the functions $H_{\mathrm{hor}}$ and $G_{\mathrm{hor}}$ have to match the flat metric at $\varrho=\varrho_{1}$.

Recall that we construct the solution numerically by means of a pseudo-spectral multi-patch scheme, i.e. using essentially two-dimensional splines of Chebyshev polynomials to cover the entire solution domain. Here, we choose Poincare-Steklov conditions ($\mathcal{C}^{1}$ continuity) at the inner domain boundaries for coupling the patches together. Based on this we now can gauge the full solution by either choosing some $\mathcal{C}^{\infty}$ matching of the reference metric at $\varrho=\varrho_{1}$, or some simplified $\mathcal{C}^{k}$ ansatz. In case of $\mathcal{C}^{\infty}$ matching we use161616We incorporated this ansatz in $D=5$ and $D=6$ with the parameter values specified in section A.3. In any case, it is crucial that the function $G_{\mathrm{hor}}$ stays positive.

[TABLE]

where the auxiliary function $E(\varrho)$ is given by

[TABLE]

Obviously, the function $E(\varrho)$ is one at $\varrho=\varrho_{0}$ and it decays exponentially fast to zero for $\varrho\rightarrow\varrho_{1}$. This $\mathcal{C}^{\infty}$ matching is similar to the approach used in [9].

In case of $\mathcal{C}^{k}$ matching with $k=2$, we use the following simplified ansatz

[TABLE]

The coefficients $h_{1}$, $h_{2}$, $h_{3}$, $g_{1}$, $g_{2}$ and $g_{3}$ can be calculated straightforwardly by matching $H_{\mathrm{hor}}$ with $H$ and $G_{\mathrm{hor}}$ with $G$ up to the second derivative at $\varrho=\varrho_{1}$.

The $\mathcal{C}^{k}$ ansatz has some additional advantages. First, we do not need to take care of the essential singularity for $\varrho\rightarrow\varrho_{1}$ in contrast to the $\mathcal{C}^{\infty}$ ansatz. Second, the auxiliary function $E(\varrho)$ is $\mathcal{C}^{\infty}$ but not analytic on the regarding patch where we use the spectral approximation. This slows the convergence rate of our spectral approximation, which is circumvented by choosing polynomials for $H_{\mathrm{hor}}$ and $G_{\mathrm{hor}}$ which are perfectly analytic within the patches. In practice we saw no difference between the two approaches (apart from numerical fluctuations) when extracting observable quantities with the two different methods with a sufficient number of collocation points, but we saw a faster convergence of the spectral coefficients for the $\mathcal{C}^{k}$ ansatz.171717Note that this matching procedure favors the inner boundary at the coordinate line $\varrho=\varrho_{1}$, see figure 2, which is now justified with hindsight.

The $\mathcal{C}^{k}$ procedure rises the question about the smoothness of the target metric. First, we should keep in mind that we numerically construct a spectral spline approximation and can only monitor its convergence to the desired $\mathcal{C}^{\infty}$ solution (if existing), i.e. also the reference does not necessarily need to be of such a high regularity class. Furthermore, recall that the detailed form of the reference metric is a gauge choice (apart from its boundary behavior). We only need to ensure that it gives rise to some reasonable cover of the underlying manifold and is of sufficient regularity for extracting the regarding observables.

A.2 Decomposition of the domain of integration

We introduce two different charts, the asymptotic chart (3), which is suitable to describe the spacetime in the asymptotic region $x>L/2$, and the near horizon chart (8), which is appropriate for the spacetime in the near horizon region $x<L/2$. The basic domain setup is shown in figure 2.

The asymptotic region

To cover the entire domain five (see figure 2) up to $x\to\infty$, we introduce an appropriate compactification of infinity using the new coordinate $\xi\in[-1,1]$

[TABLE]

where $x=L/2$ corresponds to $\xi=-1$ and $x\to\infty$ corresponds to $\xi=1$. Near infinity we have to take care of the specific behavior of the metric functions. For $\xi\to 1$, there exist exponentially decaying $y$-dependent modes and polynomially decaying $y$-independent modes. The latter ones may even contain logarithms. The bottom line is that the fall-off of the spectral coefficients with respect to $\xi$ will be rather slow. In contrast, when $\xi$ is close to one, the fall-off of the spectral coefficients with respect to $y$ will be very rapid. Therefore, we use the following trick already employed in [12]: The domain five will be further divided into several (e.g. three) linearly connected subdomains, see figure 11. By choosing narrow domains in the vicinity of $\xi=1$ this takes into account the non-analytic behavior of the functions. Moreover, this allows us to save a lot of grid points by adapting the resolution in $y$-direction in each of these domains. Alltogether, with this domain setup in the asymptotic region we achieve more accurate results, while simultaneously the time and memory consumption of the algorithm is reduced.

For the calculation of mass and relative tension we have to read off the asymptotic values $c_{t}$ and $c_{y}$ by differentiating the metric functions $T_{\mathrm{a}}$ and $B_{\mathrm{a}}$ (cf. (3) and (16)) with respect to the compactified coordinate $\xi$ and then taking their values at $\xi=1$. According to (16) we have to differentiate $(D-4)$ times. Therefore, the accuracy of $M$ and $n$ will decrease if the number of spacetime dimensions $D$ increases, since each numerical differentiation is accompanied by small errors. Nevertheless, one could naively think that it is not too difficult to avoid those errors caused by differentiation. For example the function $T_{\mathrm{a}}$ could be expressed as $T_{\mathrm{a}}=1-\tilde{T_{\mathrm{a}}}/x^{D-4}$. The numerical scheme now could solve for the new function $\tilde{T_{\mathrm{a}}}$ from which the value of $c_{t}$ can be read off directly (and without differentiation). Unfortunately, for the problem at hand things are not that simple. We refer to [12] where the non-uniform black string case was treated with a sophisticated ansatz for the metric functions near infinity including a detailed analysis of their asymptotic behavior, which led to highly accurate values for $M$ and $n$. As such an ansatz goes along with some subtle technical difficulties, we decided not doing this effort in this work. However, for the results presented here, which only concern the cases $D=5$ and $D=6$, the accuracy of $M$ and $n$ is reasonably good for all constructed solutions.

**The near horizon region **

As evident from figure 2 the geometries of the domains one and two favor the polar coordinates (7). However, the coordinates (7) are not appropriate to cover the domains three and four. Therefore, we modify the polar coordinates in those two domains such that

[TABLE]

where $v\in[\varrho_{1},L/2]$ is a new radial coordinate and the function $r(v,\varphi)$ takes the form

[TABLE]

The coordinate value $v=\varrho_{1}$ corresponds to the contour $\varrho=\varrho_{1}$, while in domain 3 the coordinate value $v=L/2$ corresponds to $y=L/2$ and in domain 4 the coordinate value $v=L/2$ corresponds to $x=L/2$.

Let us now briefly discuss the functions’ behavior in the near horizon region, especially near the horizon, $\varrho=\varrho_{0}$, and the exposed axis, $\varphi=0$. If we consider solutions of rather tiny localized black holes, or in other words solutions with large $\kappa$ and hence low mass, the functions exhibit steep gradients near the horizon. For a better resolution in this region, we further divide the domains one and two along a contour $\varrho=\varrho_{i}$ with $\varrho_{0}<\varrho_{i}<\varrho_{1}$ into four subdomains. Consequently, we are able to increase the resolution especially in the vicinity of the horizon. Furthermore, by moving along the localized black hole branch we observed that near the exposed axis the functions develop higher and higher gradients.181818This observations concerns the dimensions $D=5$ and $D=6$. We did not consider higher dimensions, but we believe that this property will be qualitatively the same. Therefore, we use the same trick as before by choosing a $\varphi_{i}$ with $0<\varphi_{i}<\pi/4$ and by splitting all domains along this contour. Altogether, instead of the initial four subdomains in the near horizon region (see figure 2) we now have nine subdomains (see figure 12).

There are two further adaptions necessary which are important for the construction of localized black hole solutions close to the critical solution. First, while approaching this solution, all of the functions develop steep gradients near $\varphi=0$ and also the values of the functions $B_{\mathrm{h}}$ and $S_{\mathrm{h}}$ become exceedingly high. To ensure bounded values, we redefine the functions as

[TABLE]

Within our numerics we solve for the new functions $\tilde{B}_{\mathrm{h}}$ and $\tilde{S}_{\mathrm{h}}$. Second, we employ an analytic mesh refinement [44, 45] near $\varphi=0$ to flatten the occurring gradients considerably, which may lead to a much more rapid fall-off of the spectral coefficients. In our setup, we apply this trick by introducing a new azimuthal coordinate $\bar{\varphi}\in[0,\varphi_{i}]$ by

[TABLE]

As the definition suggests, this transformation is only used for $\varphi\leq\varphi_{i}$. The parameter $\lambda$ controls, how strong the flattening is, i.e. the higher $\lambda$ the stronger is the flattening. Usually, there is an optimal $\lambda$ of $\mathcal{O}(1)$. If this is the case, the number of coefficients to be taken to reach a certain accuracy is minimized.

A.3 Parameter values

We defined a set of parameters that enter the numerical scheme. Some of them have an explicit physical meaning ($L$ and $\kappa$), while the other ones are useful to fix the gauge ($\varrho_{0}$ and $\varrho_{1}$) and to control the numerical grid ($\varrho_{0}$, $\varrho_{1}$, $\varrho_{i}$, $\varphi_{i}$, $\xi_{1}$, $\xi_{2}$ and $\lambda$). In this section we provide the parameter values which were convenient for our implementation.

First, throughout the numerical calculation we chose the asymptotic circle size to be $L=8$, which sets a scale for all other parameters.191919Instead, we may scale each parameter by appropriate powers of $L$ in order to obtain only dimensionless quantities. In this way, $L$ completely drops out of the numerical algorithm. In table 3 we list the corresponding values of the relevant parameters. We note that these values correspond to the reference metric ansatz (A.1), which we incorporated to approach the critical regime.

However, these values are not necessarily appropriate to construct a first solution, since the corresponding initial guess (in our case the reference metric) for the Newton method may not be close enough to the actual solution. In fact, we are rather flexible to construct a good initial guess with trial and error by fixing $L=8$ and by changing $\varrho_{0}$ and $\varrho_{1}$ accordingly. Moreover, we used $\kappa\approx 2$ as a starting point. Note that the other parameters do not change the reference metric and hence the initial guess. Once we find a first solution, we slightly modify the surface gravity $\kappa$ to obtain another physically inequivalent solution. In this step we use the former solution as the initial guess. This procedure works fine until we reach a turning point, i.e. an extreme point in $\kappa$. To overcome this point we utilize the trick presented in [27] section VII.B.

Finally, we note that far from the critical regime we used small values for the parameter $\lambda$, i.e. $\lambda<1$, see equation (42). When approaching the transition, we increased this parameter accordingly to flatten the steep gradients near the exposed axis. For our solutions which are closest to the transition we incorporated a value of $\lambda\approx 10$.

Appendix B Short review of the double-cone metric

In this section we review essential features of the double-cone metric. Kol conjectured [22] that the double-cone is a local model of the critical transit solution at the point where the poles of the localized black hole merge or the horizon of the non-uniform black string pinches off, respectively. The Ricci-flat metric of a cone over $\mathbb{S}^{2}\times\mathbb{S}^{D-3}$ reads

[TABLE]

Here, $r$ is the distance from the singular tip of the cone where both the $2$-sphere and the $(D-3)$-sphere have zero size. The question arises how this metric is related to the black hole/black string context. Obviously, the $\mathbb{S}^{D-3}$ represents the inherent spherical symmetry of the setup, while the origin of the $\mathbb{S}^{2}$ is more subtle. It is the fibration of the Euclidean time circles on a path that connects equivalent points on the horizon and crosses the periodic boundary. For example, in the localized black hole case, such a path could go along the exposed axis of spherical symmetry (see for instance figure 7). Clearly, as one starts at the horizon, the Euclidean time circle has zero size, while its size grows when moving away from the horizon and shrinks back to zero size when the endpoint (again on the horizon) is reached. The fibration of these circles produces a topological $\mathbb{S}^{2}$. At the merger or pinch-off point of the transit solution between localized black holes and non-uniform black strings the $\mathbb{S}^{D-3}$ has zero size, since the exposed axis is touched, and the $\mathbb{S}^{2}$ has zero size, since the described path has zero length. This is nicely modelled by the double-cone metric (43). See [46] for a more detailed and pictorial discussion.

Geometrically, the double-cone metric dictates the shape of the horizon of the transit solution near the singular point. The embedding of this horizon into $(D-1)$-dimensional flat space (23) gives

[TABLE]

with an arbitrary constant $Y_{L}$. The tip of the cone is located at $Y=Y_{L}/2$. To compare with the actual localized black hole or non-uniform black string solutions, $Y_{L}$ is chosen to correspond to the length of the compact dimension measured in the embedding coordinate $Y$, see figure 7.

Furthermore, in references [22, 47] perturbations from the double-cone metric of the form

[TABLE]

were analyzed. In linear order of perturbation theory the equations of motion have the following solutions

[TABLE]

Obviously, for $D\geq 10$ the exponents $s_{\pm}$ become purely real, while for $D<10$ the imaginary part of the exponents produces oscillations in $\epsilon(r)$.

In [24] the exponents $s_{\pm}$ were interpreted in the following sense: Suppose $\epsilon$ is some quantity $\delta p:=p-p_{\text{c}}$ which measures the deviation from the double-cone. Here, $p_{\text{c}}$ denotes the critical value of the quantity $p$ of the transit solution between localized black holes and non-uniform black strings. In addition, the coordinate $r$ is rescaled by a characteristic length scale of the perturbed cone, $r_{0}$. For instance, we may choose $r_{0}$ in such a way, that $r_{0}^{-2}$ is a measure of the maximal curvature of the perturbed cone. According to [24] this leads to

[TABLE]

For $D<10$, we obtain after a short algebra

[TABLE]

with $b=-\mathfrak{Re}(s_{+})$, $c=\mathfrak{Im}(s_{+})$ and real constants $a$ and $d$.