Replica wormholes for an evaporating 2D black hole

Quantum extremal islands reproduce the unitary Page curve of an evaporating black hole. This has been derived by including replica wormholes in the gravitational path integral, but for the transient, evaporating black holes most relevant to Hawking’s paradox, these wormholes have not been analyzed in any detail. In this paper we study replica wormholes for black holes formed by gravitational collapse in Jackiw-Teitelboim gravity, and confirm that they lead to the island rule for the entropy. The main technical challenge is that replica wormholes rely on a Euclidean path integral, while the quantum extremal islands of an evaporating black hole exist only in Lorentzian signature. Furthermore, the Euclidean equations for the Schwarzian mode are non-local, so it is unclear how to connect to the local, Lorentzian dynamics of an evaporating black hole. We address these issues with Schwinger-Keldysh techniques and show how the non-local equations reduce to the local ‘boundary particle’ description in special cases.


Introduction
The entropy of Hawking radiation is a diagnostic of information loss. It was long believed that this entropy could only be computed in the ultraviolet theory. Recently, however, the entropy was calculated in the low-energy theory by Almheiri, Engelhardt, Marolf and Maxfield [1] and simultaneously by Penington [2] using an extension of the gravitational entropy formula [1][2][3][4][5][6][7][8]. They discovered that after the Page time, a quantum extremal island appears in the black hole interior. The gravitational entropy formula in this context, called the 'island formula' [9], is where R is the radiation, I is the island, and S QFT is the entropy of quantum fields calculated by the traditional methods of quantum field theory in curved spacetime. In two-dimensional gravity, 'Area' means the value of the dilaton. Almheiri et al. [1] did the calculation in a two-dimensional model of Jackiw-Teitelboim (JT) gravity [10,11] glued to non-dynamical flat spacetime. Penington [2] also gave some general arguments that the quantum extremal surface should exist for higher-dimensional black holes. In further work, the entropy formula has been applied to various other setups and it consistently produces a unitary Page curve [9,[12][13][14][15][16][17][18][19][20][21][22][23][24][25][26][27][28] (see [29] for an introductory review).
The island formula (1.1) was originally proposed based on holographic reasoning, and therefore assuming unitarity. It was then derived without holography from the gravitational path integral in [13,14], using the gravitational replica method developed earlier in [5,30]. This entails calculating Tr (ρ R ) n by a Euclidean path integral, and analytically continuing n → 1. At finite n, there are 'replica wormholes' connecting the different copies of the black hole through their interiors. In the limit n → 1, the mouth of the wormhole degenerates to the island region, I.
In this paper we undertake a detailed analysis of replica wormholes for a black hole in two-dimensional JT gravity that begins near the vacuum state, forms by gravitational collapse, and then evaporates. We focus on n ∼ 1. In this example we can be much more explicit about the global, real-time solutions in Lorentzian signature, as compared to [5,30] and the general arguments in [13,14]. Along the way, we will gain a better understanding of how to apply gravitational path integral techniques to collapsing black holes. This is important because the replica wormhole derivation of the island formula requires a state prepared by a Euclidean path integral -but this seems to exclude ordinary black holes formed by collapse, because they have no time-reflection symmetry and therefore no real Euclidean continuation. We will overcome this by viewing the collapsing black hole as a limit of a Euclidean solution that scales away the time-reversed process.
We will also repeat the derivation of the island formula in this explicit example. That is, we will fully define the global replica equations in this Lorentzian setup, analyze them as n → 1, and use the Schwarzian equations to show that the gravitational path integral in the low energy theory reproduces (1.1). We assume a large-N matter sector. As with all JHEP04(2021)289 derivations of gravitational entropy from the replica method, this is on a similar conceptual footing to the Gibbons-Hawking derivation of the area law -it computes the entropy without telling us what microstates are responsible, or whether such microstates really exist.
Finally, this example illustrates some conceptual points that are hidden in the derivation of the island rule based on the action [13,14], such as the role of conformal welding, the relation between the local analysis [5] and the Schwarzian theory, and other subtleties that arise in Lorentzian signature. None of these subtleties lead to any surprises in the end, but there could be other contexts, like cosmology [16,17,[31][32][33][34] or applications to other observables [35][36][37][38][39][40], where the answers are less clear a priori and these subtleties come into play.
Our setup is similar to Almheiri, Engelhardt, Marolf, and Maxfield (AEMM) [1]: an evaporating black hole in AdS 2 glued to non-dynamical flat spacetime. (One important difference is that we create the black hole with an operator insertion, instead of a joining quench.) The Lorentzian theory is described by the Schwarzian action coupled to an external system, which adds a source term in the Schwarzian equation of motion [41][42][43].
Our main new results are as follows: • We construct an evaporating black hole from the Euclidean path integral. This is nontrivial because the Euclidean equations in the Schwarzian theory glued to flat space involve the non-local (and generally intractable) conformal welding problem. We start with a shockwave created by operators inserted in the Euclidean path integral. Then we take a scaling limit where the nonlinear conformal welding problem becomes exactly solvable. This turns out to reproduce the local Lorentzian equations for the Schwarzian 'boundary particle' derived in [41,42] and studied in AEMM [1]. (See section 3).
• We derive the replica equations for the Schwarzian theory defined on a Schwinger-Keldysh contour. This is necessary to study replica wormholes for evaporating black holes because the nontrivial island we seek does not exist in Euclidean signature.
(See section 4.) • We derive the extremality conditions for the island from the replica equations of motion in the Schwarzian theory as n → 1. (See section 5.1.) The derivation applies to any state in this theory created by a Euclidean path integral with local operator insertions, generalizing the eternal black hole analysis in [13]. In the above-mentioned scaling limit, this includes black holes that form by collapse, then evaporate.
• Using the Ward identity for CFT coupled to gravity, we show that the extremality conditions can be integrated to yield the entropy formula (1.1). This argument is equivalent to the derivation of the entropy in [13,14] from the defect action, but easier in practice. (See section 5.3.) • We also analyze the two-interval replica wormhole for an eternal black hole at late times, elaborating on a calculation in [13]. We show that at late times the wormhole factorizes into two copies of the one-interval solution. This result, which is logically JHEP04(2021)289 Figure 1. A black hole in AdS 2 that is created by a shockwave, then evaporates. The AdS 2 region is initially in vacuum. We will calculate the von Neumann entropy of region R.
independent from the rest of the paper, confirms a physical argument made in [13] that wormholes factorize in the OPE limit of the twist operators. This assumption was necessary in order to apply replica wormholes to the information paradox of an eternal black hole. (See section 6.) Our results give an explicit realization of the Schwinger-Keldysh approach to the gravitational entropy formula developed in [30], incorporating quantum effects. In the rest of this introduction, we describe our setup in more detail, highlight some technical challenges, and summarize how things unfold.

Summary
We are primarily interested in the Page curve for the black hole in figure 1. The theory is JT gravity in AdS 2 coupled to a large-N CFT, glued to a flat spacetime as in [1,42]. The CFT lives everywhere while gravity lives only in AdS 2 , so the flat spacetime is nondynamical. A delta-function shockwave is sent from I − to form a black hole. The black hole evaporates and Hawking radiation escapes toward I + . We will calculate the entropy of region R. The final answer [1] is given by the gravitational entropy formula (1.1), with an island that appears after the Page time.
For t > 0, this solution is nearly identical to the evaporating black hole studied in [1] (AEMM). In the AEMM setup, the AdS 2 and flat regions are initially disconnected, and the shockwave is created when these two regions are suddenly joined at t = 0. We will not take this route because we found it difficult to study this setup with Euclidean path integrals. Fortunately our t > 0 solution for the Schwarzian mode in Lorentzian signature is identical, after taking a limit described below, so for the island analysis we can borrow the results of AEMM (and its extension to finite β in [19,22,23]). Instead of figure 1, we will start with the more general setup in figure 2. At t = 0 the gravity region is a black hole at temperature β. The solution is time-symmetric, so there are ingoing and outgoing shockwaves. The shockwaves are not delta functions because we want the observables, and in particular the entropy, to be analytic functions of position (up to the usual singularities associated to coincident points in Euclidean and lightcones in Lorentzian which are handled by the i prescription). So the shockwaves are smeared over a width δ. This state is defined by inserting the CFT operators ψ(y 1 )ψ(y 2 ) (1.2)

JHEP04(2021)289
into the Hartle-Hawking path integral, with Here y is a complex coordinate defined below. L is the distance from the AdS boundary, and the shift by iδ is an offset in imaginary time that smears out the shockwave. The shockwave operator ψ is a scalar primary with conformal weights h ψ =h ψ , related to the energy of the shockwave by The CFT state created this way has been studied extensively in [44][45][46][47][48][49][50], mostly for applications to the AdS 3 /CFT 2 correspondence, with the shockwave operators inserted in the boundary CFT. Here we will use similar techniques but there is no holography; the operators are inserted in the matter CFT, which is directly coupled to gravity in the AdS 2 region. The operators are inserted in the background of the Euclidean eternal black hole. But even at t = 0, the geometry backreacts, so the gravitational solution is not exactly an eternal black hole. The reason is simply that the CFT stress tensor produced by the operator insertions is everywhere non-vanishing. This backreacts on the geometry, which JHEP04(2021)289 feeds back into the matter stress tensor and thus leads to a complicated set of non-local equations coupling the geometry to the matter stress tensor. The equations involve an implicit solution to the conformal welding problem, and cannot be solved analytically at finite δ. We will not solve these equations, but we will write them down, and solve them in the limit δ → 0.
Even before getting to the replica problem, this addresses a puzzle raised by comparing the two papers [1] and [13]. In [1], the Lorentzian geometry associated to a shockwave was found by solving a relatively simple local equation for the Schwarzian mode [42] at the boundary of AdS 2 . By contrast the Euclidean equations of [13] are non-local due to conformal welding whenever there are nontrivial operators inserted, including the operators that produce the shockwave. How does conformal welding in Euclidean signature lead to nice, local equations in Lorentzian signature? The answer we will find is that as δ → 0, the conformal welding equations have an exact nonlinear solution in Lorentzian signature that agrees with [1]. This is somewhat surprising, since it is impossible to find exact nonlinear solutions to conformal welding in Euclidean signature, but it had to be the case for consistency between [1] and [13]. We will postpone the more technical discussion until later, but the basic observation is that nonlinear conformal welding is tractable when the Schwarzian mode is purely positive or purely negative frequency, a situation that is possible only in Lorentzian signature. In the δ → 0 limit the outgoing shockwave decouples from the ingoing shockwave -this simplifies the conformal welding problem, and then the Euclidean equations reduce exactly to the 'boundary particle' equations for the Schwarzian theory studied in [1,41,42]. 1 The boundary particle equation does not apply in the more general case with simultaneous ingoing and outgoing matter.
At finite δ, the black holes in figures 1-2 are described as solutions of the Schwarzian equation on a Schwinger-Keldysh contour. The equations can be solved as δ → 0. The one-sided shockwave in figure 1 is obtained by sending the operator insertions to I − . (Alternatively we could smear the insertion against a wavepacket, but we will not do this in detail.) Another puzzle in figure 1 is the role of the left endpoint of the island. It appears to sit on the left boundary of AdS 2 . However this is potentially problematic because in calculating the entropy of the radiation region R, we should not include the boundary point on the left side of the Penrose diagram. We could therefore question whether the replica manifolds associated to figure 1 obey the correct boundary conditions on the left boundary. This puzzle is eliminated by going to finite temperature in figure 2. The left endpoint of the island is now a second quantum extremal surface that sits near the bifurcation point of the original black hole. Now it is manifest that the correct boundary conditions are obeyed. So far we have described the background solution, or the n = 1 replica manifold. We now turn to finite n, which is relatively straightforward now that we have framed the problem in a convenient way. Our strategy (see section 4) will be to write the gravitational JHEP04(2021)289 equation of motion (i.e., the Schwarzian equation) in Lorentzian signature and solve it in the limit δ → 0, n ∼ 1. This gives the extremality condition for the quantum extremal surface. This derivation (see section 5) is very general for JT gravity, not limited to the shockwave state. It combines elements of [13] and [14]. Conformal welding enters this calculation in the intermediate steps, but it ultimately drops out of the extremality condition. We will also discuss how this derivation relates to the local analysis of the dilaton equations of motion around the defect that was used by Dong and Lewkowycz [51] to derive the quantum extremality condition.
Finally, we must calculate the entropy from the gravitational action plus the 1-loop effective action of the matter fields. There is a shortcut: we demonstrate using the Ward identity that once the extremality conditions have been derived, the entropy automatically agrees with (1.1) (see section 5.3). This is similar in spirit to Cardy and Calabrese's derivation of entanglement entropy in CFT from the conformal Ward identity [52], but the argument is modified in the presence of gravity.

Evaporating black holes in JT gravity plus a CFT
In this section we will review the formulation of the information paradox in JT gravity along the lines of [1], in Lorentzian signature. Readers familiar with [1] can skip to section 3.
There are two differences in our setup compared to [1]. The first is that we retain a finite temperature parameter β for the initial black hole; in [1] the temperature is taken to zero, but the case of finite β has been studied in [19,22]. The second difference is that our shockwave is produced by local operator insertions, rather than a joining quench. The Lorentzian solution is identical for t > 0, so this will not make any difference until the next section.

Jackiw-Teitelboim gravity theory plus a CFT
We begin with the Lorentzian theory in AdS 2 coupled to a CFT on just one side of AdS. We will generalize to the 2-sided gluing below. The action of Jackiw-Teitelboim gravity in AdS 2 coupled to a CFT is where (in Lorentzian signature) where we set AdS = 1. The first line is topological. We take the matter action to be independent of the dilaton. We couple this system to the same CFT living in a portion of rigid Minkowski spacetime and impose transparent boundary conditions for the CFT at the interface [1,13]. Variation with respect φ yields the equation R + 2 = 0 which requires JHEP04(2021)289 the metric to be locally AdS 2 . We take Poincare coordinates x ± for the AdS 2 geometry in the interior region and y ± = t ± σ for the flat space in the exterior region, The two spacetimes are glued together near the AdS 2 boundary according to a map x(t) that will be determined dynamically. Given x(t), the gluing identifies the curve y + = y − = t in the center of Minkowski spacetime to the curve in AdS. In the limit → 0 the gluing is simply along x ± = x(t). At the interface, the boundary conditions for the metric and the dilaton are The variation with respect to the metric yields the dilaton equations is set by the anomaly and can be absorbed into the definition of φ 0 , so we will adopt this convention and drop it from the first equation. Then the general solution of these equations is [53] φ(x + , x − ) = − 2πφ r β up to SL(2) transformations. As argued in [41], it is convenient to express the dynamics of JT gravity in terms of the shape of the boundary curve, or gluing map, x(t). The dynamics of x(t) are governed by the Schwarzian action. After setting R = −2, the JT action is a pure boundary term that evaluates to The ADM energy of the gravity region is given by Conservation of energy relates dM/dt to the net flux of energy across the interface, This is the equation of motion for the gluing map x(t).

JHEP04(2021)289
We will be particularly interested in situations where T x − x − = 0 at the interface. In this case we can use the equation of motion to re-express the dilaton (2.8) in terms of the gluing map as Details about the derivation of this expression can be found in appendix C.2. If AdS 2 is glued to Minkowski spacetime at both of the AdS 2 boundaries, then this expression for the dilaton holds in a region on the right side of AdS 2 that is spacelike separated from the left interface.

Eternal black hole coupled to a bath
We now review the eternal black hole coupled to an external bath system which is analyzed in [1]. We couple an eternal black hole at inverse temperature β to Minkowski spacetime at the same temperature. In the Lorentzian picture, the geometry is an eternal black hole in AdS 2 connected to two non-gravitating half-spaces, one on each side as drawn in figure 3.
The gravitational system is in equilibrium with the bath. Therefore the net flux is zero, and the Schwarzian equation of motion is The solution corresponding to inverse temperature β is We can also act on this with an SL(2) transformation, which would correspond to a different coordinate system in AdS 2 . A priori, the gluing map x(t) is defined only at the interface, so it is a real function of a real variable, and the y ± coordinates are defined only in the Minkowski region. But we can extend the y ± coordinate system into the gravity region by the coordinate change This can be done for any gluing map, but it is only useful in certain cases (in particular we will find this does not suffice to describe the replica wormholes or finite-δ shockwaves). In the eternal black hole the reason this is useful is that the gluing map can be analytically continued to a holomorphic function on the Euclidean disk (as we will see in more detail below), and this provides a natural Hartle-Hawking state for the CFT in which T x ± x ± is given entirely by the conformal anomaly. In this state, {x(y ± ), y ± } = πc 12β 2 .  Putting the solution (2.14) into (2.3) and (2.12), we obtain the following physical metric for the inside region and the dilaton profile . (2.17) An SL(2) transformation acting on x ± relates different choices for how to embed the (y + , y − ) coordinates into the Poincare patch. Our choice (2.15) is depicted in figure 3. 2 The coordinates (y + , y − ) cover the right wedge of the eternal black hole. The left wedge corresponds to the region x ± < 0 which can be obtained by the analytic continuation y ± → y ± + iβ/2. The Euclidean continuation is where w is a complex coordinate. The gravity region is the interior of the unit disk, with the hyperbolic metric The hyperbolic disk is glued at the boundary to a flat spacetime with the metric for |w| > 1. See figure 3. The Euclidean continuation of the solution (2.14) is given by

JHEP04(2021)289
which is the map between the Euclidean cylinder and the plane coordinate. The black hole mass and the Bekenstein-Hawking entropy are given by respectively. Here S 0 is the extremal entropy given by S 0 = φ 0 4G N , and φ h represents the dilaton value at the horizon.

Evaporating black hole
We will now add a shockwave to the eternal black hole, leading to an evaporating solution. At times |t| < L, we have the eternal black hole at temperature β in equilibrium with the bath. The evaporating black hole is created by injecting a shockwave from the flat spacetime region at time t = L along the AdS boundary. To make the solution timesymmetric, we also have a shockwave exiting the AdS region at t = −L. For now, these shockwaves will be treated as delta function sources, but below we will describe how they can be obtained by the singular limit of a smooth state with operator insertions in Euclidean signature. Since the solution is time-symmetric, we will focus on t > 0.
Injection of the ingoing shockwave increases the temperature of the black hole and also changes the gluing function x(t) as compared to the eternal black hole. While the ingoing stress tensor remains thermal with the original temperature β away from the shockwave, the outgoing stress tensor is changed according to the new gluing function: where E ψ is the energy of the shockwave. The Schwarzian equation of motion for t > 0 is Integrating once gives We have set the integration constant to agree with the eternal black hole Schwarzian at early times, {e 2πt/β , t} = − 2π 2 β 2 , and defined a parameter which controls the evaporation rate. This represents a black hole that increases in mass when the shockwave enters, then decays back toward the original mass,

JHEP04(2021)289
From the above expression, we can define the temperature of the evaporating black hole as where we introduced a dimensionless parameter which controls the increase in temperature of the black hole due to the shockwave with energy E ψ . In this paper, we will consider the weak gravity limit κ 1 and c 1 with κE ψ /c fixed , (2.31) such that the change in temperature is finite.
The solution x(t) to the equation of motion (2.25) is explicitly given by which indicates that the new horizon sits outside the old horizon of the original black hole before the shockwave as depicted in figure 2. For the arguments in the following subsections, we need to know the late-time behavior of the solution x(t). At late times t ∼ O(1/κ), the solution can be approximated as where we used the asymptotic behavior of the modified Bessel functions by fixing κ(t − L) and κE ψ in the κ → 0 limit. Here we introduced a function η(u) (2.36)

JHEP04(2021)289
The low-temperature limit β → ∞ corresponds to η(u) → u and ν → 0, and in this limit we recover the solution studied in AEMM [1]. Taking derivatives with respect to t, we also have the following expressions . From these expressions, we can check the solution satisfies (2.26) at the level of this late-time approximation.

Hawking calculation of the entropy
In the metric −Ω −2 dz + dz − , the vacuum entropy of a CFT on an interval [z 1 , z 2 ] is . (2.38) This is the entropy in the vacuum state with respect to the z coordinate. To measure the entropy of the Hawking radiation of the evaporating black hole, we will apply this to an interval R : [σ B , σ B ] placed in the right flat region at time t. The endpoint σ B is an IR regulator which is taken to be large.
In the vicinity of point B, as we can see from the stress tensor (2.24), the CFT is in the vacuum state with respect to the coordinate That is, left-movers are in a thermal state at inverse temperature β, while right-movers have a thermal contribution at the same temperature plus a contribution from Hawking radiation. Near the point B , both left-and right-movers are thermal at inverse temperature β, so The conformal factors in the metric −dy + dy − are Ω 2 = z (y + )z (y − ). Therefore, applying (2.38), the von Neumann entropy of the radiation in the standard semi-classical approximation is is the entropy of the Hawking radiation, with x(t) the gluing map obtained in (2.32). The second term in (2.41) πc 3β (σ B − σ B ) is the thermal entropy of the initial equilibrium state. S shock is here because we are not in the vacuum state -it is the entanglement of the matter in the ingoing shockwave with the matter in the outgoing shockwave. In a coherent state, this would vanish, S shock = 0 [54]. In a state created by local operator insertions, S shock is non-zero, but assuming a large-c limit (with h ψ /c held fixed but small), it is a constant JHEP04(2021)289 that grows only logarithmically with the energy [45,48]. This is subleading compared to other contributions to the entropy and we will therefore neglect it.
At late times t ∼ O(1/κ), the entropy of the Hawking radiation increases as and finally asymptotes to a finite value 3 As a diagnostic for information loss we will compare this to the entropy of the black hole. Actually, this setup is slightly more complicated than an ordinary evaporating black hole, because the Hawking radiation is entangled with both the black hole and the left Minkowski wedge. To allow for this additional entanglement we should compare to the generalized entropy of the complementary region R C , which extends from spatial infinity in the left Minkowski region to the point B. That is, we consider it a 'paradox' if The generalized entropy of R C is calculated in a quasi-static approximation, where we sum the gravitational entropy at the two event horizons and the entropy of the quantum fields in the union of the two exterior regions (one extending from left infinity to the left event horizon, and the other extending from the right event horizon to B). The important piece of this entropy is the gravitational contributions at the left and right event horizons, yielding at late times t ∼ O(1/κ). This asymptotes to twice the entropy of the original black hole entropy before the shockwave. 4 For large u 0 , we find S QFT (R) > S gen (R C ), so the black hole does not have enough entropy to purify the Hawking radiation. This is the information paradox.

Island calculation of the entropy
The island formula [1,2,9,14,29] states that the true von Neumann entropy of region R is given not by Hawking's calculation, but by an extremum of the generalized entropy, Plus terms from 'const.' in (2.42), including UV divergences from the endpoints. These terms cancel in a more careful comparison to the generalized entropy of R C . 4 We have dropped subleading matter contributions, which are suppressed by κ, and UV and IR divergences from the matter entropy, which cancel against identical divergences in S QFT (R). We have done this calculation by the old-fashioned, pre-Ryu-Takayanagi method, using the quasi-static approximation in order to avoid using the island rule in the statement of the paradox. The more accurate semiclassical calculation of the entropy of R C would actually use the island rule, i.e., the gravitational contributions would be evaluated at the quantum extremal surfaces. The difference between these two calculations is subleading at small κ, so it is sufficient to use the old-fashioned calculation of Sgen(R C ) to state the paradox. The island has two ends A and A determined by the extremal conditions. As we will see, the island is placed outside the original horizon of the right black hole, which is depicted as the dashed null lines.
The term S div is a counterterm subtracting the UV divergences in S QFT coming from the boundary of the island. 5 In two dimensions, 'Area' means the value of the dilaton. At early times I is trivial and this formula leads to the usual Hawking entropy. After the Page time, a non-trivial island appears inside the black hole and the entropy starts to decrease, leading eventually to the same entropy as the original black hole, as required by unitarity. Now let us briefly explain how to compute the entropy of the Hawking radiation using the island formula. We assume the energy of the shockwave is large u 0 1. To find a non-trivial QES, we also assume late enough times t ∼ O(1/κ) when the black hole has lost an order one fraction of its mass (see [1]), but no too lates times u 0 e −κ(t−L)/2 1 so that the mass is still much bigger than the original one. We take a single interval R which lies in the right flat region. In this case, the island has two ends; we call the right end A and the left end A as depicted in figure 4. The area terms associated to A and A can be easily computed from the formula (2.12) by using the asymptotic expansion of the solution x(t) (2.34), (2.37). For the variations of the matter entanglement entropy S QFT (I ∪ R), we assume the factorization up to small corrections from the non-factorizable part of the entropy which we can neglect when evaluated at the QES. We will check that a posteriori once we have solved for the QES. The extremization for each endpoint of the island is then independent of the other. We also assume the points A and B are both placed before the shockwave, where the

JHEP04(2021)289
background solution is locally equivalent to the black hole at inverse temperature β. In the limit σ B → ∞, the entropy S QFT (A , B ) is extremized at σ A = −∞, i.e., the bifurcation point of the original black hole. 6 Therefore the generalized entropy for the region The term cπ 3β σ B is the IR divergent thermal contribution, far from the black hole. Next, we will extremize the right QES, A. The calculation of the matter entropy for [A, B] is similar to the calculation for [B, B ] around (2.41), except that now both endpoints are after the shockwave, so in (2.38) we use the relations for both endpoints, with v(y) = e 2πy/β . The other difference is that because point A is inside AdS, it has the conformal factor Ω = 1 This leads to the generalized entropy . (2.51) The derivative of the generalized entropy with respect A can be expressed as at leading order for large u 0 . Here we defined . Thus the location of the right QES satisfies As observed in AEMM, by sending the point B to the AdS boundary y − B = t, the outgoing coordinate of the quantum extremal surface y + A is associated to the boundary time t as With R in the right flat region, the island should be outside the original eternal black hole, since R cannot be used to reconstruct information on the left. This follows from the explicit calculation, and is also guaranteed by the quantum normal conditions derived in [1,16]. It also follows from the quantum focusing conjecture [12]: if the left end of the island A is placed behind the original horizon of the right black hole, one can send a light ray from A to the left physical boundary where φ ∼ φr/ . The quantum focusing conjecture, taking into account of ∂Sgen(I ∪ R) = 0 at the QES A , implies that Sgen would monotonically decrease toward the left boundary along the light ray. This contradicts with the fact that the dilaton, and so Sgen as well, diverge at the left boundary.

JHEP04(2021)289
. This time delay can be interpreted as the scrambling time in the thermal system with temperature T (t) and gives a nice geometric realization of the Hayden-Preskill protocol [1,2,55].
Plugging the solution for the quantum extremal surface back into the entropy, after the Page time we have S(R) ≈ S BH where the black hole entropy S BH was given in (2.46). So with the island prescription, the entropy of the radiation obeys the unitary Page curve.
We assumed factorization of the entanglement entropy to determine the QES. Now, we must go back and check that this assumption is self-consistent, given the resulting QES. This analysis is done in appendix A. Actually we find that the entropy does not quite factorize at the position of the QES, but the additional terms do not affect the extrema at leading order in κ. That is, any large non-factorized contributions to the entropy drop out when we take the derivatives with respect A or A .

Shockwave in Euclidean signature
Our goal in this section is to obtain the evaporating black hole, joined to an external bath, from the Euclidean path integral. As usual, once the quantum state and equations of motion for geometry and matter are constructed in Euclidean signature, real time observables are given by analytic continuation. We will focus on the state created by a local operator insertion, but it should be straightforward to generalize the methods to arbitrary states created by a Euclidean path integral with operators or sources inserted in the non-gravitating region.

Euclidean setup
At finite temperature we prepare the Euclidean path integral by imposing periodic boundary conditions in Euclidean time for quantum fields and gravity. A semi-classical picture JHEP04(2021)289 of the path integral is shown in figure 6. We require that matter fields and the metric be continuous at the interface between the two regions.
In the flat region, we use coordinates y = σ + iτ,ȳ = σ − iτ , We set the inverse temperature to β = 2π. The shockwave is created by a scalar primary operator ψ inserted in the flat region. In the path integral we have two operator insertions (for the bra and ket), The shockwave operator ψ has conformal weights (h ψ , h ψ ) and scaling dimension 2h ψ . This is similar to the local operator quench [45,56] studied in the literature of 2d CFT. However, in the present setup gravity is dynamical.
In the interior of the disk corresponding to the gravity region, the coordinate is denoted by w. In JT gravity, the metric is hyperbolic, so for the disk topology we have The boundary of the interior region in the Schwarzian limit → 0 is denoted by w| |w|=1 = e iθ . The dynamics is encoded via the Schwarzian action in a non-trivial diffeomorphism between the boundary of the disk θ and the time τ along the boundary of the flat space region. The equation of motion for the gluing function θ(τ ) involves the energy flux across the interface, which in turns depends on the manifold and therefore on θ(τ ). So far two different coordinates w and y corresponding to the interior and the exterior of the disk were introduced. Gluing the two regions amounts to finding a coordinate system which covers both regions. The coordinate system is not unique and could be explicitly written for non-holomorphic extensions of the form (w(y,ȳ),w(y,ȳ)). However, a nonholomorphic coordinate extension would not allow us to compute the CFT stress tensor. Instead, we look for a holomorphic coordinate which covers both regions. In practice, this JHEP04(2021)289 means given e iθ(τ ) , we implicitly define a coordinate z that covers both gravity and flat regions where the metric can be written globally on a plane as ds 2 = Ω −2 (z,z)dzdz. The maps y → z and w → z must be holomorphic on their respective domains. Once the metric is known in this form, the stress tensor in the original coordinate system is calculated by the Weyl anomaly.
Finding the coordinate z given θ(τ ) is known as the conformal welding problem [13,57]. We must find two analytic maps G, F from the interior and exterior of unit disk respectively, to the z plane such that and where v = e y . See figure 7.
Here θ(τ ) is a bijection and in particular θ(τ + 2π) = θ(τ ) + 2π. Due to the Riemann mapping theorem, the solution for bijective maps G, F always exists and each map is unique up to a PSL(2, R) transformation. However, it is generally not possible to write a solution in closed form.
From the stress tensor on the z-plane, the original stress tensor is determined by transforming back to the y-plane. Assuming the functions F, G are known, and using v = e y , we have {F (e y ), y}, where T zz is the stress tensor in the z-plane with the Weyl-rescaled metric ds 2 = dzdz.

JHEP04(2021)289
For the insertion of two operators, the stress tensor on the z-plane is determined by conformal invariance and given by with a similar expression for the anti-holomorphic stress tensor. Here z 1 = F (e y 1 ), Finally, using the flux of energy (3.8) evaluated at the interface, we can write the Schwarzian equation of motion for e iθ(τ ) , Note that (3.9) is not simply a differential equation for θ(τ ), because the welding function F depends non-locally on θ(τ ).
However, we will see that in the limit of a delta-function localized shockwave, the equations become local and (3.9) can be solved. This is the limit δ → 0 with held fixed. In this limit, the stress tensor (3.8) vanishes away from the insertion points. Therefore the Euclidean gluing map becomes trivial, θ(τ ) = τ for real τ ∈ (0, 2π), and the operators are inserted on the z-plane at The stress tensor has delta function contributions supported on the ingoing and outgoing lightcones emanating from the point y = L, Our strategy will be to setup all of our calculations at finite δ, analytically continue to Lorentzian signature, then take δ → 0. Observables in the localized shockwave state are not analytic, but they are analytic at any finite δ, so it is important to do things in this order.

The welding solution for small E ψ
As a warmup, let us first consider the case that E ψ is small and fixed as δ → 0. This was treated in appendix B and C of reference [13] in detail, and here we review that argument. We need to solve the welding problem and equation (3.9) to first order in E ψ . The holomorphicity conditions on F, G in their respective domains requires them to have the series expansions

JHEP04(2021)289
The PSL(2, R)×PSL(2, R) ambiguity in the welding problem is fixed by setting The zeroth order solution in E ψ to the welding problem and (3.9) is given by the trivial solution, with F and G the identity maps, (3.14) Expanding the matching condition (3.6) to first order around the trivial solution θ(τ ) = τ + δθ(τ ), we find Note that as a function of z = e iτ , δG has a series expansion with only positive powers of z, whereas δF contains z m , m ≤ 0. Therefore δG and δF are given by where +, − indicate positive and negative frequency projections defined with respect to the background welding coordinate, z = e iτ . That is, for any meromorphic function K = ∞ n=−∞ a n z n , A convenient way to calculate these projections is by the contour integrals . (3.18) This expression is valid for K(z) + when |z| < 1 and it is valid for K(z) − when |z| > 1.
Using (3.15), the Schwarzian equation of motion to first order becomes where δS = δ{e iθ , τ } and F shock = T shock Decomposing the flux into positive and negative modes, (3.19) decouples into two separate equations where the flux from (3.8) has positive and negative projections In the limit δ → 0, the flux vanishes in Euclidean signature, but continuing to Lorentzian with τ = it we find where E ψ = h ψ /δ. The solution is where we chose decaying boundary conditions δS → 0 as |t| → ∞. These solutions are plotted in figure 8. Recall from section 2 that S + + S − is proportional to the time-dependent black hole mass, M (t). From the figure we see that δS − describes the formation and decay of a black hole from the infalling shockwave, and δS + is the time-reversed process.
There is no information paradox for a single interval when the geometry is near the eternal black hole solution (δθ small) -the sharp paradox of AEMM discussed in section 2 required the entropy of the black hole created by the shockwave to be large. So we will now turn our attention to non-linear geometries corresponding to fixed E ψ which is not assumed to be small.

Nonlinear solution of welding
The discussion of section 3.2 showed that to first order in E ψ , at small δ, the solution for δS = δS + + δS − consists of only positive or only negative modes at any given t. We call this the no-mixing condition:

JHEP04(2021)289
At the linearized level of the time-symmetric shockwave, this was a consequence of taking δ → 0. For finite δ, the positive and negative solutions have overlapping support (see appendix B for a numerical analysis). In Euclidean signature, the reality condition for θ(τ ) makes it impossible to satisfy the no-mixing condition except in the trivial case θ = τ . However, as evidenced by section 3.2, the no-mixing condition has non-trivial solutions in Lorentzian signature.
For the shockwave configuration, the no-mixing condition is the key to solving the nonlinear welding problem at finite E ψ . For example, if we expand θ(τ ) ≈ τ + δ (1) (2) θ are first order and second order corrections in E ψ , and solve the matching condition (3.6) to the second order, we find where equation (3.16) is used in the first line, and primes denote derivatives with respect to the background coordinate z = e iτ . In order to solve for δ (2) F , we project onto negative modes and similarly for δ (2) G. If the support of positive and negative solutions have an empty intersection, the second term in equation (3.28) is zero. In the same way, the no-mixing condition can be generalized to all orders in perturbation theory. For an expansion to higher These are identical to the linearized equations (3.16) but now the variations are nonlinear, Here is another perspective. Let us assume the ansatz [δθ] + = 0 for t > 0, and [δθ] − = 0 for t < 0. We need to solve the matching condition G(e iθ ) = F (e iτ ). The nonlinear solution we have just constructed is simply

JHEP04(2021)289
In Lorentzian signature, τ = it, t ∈ R, and θ(it) is purely imaginary. Therefore the solutions (3.29) are real for real t, and the reality conditions becomē We can now restate the welding solution for t > 0 as where we introduced the gluing function which is real for t ∈ R, and functions F + (e t ) =F (e t ) , F − (e t ) = 1/F (e −t ) , G + (e −iθ(it) ) = G(e −iθ(it) ) , G − (e iθ(it) ) = 1/G(e iθ(it) ) in Lorentzian signature. 7 Here θ(it) has only negative modes when t > 0 and in particular it does not have to be a small deformation of the eternal black hole solution. Given this ansatz for the welding solution, the Schwarzian equation of motion in (3.9) becomes where on the right-hand side we have plugged in the shockwave stress tensor obtained in (3.12). These are identical to the equations of motion (2.25) obtained by Lorentzian methods.
At this point we have resolved one of the puzzles described in the introduction: the Lorentzian methods of [1,41,42] led to relatively simple, local dynamics, whereas the Euclidean approach of [13] requires non-local solutions of the conformal welding problem. We have shown that these two methods are in perfect agreement. The boundary particle equations of motion used in [1,41,42] correspond to a limit where the no-mixing condition is satisfied, so that the welding problem can be solved exactly and used to write simple Lorentzian equations. In other words, we recover the Lorentzian equation of motion from [41,42] when the only outgoing matter is that coming from Hawking radiation, as opposed to the explicit outgoing shockwave. The equations at finite δ are more complicated because the outgoing shock does not decouple and the equations remain non-local. 8 The no-mixing condition is what allowed us to find an exact solution of conformal welding. This condition has a nice interpretation in the x-coordinate: it can be satisfied if and only if the stress tensor is chiral, T x − x − = 0. In the y-coordinate, this is equivalent to 7 The inversion appearing in these definitions, for example F − (e t ) = 1/F (e −t ), is for comparison to section 2. If we continued to Lorentzian signature without this inversion we would find coordinates on the Poincare patch which has t = 0 at the point of time reflection symmetry. In section 2, we instead used the Poincare patch pictured in figure 3, and these two patches differ by the null inversion x − → −1/x − (see equation (2.18)). 8 Note that there is nothing inconsistent about having non-local equations of motion for the boundary mode. The theory is still local, since it comes from a local two-dimensional action. saying that the only contributions to T y − y − come from the Schwarzian, i.e. from Hawking radiation.

JHEP04(2021)289
For the purposes of studying the black hole background solution, there is no obvious advantage to our Euclidean setup. The advantage comes when we want to apply Euclidean path integral methods, including replica wormholes, to the evaporating black hole, which we will come to in the next section.

Schwinger-Keldysh and asymmetric shockwaves
We can view the shockwave geometry as a solution defined on a Schwinger-Keldysh contour [58]. The contour in the complex τ plane is illustrated in figure 9. This is the contour on which the gluing function θ(τ ) is defined. In the Schwarzian equation of motion, there are singularities in the complex τ plane from the insertions of the ψ operators at iτ = y 1 , y 2 and −iτ =ȳ 1 ,ȳ 2 . At finite δ, the solution is smooth along the Schwinger-Keldysh contour. As δ → 0, there is a pinch singularity where the two ψ insertions move into the real-t line, which leads to non-analytic behavior in θ(τ ) in this limit.
The solution is obviously symmetric under t → −t. However, by moving the shockwave insertion to I − , we can move the outgoing shockwave to past infinity, and the solution at finite t becomes identical to an asymmetric solution with only an ingoing shock. This simply moves the singularities in the lower-half τ plane to τ = −i∞; as δ → 0, the solution for t > 0 is unchanged. Therefore we can expect to reproduce all of the physics of the asymmetric shockwave, including the Page curve, by taking this limit of the Euclidean solution. 9 9 The ingoing and outgoing shocks are entangled, so the purely-ingoing shockwave produced in this way is in a mixed state. We are assuming the entanglement entropy of ingoing and outgoing shocks is subleading compared to the entropy of the black hole that forms, as discussed below (2.42).

Replica wormhole equations
In section 3, we discussed the construction of the shockwave background from the Euclidean path integral. In this section, we write the equations for replica geometries at n = 1. The equations are implicit due to the non-local welding contributions, but can be simplified somewhat near n = 1. We set β = 2π. We will not use explicit features of the shockwave solution, so this analysis applies to any state created by a similar Euclidean path integral.

Replica geometry setup
For now, we will focus on a region R which has a single endpoint in the non-gravitating region. That is, on the Euclidean y-plane, including the right boundary point at y = 0. Equivalently we can take the complement region, Here b is complex so at the end we can continue to Lorentzian time. See figure 10.
To compute entropy by the replica method, one needs to first compute the replica partition function Z n = Tr (ρ R ) n for n copies of the system. The entanglement entropy is then computed by analytic continuation in n, In a QFT without gravity, Z n is given by a path integral on an n-sheeted Riemann surface M n . The original manifold is the quotient M = M n /Z n . The replica partition function can be viewed as a correlation function of twist operators Ψ|T n (y 1 )T −n (y 2 )|Ψ , where |Ψ is a replicated version of the original state, created by inserting the same operators in each copy. The dimension of the twist field in 2d CFT is ∆ T = c 12 (n − 1/n). In a theory with dynamical gravity, there are two crucial differences in the replica method. First, the geometry backreacts so that the replica manifold M n solves the gravitational equations of motion. This means that the quotient geometry M n = M n /Z n now depends on n, and it has conical defects at the twist insertions with angular identification 2π/n (so that M n is smooth). For n ∼ 1 this leads to the area term in the entropy. The second crucial difference is that higher topologies can contribute to the gravitational path integral. The topologies responsible for the island rule are Z n -symmetric replica wormholes [13,14]. On the quotient manifold, these solutions are realized by the insertion of extra, dynamical twist fields. Any number of dynamical twist fields can appear in the gravity region, so long as the total twist correlator (including the non-dynamical twist fields in the exterior region) is neutral.

JHEP04(2021)289
For n ∼ 1, the dynamical twist fields obey an equation of motion that extremizes the generalized entropy. Thus the twist fields become the endpoints of the island. A general argument for this based on the effective action for twist defects appears in [5,13,14,51], but in the island context, only the eternal black hole has been treated in detail [13,14]. Here we will consider the effects of additional operator insertions in the exterior region in order to see in detail how the replica equations lead to the island rule for the evaporating black hole.
For an interval with one endpoint, the important contribution to the path integral comes from a saddlepoint with one dynamical twist field at y = −a. See figure 11 for the replica manifold with n = 2. The goal is to write the equations of motion for an n-fold replica, then take n → 1 to determine the quantum extremal surface, a. As in section 3, we use a coordinate w in the gravity region and v in the flat region. The twist points are at and the quantum extremal surface is A = e −a as n → 1.

Finite n equation of motion
Notation. The background solution is the n = 1 replica geometry, so in this section the notation F 1 , G 1 , θ 1 , etc., denotes the solutions obtained in section 2. The Euclidean coordinatew is a global coordinate on the interior (gravitating) part of the replica manifold M n , and w is a coordinate on the interior part of the quotient manifold M n = M n /Z n . The function w n (τ ) is the gluing function, so the gluing in Euclidean signature is w = w n (τ ). The Euclidean coordinates used here are related to the Lorentzian coordinates of section 2 by w = 1/x − ,w = x + , so the background gluing function is w 1 (τ ) = 1/x 1 (τ ) with x 1 (τ ) the solution found in section 2. The stress tensor is denoted by T for general n, and T (1) for the n = 1 background. For a single interval, there is a uniformization map from the n-fold cover of the unit disk, with a twist field at w = A = e −a , to the unit disk itself. The map is an SU(1, 1) isometry taking A to the origin, followed by z → z 1/n . Denoting the uniformized coordinate byw, it is given byw In order to solve Einstein's equation, the curvature in the gravity region should satisfy R = −2. Therefore, the metric in the unit disk in the covering manifold is the standard hyperbolic metric 6) and the metric in the quotient unit disk with coordinate w is easily found from the coordinate transformation (4.5). It is in thew coordinate that the equation of motion is given by (2.11). That is, if we denote the boundary curve in thew plane byw n (τ ) and the boundary curve in the w plane by w n (τ ), the extrinsic curvature in the Schwarzian limit is where (4.8) The equation of motion that follows from varying the Schwarzian action is therefore In order to have a complete set of equations, we also need to determine the stress tensor. This is done in two steps. First we map to the z plane, using the conformal welding map z = F n (e y ).

JHEP04(2021)289
Next we need to find the stress tensor on the z-plane. Since z covers the whole plane with a single coordinate system, this can be calculated in principle by the usual CFT methods. There are contributions to T zz from the ψ insertions and from the twist operators. In general these do not decouple, so we cannot calculate T zz exactly in an interacting theory. Let us write where the first terms is the universal contribution given by equation (3.8), with the insertion points at The second term in (4.12) is non-universal -it encodes the Renyi entropy of the shockwave and therefore depends on the details of the CFT. Also, the metric in the z-plane is not flat and the Weyl factor gives an additional contribution to the stress tensor. The metric in the gravitational region in terms of z coordinate is given by (4.6), wherew is given by (4.5). Therefore, in this background the stress tensor is related to stress tensor in the flat metric dzdz as In the vicinity of the twist operator at z A = G n (w = A) in the flat metric, the stress tensor takes the form − 2πT twist,flat where and A = G −1 1 (z A ). Note that the double pole term in (4.16) vanishes in (4.17) due to the non-trivial Weyl factor in (4.14). Although we have focused on the shockwave state created JHEP04(2021)289 by a local ψ insertion, this discussion is general. The only difference in an arbitrary state is that T ψ zz is given by the expectation value in the background state, T ψ zz = ψ|T zz |ψ . To summarize, the gravitational equation of motion at finite n is determined as follows: find T zz , transform to the y-variable using (4.10)-(4.11), and plug into (4.9). Needless to say, even in situations where T zz is known, this equation is not easy to solve because the welding map depends on w n (τ ). It could perhaps be solved numerically, as was done for the eternal black hole in [63].
Assuming that T zz can be expressed as an analytic function of n, the final equations make sense at non-integer n. This is because we assumed replica symmetry and formulated the equations on the quotient manifold M n /Z n . These are believed to be the dominant contributions in the shockwave state away from the Page transition time, but there are other cases where non-replica-symmetric contributions to the path integral are important [14,38,64,65].
For the Lorentzian problem, we must allow B,B to be independent, and integrate the equation of motion on the Schwinger-Keldysh contour shown in figure 12. This is the same contour as the background geometry in figure 9 but we have also shown the singular sources at the locations of the twist operators. These will lead to branch cuts in w n (τ ), so that at finite n, the gluing function θ n (it) has both real and imaginary parts along the real-t part of the contour.

Equations for n → 1
The equations simplify for n ∼ 1 so we can be more explicit. The twist field stress tensor and R n terms have explicit factor of n − 1 as their coefficients. Hence, to the first order, they are evaluated on the background solution w n=1 = e iθ(τ ) . Moreover, given the solution
Putting everything together, with w n (τ ) ≈ w 1 (τ ) + δw where δw is O(n − 1), we find to first order (4.20) where (4.21) and z is the solution of the background welding problem. The right hand side in (4.20) could be considered as the source term which depends only on the background solution whereas the left hand side is only the function of δw. However, even near n = 1, the equations for a general background are non-local and they are not reduced to differential equations.

Derivation of QES from replica equations
In this section we derive the QES condition from the replica equations of motion for n ∼ 1. We will do this two different ways. First, we will derive the QES directly from the Schwarzian replica equations discussed in section 5.1. This was done for eternal black holes in [13,14]; here we will do it for a general background, then specialize to the shockwave. The other approach is to solve the equations of motion locally near the defect, as in [5,51]. We review this derivation (specialized to JT gravity) in section 5.2.
To derive the island rule, we also need to show that the entropy derived from the replica method is equal to the generalized entropy associated to the QES. This was demonstrated from the effective action for twist defects in [13,14]. In section 5.3 we give a simpler (but ultimately equivalent) argument based on the Ward identity for a CFT coupled to gravity.

QES from the Schwarzian equations
Our starting point is the replica equation of motion for w n (τ ) as n → 1, given in (4.20). Working around a general background w 1 , with w n = w 1 + δw, the equation is The flux i(T yy − Tȳȳ) has contributions from the background state, the twist operators, and the Schwarzian of the conformal welding map.

JHEP04(2021)289
It is well known that the Schwarzian theory has an SL(2, R) symmetry. Using this symmetry, the left hand side of (5.1) vanishes if we integrate it against the following SL(2) generators: where α ∈ {0, 1, 2}. The equations (5.2) follow from the relation ∂τ {wn(τ ),τ } w n (τ ) which can be used to show that integrand is a total derivative for α = 0, 1, 2. The strategy is to show that these identities, applied to the equation of motion, give the extremality condition. In Euclidean signature, the contour C of integration is the boundary of a unit disk, while in the Lorentzian setup, the integral is taken over the interface, i.e., over the Schwinger-Keldysh contour depicted in figure 12. We will give the argument in Euclidean signature and generalize to Lorentzian at the end. Energy conservation relates the stress tensor across the interface as Therefore, integrating matter flux against SL(2) generators, we have where in the second line we use the relationw n = 1/w n which holds along the interface. Since α ∈ {0, 1, 2}, the equation (5.4) shows that any source contribution for the holomorphic and anti-holomorphic stress tensors that are analytic inside the disk vanishes when they are integrated against SL(2) kernels. In particular, the conformal welding stress tensors proportional to {G n (w), w}| w=wn , and {Ḡ n (w),w} w=wn , drop out from (5.4). Similarly, matter sources such as operator insertions outside the unit disk do not contribute to (5.4). This matches with the fact that in the local argument for finding the QES, reviewed below, the only important terms in the stress tensor are the residues near the dynamical twist defects.
The QES condition is obtained by integrating an appropriate linear combination of these kernels. One way to guess the correct kernel is as follows. We expect the R-term in (5.1) to integrate to give the dilaton term in the QES equation, ∂ A φ, because of the relation derived in appendix C: Therefore we look for a kernel that satisfies φ r 2π Starting with a general linear combination of the integrals for α = 0, 1, 2 we find that this holds for the kernel There is a similar kernelKĀ obtained by taking conjugates that integrates to ∂Āφ. Here r is an arbitrary complex number. Another way to determine the correct kernel is starting from the Ward identity. We will try to design a kernel which, upon integrating against the flux T yy − Tȳȳ, produces the entropy term ∂ A S CFT in the extremality equation. The Ward identity in the CFT on the replica manifold relates the derivative of the entropy to the stress tensor [52,67]. The relation is Here primes are τ -derivatives. The first line is the replica calculation of the entanglement entropy; the second line is the usual Ward identity; in the third line we have gone from the global replica manifoldw to the quotient w (there is no Schwarzian because we are not doing a Weyl transformation to remove the conformal factor in the w-metric); the fourth line is the coordinate change w = e iθ . In the last line we used the background equation of motion to drop the (∂ n | n=1 n) term, and used the fact that w = 1/w along the contour of integration to combine the two integrals, and used the conservation of flux, discussed above, to rewrite the stress tensors in the y-coordinate.
Using the explicit coordinate change (4.5), assuming |w| = 1, and setting w = w 1 (τ ), the kernel in the last line is This is equal to K A if we assign r = 1 A(|A| 2 −1) . Therefore we have 14) The ambiguous term proportional to r in the kernel corresponds to the Ward identity for rescaling δw ∝w, and does not affect this integral. The identity (5.14) can also be checked by explicit evaluation of the integral. 11 11 In more detail: the stress tensor Tyy has three contributions, discussed around (4.17). The first is T flat zz , the stress tensor in the metric dzdz on the welding plane. The expansion of this stress tensor near the twist point is known [52,67], and its residue is proportional to ∂AS flat CFT . The second contribution is from welding; this drops out of the integral (5.14) because it is analytic inside the disk. The third contribution comes from the Weyl factor and accounts for the extra term in ∂AS CFT = ∂A S flat CFT − c 6 log(Ω(A)) .

JHEP04(2021)289
Combining everything, we found that integrating the equation (5.1) with kernels K A and KĀ yields the extremality conditions, For the Lorentzian problem where we have independent A,Ā and B,B, all of the integrals are done over the Schwinger-Keldysh contour shown in figure 12 rather than the unit disk, and the final equations are the same. For the entropy term, this follows by deforming equation (5.14) into Lorentzian signature by moving the point (A,Ā) and simultaneously deforming the contour to prevent any singularities from crossing the contour as we smoothly move (A,Ā). The same argument applies to the dilaton term; see appendix C.3 for a more explicit calculation.

Local derivation of the QES
We will now apply the methods of [5,51] to JT gravity to re-derive the QES condition by locally solving the equations of motion near the defect (see also the appendix of [15]). In this approach we directly use the dilaton equation of motion. In conformal gauge, the dilaton equations of motion in JT gravity are [41] where we set 4G N = 1. At n = 1, we can take the metric of a standard hyperbolic disk, (5.18) Expanding around n ∼ 1, we have dynamical twist defects inserted in the gravity region. We focus on a single dynamical branch point and determine its position by locally solving the equation of motion. We choose the coordinate (w,w) so that the branch point A is placed at the origin w =w = 0, and solve the equation (5.17) and its barred version near the origin to leading order in n − 1. 12 The quantities appearing in (5.17) are expanded around their background values as 19) and the first order equations are (5.20)

JHEP04(2021)289
The metric near the branch point on a single sheet of the replica manifold may be locally written as (5.21) since the new coordinate w = w 1/n that covers a neighborhood of the defect in the full manifold gives a smooth metric ds 2 ∝ d wd w. Accordingly, the conformal factor at first order is The dilaton can be generally expanded near the branch point as δφ ≈ a 00 + a 10 w + a 01w + a 11 ww + a 20 w 2 + a 02w 2 + · · · (5.23) We can relate the residue of the stress tensor δT ww to the derivative of the entanglement entropy as Comparing the residues of the simple poles in (5.24) and (5.25), we obtain the conditions for the QES where we have restored G N = 1/4. Let us make a brief comment on the relation with the arguments in the previous subsection using the contour integral. The derivative of the entanglement entropy can be expressed by a contour integral derived from the Ward identity as (5.9). By using (5.13), we can rewrite it as (5.27) in our gauge: A =Ā = 0. Now we plug the relation (5.25) into this integral and pick the residue at A, then we obtain the condition for the QES (5.26). The Ward identity used here is that for CFT in a fixed spacetime. The QES equation (5.26) itself can be viewed as the full Ward identity in the CFT coupled to gravity; the vanishing of the derivatives is the statement of diffeomorphism invariance. This is similar to the point of view in [51]. We will discuss the full gravitational Ward identity further in the next subsection.

Island entropy from the gravitational Ward identity
So far, we have derived the extremality condition, but we have not yet derived the island formula. The last step is to calculate the entropy from the on-shell action and verify (1.1). The entropy is obtained from the replica partition function by (5.28) This can be calculated directly from the action of the defect, and the result [13,14] is the island formula (2.47). Here we will give a different, simpler derivation from the Ward identity. The point is that once we have derived the extremality condition from the equations of motion, the entropy is automatically correct as well.
Let us first recall how the Ward identity is used to calculate the entropy in a CFT without gravity [52,67]. Consider the entropy of an interval with endpoints (A,Ā) and (B,B). This obeys the Ward identities 13 and similarly for ∂ B and ∂B. We used these relations above in writing the stress tensor near the defects in equations (5.9) The Ward identities can be integrated to find S. For example, for a CFT (not coupled to gravity) in vacuum on the flat dwdw plane, . This integrates to the well known formula S = c 6 log |A − B| 2 + const [52,67]. In a CFT coupled to gravity, S is only a function of the point (B,B) in the asymptotic region, since (A,Ā) is determined dynamically by the equations of motion. So there is no Ward identity for ∂ A , ∂Ā, but the relations for ∂ B , ∂B are unchanged. Thus (5.31) and similarly for ∂B. The matter stress tensor on the replica manifold is the same in gravity+CFT as it is in CFT, up to O((n − 1) 2 ). Therefore the right-hand side can be evaluated by the CFT Ward identity, and we have note that the kernel 1−Āw 1−|A| 2 vanishes at w = 1/Ā and is equal to one at w = A, so the effect of the kernel on the last line of (5.8) is to remove the anti-holomorphic contribution to the integral. This argument also uses the fact that once the conformal factors are included, the stress tensor has no double pole near the defect, as in (4.17).

JHEP04(2021)289
CFT calculation is that now A andĀ are not independent variables -they are functions of (B,B) determined by the QES condition. Note that for CFT, we had four Ward identities associated to A,Ā, B,B. For gravity we also have four Ward identities, two are given by (5.32) and the other two are the QES conditions. The last step is to integrate this equation. It is easy to check that the solution is the generalized entropy at the extremum, as in the island rule: To check this we simply take the derivative as needed. The other terms dropped out by the extremality condition. (This argument actually fixes S only up to an overall constant, which requires a direct evaluation of the action in one example). This derivation can be illustrated by a simple example: the QES in empty AdS 2 at zero temperature [12,13]. We refer to [13, section 4.1] for details of the island in this example. For the interval y ∈ [−a, b], with a, b > 0 so the left endpoint is in the gravity region and the right endpoint is outside, the generalized entropy is We are working at t = 0, so a =ā, b =b. The extremality condition ∂ a S gen = 0 gives the QES as Now there are two ways to calculate the entropy. The first is to apply the island formula, i.e., plug this value of a = a(b) into (5.35), The second is to integrate the Ward identity. We want to check that this gives the same answer. The Ward identity for gravity+CFT is (5.32), which in this example states with S CFT = c 6 log (a+b) 2 a . Integrating this equation gives (5.37) up to an integration constant, which can be fixed from the b → 0 limit.

Zero temperature limit
All of our calculations thus far have been at finite β, which is the inverse temperature of the black hole in the initial state. As emphasized in the introduction, this is only for generality and technical convenience. We can smoothly take β → ∞ to recover the island rule for an evaporating black hole that starts in vacuum. In this limit, the Schwarzschild-like y ± coordinates of section 2 become Poincare coordinates. In our fundamental equations for the background geometry (3.9) and replica geometries (4.9), the dependence on β can be reintroduced by dimensional analysis. Subsequently sending β → ∞ gives a Schwarzian equation that governs the gluing of a planar interface. We define the welding conditions on the plane by taking the β → ∞ limit of the welding conditions on a circular interface, so this ensures that all of the results on replica wormholes in sections 3-4 and the island formula in section 5 also carry over to this limit.

Factorization of the two-interval solution
So far we have discussed the setting of evaporating black holes. There is also an information paradox for eternal black holes [12], and replica wormholes can be applied in this context [13]. This setup is simpler than an evaporating black hole because it has no operator insertions in the definition of the state, but it has an extra complication: to see the information paradox, one must consider the entropy of two disjoint intervals, with one on each side of the black hole [12]. Otherwise the entropy is time-independent and there is no paradox.
It is reasonable to expect that at late times, the replica wormhole for this setup factorizes into two independent wormholes for the left and right regions, since this is a regime where the twist operators are in an OPE limit. This factorization property was advocated on physical grounds in [13, section 5] in order to make contact with the information paradox. In this section we will confirm factorization of the wormhole geometry explicitly for n ∼ 1 by explaining how to patch together two separate solutions of the Schwarzian equation, and recover the expected QES's.
When the replica manifold is branched along a single interval, with only one branch point in the gravity region, the topology of the gravity region in the replica manifold is unchanged -it is a disk for all n. By contrast when there are two branch points in the gravity region, the replica manifold has wormholes with higher topology. This introduces moduli, and therefore new equations of motion, beyond the Schwarzian equation. We will see however that in the factorized limit these extra equations are not necessary. We will start with the one-sided problem, then describe how at late times two copies of the solution can be patched together.

Single interval geometry in the eternal black hole
All equations discussed in sections 4-5 for the one-interval case also apply to the situation with no shockwave operator insertions. This is the eternal black hole. The geometry is given by w 1 (t) = e −t (we set β = 2π), F = G = id.
The endpoints of the interval in (σ, τ ) coordinates are given by (σ a , τ a ), where −σ a , σ b > 0 and τ a is a complex number. The problem has a boost symmetry and therefore we can always set τ b = 0. See figure 13. The Lorentzian time and Euclidean time are related by τ = it. The equation of motion for the replica manifold is the same as equation (4.20). Defining δM by the expansion the equation for the perturbation around n = 1 is φr . H is the Hilbert transform which acts on Fourier modes as H[e imτ ] = −sgn(m)e imτ . This is the equation governing the geometry of the replica wormhole solution for a single interval near n = 1. Solving (6.3) by matching the three Fourier modes e iτ , 1, e −iτ , one finds the position of point σ a according to the island rule which gives τ a = 0 and the condition [12,13] 2κ . . The full answer is δM = δM + + δM − . An example is plotted in figure 14 in Lorentzian signature.
There are several features of this solution in Lorentzian signature that we want to highlight: • The series expansion as a sum over Euclidean Fourier modes e imτ diverges at the lightcones of the twist operators. This leads to singularities in Lorentzian signature, which have been regulated by an i shift in the plot. The i prescription corresponds to the contour in figure 12 (without any ψ insertions).
• Since the series diverges, it must be summed before continuing past the lightcone. This changes the naive behavior in an essential way. For example, positive Fourier modes e imτ with m > 0 decay as e −mt as t → +∞. After doing the sum and analytic continuation around the lightcone singularity, we find that the sum over positive modes also decays in the past: δM + ∼ e κt as t → −∞ . (6.6) This is the crucial fact we will use below to glue two one-interval solutions into a two-interval solution.
• To obtain this solution we first set a to the QES. In other words, we first solve for the m = 0, ±1 Fourier modes. If we do not impose this condition, keeping a general, then we can still solve the equation of motion for all of the other Fourier modes, but we find δM diverges as t → ±∞. Therefore imposing regularity at early/late Lorentzian times is equivalent to the extremality condition.

Two interval geometry in the eternal black hole
Now we consider the geometry near n = 1 for two intervals in the eternal black hole. The coordinates in terms of (σ, τ ) are given by and the setup is shown in figure 15. The entropy of region R without considering the island grows linearly with time t b , S Hawking ∼ πc 6β t b . This follows from the fact that the wormhole grows linearly with time [68,69]. This indefinite entropy growth of a region is the Hawking paradox for the eternal black hole. However, an island dominates after the Page time, and the island prescription S(R) = min{S Hawking (R), S island (R)}, (6.8) gives the unitary Page curve. At late times, the twist operators are in an OPE limit, so the matter entropy factorizes as In addition the island prescripion sets t a = t b . The dilaton contribute to the entropy obviously just adds the two endpoints, so the upshot is that the generalized entropy with the island at late times is given by twice the single-interval answer. We want to argue that factorization at the level of geometries also holds. Intuitively this is expected because the points P 1 , P 2 are far from points P 3 , P 4 . Therefore, the matter stress tensor and the effects of conical defects factorizes to two single interval answers.
To be more concrete, let us denote the perturbative Schwarzian of the single interval case due to points [P 1 , P 2 ] and [P 3 , P 4 ] by δM [P 1 ,P 2 ] and δM [P 3 ,P 4 ] . We are now solving the Schwarzian equation along a contour similar to figure 12, but with two Lorentzian pieces -one on the right side of the black hole, and the other on the left. We must show that along this contour the perturbation is additive, This is indeed a solution to the equation of motion. It was shown around (6.6) that in the single-interval case, on the Lorentzian part of the contour, δM [P 1 ,P 2 ] is exponentially small for t t a . Therefore for t a > t Page , The same is true for δM [P 3 ,P 4 ] . Therefore we can add these two solutions, gluing them together in the Euclidean regime where both perturbations are negligible. Practically, this factorization implies that for evaluating observables like the gravitational action or the matter effective action on a Schwinger-Keldysh contour, one could alternatively approximate the exact two intervals answer by the sum of two Schwinger-Keldysh contours for single interval geometries as shown in figure 16. Let us make a comment about number of equations needed to find all QES conditions. In an exact analysis of two intervals geometries, only a subset of QES conditions at points P 1 , P 3 are expected to be given by integrating the boundary curve using three SL(2) kernels constructed in section 5. The other set of equations are supposed to follow from variation with respect to moduli of the replica manifold. However, when we consider the limit corresponding to the factorized two interval geometry to single interval solutions, we find all QES conditions from integrating the boundary curve. The reason for finding extra equations is that for each single interval solution, we impose the boundary conditions at t → ±∞ separately from the other single-interval solution. Therefore the extra equations, which should in principle come from the equations of motion for the moduli, came for free from regularity, i.e. from requiring the solution to factorize. A better understanding the two interval replica manifold is an interesting question that we leave for future work.
To summarize, we have seen that two copies of the single-interval replica wormhole can be patched together to write the replica wormhole for the late-time, two-interval problem.

JHEP04(2021)289
It follows that we have also derived the QES's for two intervals in the eternal black hole, and therefore provided further justification for the island rule in this version of the information paradox. Let us note that the welding terms in the Schwarzian equation were essential for this to work, since these terms were responsible for the early-time decay of the Schwarzian perturbation.

JHEP04(2021)289
Now let us evaluate the cross-ratio η, which is related to the z − coordinate, at the positions of QESs obtained by the analyses in section 2. As we send B → ∞, i.e, v − B = 0, the left QES A is placed at the bifurcation point v − A = ∞. Moreover, by using u 0 1 we obtain where for the first approximation we used the condition for the left QES A and B → ∞, and for the second approximation we used the condition for the right QES (2.53) and the asymptotic expansion of the map x(t) in (2.34). Therefore we need to study the entanglement entropy in the regime η ≈ 0.
where the block admits a series expansion around η = 0 as Therefore the entanglement entropy can be decomposed as where S fact.
is the factorized part of the entropy we used to find the QES in the previous section, and S non-fact.
QFT is the remaining "non-factorized" part on the entropy S QFT (I ∪ R) defined respectively as where in the last line we have used η → 0, which projects onto the vacuum block. As we can see, in general the entanglement entropy itself doesn't factorize completely due to S non-fact.
QFT even in the limit η → 0 due to the Virasoro descendants of the vacuum. Nevertheless we can show the factorization of its derivative (2.48) in the regime of interest by using the QES condition obtained from S fact. QFT . First we consider the derivative with respect to x − A in the limit η → 0. As we can see from the expression (A.8), the factorized part S fact.
QFT gives the leading singularity in the limit η → 0 as

JHEP04(2021)289
while the non-factorized part gives terms subleading in η. Therefore in the regime of interest u 0 1 ⇔ η ≈ 0 we can confirm the approximation The non-factorized part S non-fact.

QFT
(η) (after taking η → 0 limit) depends only onz = z + = e y + , so it is independent of x(t) and u 0 . Thus the dependence on . By using the asymptotic expression of where S non-fact.

QFT
(η) is non-singular and independent of x(t) and u 0 . Therefore using and confirm the relation (A.14) Example: free Dirac fermion. We explicitly demonstrate the factorization of the entanglement entropy using the free Dirac fermion. The entanglement entropy of the free Dirac fermion is given by [70] S fermions ( As we did above, we divide the entropy into the factorized part end non-factorized part as First we compute the derivatives with respect to x − A as (A.18) The leading contribution comes from the singularity η ≈ 0 in the factorized part S fact. fermions . Thus to compute the QES, we can ignore the non-factorized part of the entanglement entropy.
The derivatives with respect to x + A are (A. 19) and fermions for u 0 e −κ(t−L)/2 1, thus the factorized part gives the leading contribution. and fermions ≈ 0 at the QES in the κ → 0 limit. The variation of S non-fact.
fermions with respect to v − A may be evaluated at the QES A as Thus we confirmed that S non-fact. fermions doesn't change the position of the right QES A : v − A = 0. The variation with respect to v + A is evaluated as Since we have y + A − y + B < 0 and e −y + A ∼ O(e −1/κ ) at the QES, this gives zero up to a tiny correction O(e −1/κ ) we can neglect in the κ → 0 limit. In this way we can check the factorization (2.48).

B Shockwave solution at finite δ, small E ψ
In this appendix, we describe the shockwave solution to leading order in E ψ = h ψ /δ at finite δ. In this limit, the form of welding is known and the mixing between positive and negative frequency modes is tractable. For simplicity, we consider the zero temperature black holes for this discussion. After the analytic continuation τ → it, the equation (3.21) in the limit β → ∞ becomes These equations were also derived in [13, appendix C]. The numerical solutions are plotted in figure 17 for different values of δ. The full δS is determined by δS = δS + + δS − . The region where positive and negative modes overlap gets smaller as δ decreases.
As we increase the shock energy h ψ /δ, there are corrections to the welding term in the differential equation (B.1). However, as argued in section 3, these corrections are controlled by the mixing between positive and negative modes and if take δ → 0 first, the leading welding term is given by (B.1), and for the finite shock energy we find (3.35) and (3.36).

C Dilaton from boundary curve
In this appendix, we demonstrate that on any background w 1 (τ ) created by an arbitrary number of operator insertions outside the gravity region, the dilaton φ has a simple form (C. 2) The dilaton has a well-known closed form solution (2.8) in terms of the stress tensor. The goal here is to write the dilaton only in terms of the boundary curve. In subsection C.1 this is derived in Euclidean signature. The derivation in C.1 is general and relies only on the analytic properties of the stress tensor inside the gravity region. For solutions like sharp shockwaves which satisfy the no-mixing condition in real time, we compute the dilaton directly in terms of the boundary curve x 1 (t) in C.2. Finally, using the Schwinger-Keldysh contour, we show in C.3 that the Euclidean method is equivalent to the Lorentzian calculation.

C.1 Euclidean
Here we evaluate the dilaton at an arbitrary point (A,Ā) in the w-coordinate inside the unit disk. The argument is easier if we apply a SL(2) transformation to the disk and work instead on a plane. Let us define a new coordinate Z = X + iY from the relation − Z(τ )−Z 0 Z(τ )+Z 0 = w 1 (τ )−A 1−w 1 (τ )Ā , where Z 0 is an arbitrary constant that we assumed for simplicity is real and Z 0 < 0. The metric is the hyperbolic metric in the Poincare coordinates The boundary curve is (X(τ ), Y (τ )). The intrinsic metric condition on the boundary curve gives The equations for the dilaton in the Z-plane are given by Euclidean version of (2.8) as written in [41]. The closed form solution for these equations is is the homogenous solution. Changing the lower limit of the integral in (C.5) corresponds to a change in the homogeous part of the dilaton which would be fixed by the value of the

(C.7)
If the operators creating the geometry are all inserted outside the gravity region, each components of the stress tensor T ZZ , TZZ are analytic functions in the gravity region. In the Schwarzian limit, the boundary lies at X = 0 and therefore the real part of the stress tensor can be written in terms of the imaginary part, (C.8) The kernel appearing in (C.8) is nothing but the Poisson kernel for harmonic functions in two dimensions. The imaginary part of the stress tensor is the energy-momentum flux and it is related to the boundary curve Z(τ ) = iY (τ ) + O( ) as . (C.9) Combining (C.7), (C.9) and (C.8), the dilaton can be written in terms of the derivatives of the boundary curve. Doing the integration by parts and dropping the boundary terms, we find (C.10) Note that Z(τ ) is purely imaginary. By a SL(2) transformation, the dilaton in a general point (Z 0 ,Z 0 ) is written as (C.11) To make sure there is no homogenous term in (C.11), let us take Z 0 → X(u) + iY (u), Z 0 → X(u) − iY (u). The equation (C.11) then becomes where we used lim α→0 α 2 (r 2 +α 2 ) 2 = π 2α δ(r). This shows that (C.11) has the right boundary condition and in particular, there is no extra homogeous term in the final expression. Once the transformation − Z−Z 0 is applied back to (C.11), we find (C.1) which is the dilaton written in the disk coordinates. Note that computing (C.1) by residues would require additional assumptions on analytic properties of w 1 (τ ) in the complex τ -plane.
Here it is only assumed that w 1 (τ ) is a smooth real function and therefore (C.1) cannot be simplified further. Figure 18. The integration is performed along the null ray from x + 0 to x + . We change the integration variable from x to the boundary time t ≡ x −1 (x).

C.2 Lorentzian -direct calculation in shockwave
In the Lorentzian discussion for computing the generalized entropy, we used a particular form for the dilaton expressed in terms of the boundary curve (2.12). In this subsection, we derive this equation from the expression (2.8).
As we saw in section 3, the full stress tensor for t > 0, δ → 0 has T x + x + = 0, T x − x − = 0. We will focus on such cases, for which the dilaton is given by We take x + 0 = x(0) by using the gluing map x(t), which is placed before the shockwave, and we have T x + x + (x + 0 ) = 0. We evaluate the integral by changing the integration variable from x to t ≡ x −1 (x) (see figure 18) and by using the Schwarzian equation of motion at time t: (C.14) Using the following identity ∂t{x(t),t} x (t) = ( 1 x ( x x ) ) , we can compute this integral as where in the last equality we used the initial conditions x(0) = 1, x (0) = 2π/β, x (0) = (2π/β) 2 , x (0) = (2π/β) 3 . One can check that for an eternal black hole with a boundary curve x(y + ) = e 2π β y + , the dilaton is given by the expected answer:

C.3 Schwinger-Keldysh contour
Let us show that the dilaton evaluated by (C.11) in the shockwave background is equal to the direct computation done in C.2. Since the expression for the dilaton was derived for an arbitrary background in Euclidean signature, it is expected that the Lorentzian version of the formula is also given by Schwinger-Keldysh contour (C.15). Here the boundary curve parametrized with τ . The Euclidean variables in C.1 are related to Lorentzian ones in C.2 by the analytic continuation it = τ, Z(it) = ix(t). For the shockwave background, the only singularity in the integrand R 1 comes from residue of (C.11) at τ = ix −1 (x + ). See figure 19. By computing the residue we find This is exactly (2.8), so the two expressions match. When we integrate R 1 against the SL(2) kernels to find the QES, all boundary terms from integration by parts drop out since the Schwinger-Keldysh contour is closed. By repeating the argument in section 5.1, we find (5.5) is valid when the island appears in Lorentzian signature.

D Details of replica geometry for one interval in eternal black hole
Starting from equation (6.3), we expand the right hand side in terms of Fourier modes Open Access. This article is distributed under the terms of the Creative Commons Attribution License (CC-BY 4.0), which permits any use, distribution and reproduction in any medium, provided the original author(s) and source are credited.