Superluminal chaos after a quantum quench

Thermal states holographically dual to black holes in Einstein gravity display maximal Lyapunov growth as well as"butterfly effect cones". We study these effects in highly non-equilibrium states, obtained from an initial thermal state by the sudden injection of energy. We do this by computing out-of-time-order correlators (OTOCs) in BTZ-Vaidya spacetimes, which describe transitions between black holes at different temperatures. If both pairs of boundary operators appearing in the OTOC are inserted before the energy injection, we recover standard results, with butterfly effect cones displaying a light-cone structure. But when one pair of operators is inserted before and the other pair after the energy injection, the Lyapunov growth saturates the chaos bounds set by the local temperatures and the butterfly effect cones"open up", becoming superluminal, albeit in a way that does not violate causality. In the limiting case, in which the initial state is the vacuum, Lyapunov growth only starts after the energy injection. Our computations of the OTOCs are phrased in terms of gravitationally interacting particles, where fields are treated in a geodesic approximation and the eikonal phase shift is expressed in terms of stress tensors and shock waves associated to geodesics.


Introduction
In recent years, a connection between chaotic behavior of strongly coupled holographic theories and gravitational physics in AdS has been identified, leading to the conjecture that thermal states with holographic duals are the most chaotic systems [1]. This idea was made precise in the context of a quantum mechanical bound on chaos [2], which is saturated by thermal CFT states dual to AdS black holes [3][4][5][6][7].
In thermal systems at inverse temperature β, chaos can be diagnosed by a specific type of out-of-time-order correlator (OTOC) between "simple" hermitian operators V and W (which can be localized in spatially extended systems, or involve a small number of degrees of freedom in matrix models) [8]. In a chaotic system, the OTOC is expected to vanish at late times, independently of the choice of operators V and W [7]. The OTOC is closely related to the commutator squared (CS), which is used as a diagnostic for the presence of chaos in quantum systems because of its relation to the classical Poisson bracket and the exponentially divergent trajectories expected in classical chaos [8,9]. A quantum system is said to be chaotic if the commutator squared grows exponentially in time with characteristic Lyapunov exponent λ L and where N corresponds to the number of degrees of freedom of the system (we restrict our attention to systems without spatial extent for the moment). More precisely, this exponential behavior, which is related to the phenomenon of fast scrambling, is expected for times t d t t * , where t d ∼ β is the dissipation time and t * ∼ λ −1 L log N the scrambling time [1]. The scrambling time gives the time scale needed for a perturbation involving a few degrees of freedom to spread among all the degrees of freedom of the system [1,3]. Any operations performed after the scrambling time on a subset of the degrees of freedom will be unable to reconstruct the information of the perturbation.
The connection between the OTOC (1.1) and the squared commutator becomes clear after expanding ( (1.4) and noting that the first two terms factorize in a thermal state after a dissipation time. For example, since it can be thought of as the two-point function of W (t) in the thermal state perturbed by the insertion of operator V (0) [6]. If the energy injected by V is small, the state will relax and the expectation value will approach the thermal value multiplied by the norm of the state [6], essentially because all time-ordered correlations decay to products of expectation values in a thermal state after the dissipation time [10]. This information can be used to simplify (1.2) in the regime t t d , Thus, the normalized commutator squared D(t) is simply related to the OTOC after the dissipation time in thermal systems. In addition, the characteristic late-time vanishing of the OTOC for chaotic systems is equivalent to saturation of the normalized commutator squared after the exponential Lyapunov regime. The above considerations apply for quantum systems without a notion of spatial locality. If Hamiltonian interactions are local in space, i.e., if they decay for operators that are spatially separated, then there are two scrambling time scales. The first is associated to fast scrambling among local degrees of freedom and can be characterized as described above. The second type of scrambling is associated to spreading of quantum information in space and in some cases, a butterfly velocity v B can be associated to this propagation [5], roughly speaking measuring the rate of growth of the region that must be tracked to reconstruct an original perturbation with a given accuracy. In nonrelativistic systems the butterfly velocity is bounded by Lieb-Robinson bounds constraining the growth of commutators in discrete Hamiltonian systems to be at most ballistic [11,12]. In relativistic theories, it would be natural to expect that the butterfly velocity would be bounded by the speed of light. Note though that if the region that must be measured to reconstruct an original perturbation grows faster than light, it does not imply that this can be used to send information faster than light. Consider an astronaut on the moon who keeps track of laser dots of various colors originating from Earth; these dots can move faster than light on the surface of the moon, and the area needed to keep track of all of them can also grow faster than light. But of course this does not allow astronauts on the moon to transmit information faster than light.
Unlike in classical systems, there is a bound on the Lyapunov exponent in quantum mechanics. For thermal systems, the bound is simply set by the temperature [2], while some generalizations thereof have been recently studied in [13,14]. This bound is actually saturated in conformal field theories (CFTs) by thermal states dual to AdS black holes. The OTOC of such thermal CFT states at inverse temperature β can be computed in two distinct ways. In [6], the correlator was computed using CFT methods, by taking advantage of the exponential map between a thermal state and the vacuum state and making use of vacuum block techniques. In the gravitational bulk interpretation [7], the chaotic behavior can be understood in terms of a more conventional time-ordered high-energy bulk scattering of particle excitations created by the insertion of the operators V and W . In this picture, the OTOC is the overlap of bulk in-and out-states created by the boundary insertions. This can be seen by first viewing the OTOC as the overlap Ψ |Ψ of the two states |Ψ ≡ V (t − , 0)W (t + , x)|0 and |Ψ ≡ W (t + , x)V (t − , 0)|0 , (1.8) with |0 the state dual to the background bulk geometry with no excitations of the W, V operators, and t − < t + . The state |Ψ can be understood as follows. As a first step, W is Ψ |Ψ (c) Figure 1: (a) The bulk picture of the in-state. This state has a simple interpretation at early times. (b) The bulk picture of the out-state. This state has a simple interpretation at late times. (c) The bulk picture of the OTOC is the overlap of the in-state (a) and the out-state (b), which corresponds to bulk scattering. Figure adapted from [7].
inserted at the boundary of the empty geometry at time t + . In order to apply the operator V at time t − on this state, the perturbation W is evolved backwards in time using the bulk-to-boundary propagator. Subsequently, the V perturbation is inserted. The state as a whole has a simple bulk picture at early times (Fig. 1a), and is therefore interpreted as an instate. Expressing the same state at late times would require solving an involved scattering problem. A similar picture can be obtained for |Ψ , where the operators are inserted in the opposite order. The V operator is inserted at time t − on the empty background geometry and is evolved forward in time in order to apply the W operator at a later time t + . The V perturbation creates field excitations that fall into the black hole, and the W operator can be applied without disturbing the V field excitations too much. This state has a simple form at late times (Fig. 1b), and is interpreted as an out-state. A Fourier decomposition into momentum states allowed to derive a formula for the OTOC that describes the correlator as a 2 → 2 bulk scattering of highly energetic particles created by the perturbations V and W [7]. This scattering is illustrated in Fig. 1c.
In both methods, the OTOC yields under the assumption ∆ W ∆ V and where the parameter ε has been introduced to regulate the divergences associated to the coincident insertions of operators. For a small enough separation t + −t − −|x| between pairs of operator insertions, the commutator squared (1.6) displays an exponential regime, with Lyapunov exponent λ L = 2π β . It saturates the bound (1.7) and hence supports the idea that black holes are the fastest possible scramblers.
In the holographic bulk interpretation, the Lyapunov behavior arises from the blueshift that such particles are subjected to, due to the gravitational acceleration in the black hole background. Initially small perturbations around the black hole background end up energetic enough to create strong gravitational fields known as shock waves, causing them to interact strongly with other particles travelling near the black hole. This phenomenon is at the core of the chaotic behavior of the dual holographic CFT. From (1.10), the spatial region for which D(t + −t − , x) is non-negligible can be seen to grow ballistically with butterfly velocity v B = 1, In addition, expression (1.9) also encodes the decay of the OTOC at very late times t β 2π ln G −1 N , This decay is controlled by the lowest quasinormal mode frequency ω QN = 2π β ∆ V associated to the light operator V of conformal dimension ∆ V [15]. This shows that dissipative effects eventually take over the chaotic Lyapunov growth.
In this paper, we investigate the chaotic behavior of a system in the process of thermalization by holographically computing OTOCs in a thermal state undergoing a global quench. On the gravity side, this state corresponds to an initial BTZ black hole with inverse temperature β − (BTZ − ) which, after a global energy injection at the boundary of the spacetime that creates an infalling shell, turns into a black hole with inverse temperature β + < β − (BTZ + ), as depicted in Fig. 2. These bulk states are known as BTZ-Vaidya spacetimes.
We first extend the computation of the OTOC in the eternal BTZ black hole, which was reviewed above, to the case of dynamically growing (forming) black holes. As for the eternal black hole, boundary perturbations can be propagated into the bulk and the OTOC can be interpreted as an overlap of an in-and an out-state. When all operators are inserted before (after) the shell, the shock-wave scattering process happens in the BTZ − (BTZ + ) part of the spacetime, as can be seen in Fig. 3a. In this case, we recover the standard Lyapunov behaviors in thermal states, but we are able to extend previous results to time separations between the inserted operators of the order of the dissipation time and smaller.
However, when one pair of operators is inserted before and the other after the quench, the shock-wave scattering process often takes place in BTZ − , as illustrated on Fig. 3b (this happens if the transverse separation between the operators is sufficiently small; the precise condition will be given in Section 3.3). 1 This fact modifies the (spatial) behavior of the Figure 3: (a) When all the operators are inserted after the quench, the situation is very similar to the eternal BTZ case, as the scattering happens in the BTZ + part of the spacetime. (b) When one pair of operators is inserted before and the other after the shell, the scattering typically happens in BTZ − , which modifies the behavior of the OTOC.
OTOC, as compared to the thermal equilibrium case. There is an unexpected consequence: right after the quench, the butterfly velocity is set by the ratio of the second to the first horizon radius, and is superluminal. However, as we will see, this effect relaxes before it can violate causality and the overall spread remains within the lightcone. Moreover, two exponential behaviors develop with separate Lyapunov exponents, set by each of the temperatures, that multiply the time separation between the insertion of the quench and the insertion of each of the operators in the part of the spacetime associated to that temperature. From the holographic perspective, this feature is expected in view of the blueshift that particles experience before the scattering, which is the origin of the Lyapunov exponents in black hole geometries. By taking the zero-temperature limit of the initial state (β − → ∞), we extend our results to CFT states which start in the vacuum and thermalize after a global quench, allowing us to study the onset of chaos in dynamically forming black holes. As a cross-check we also recover the vacuum CFT result obtained by [6], for all the operators inserted before the quench. The plan of the paper is as follows. In Section 2, we derive for a general spacetime a position-space WKB analogue of the overlap integral that computes the OTOC in [7]. We apply this formula to BTZ-Vaidya spacetimes in Section 3 and discuss our results in Section 4. The appendices contain technical material regarding the properties of geodesics, the WKB approximation and shock waves.
occur in different parts of the geometry depending on the details. 6

OTOC in asymptotically AdS spacetimes
Let us consider an asymptotically AdS background spacetime on which two massive Klein-Gordon real scalar fields φ V and φ W propagate. We will derive a general formula for the out-of-time-order correlator (OTOC) where the insertion points X 2 , X 4 lie in the future of X 1 , X 3 or are spacelike-separated from them. As the radial location of the insertion points X i are taken to the boundary of AdS space, (2.1) computes (after a radial rescaling) the associated correlator in the dual conformal field theory defined on the boundary [16][17][18][19]. We will carry out all our calculations in such a near-boundary limit of the insertion points so that the dependence of the correlation functions on time and the transverse space can be interpreted as also computing these dependences in the dual field theory.
In the limit that all insertion points are close to the conformal boundary, one can make use of the WKB (geodesic) approximation for propagators of the scalar fields. Indeed, in any asymptotically AdS space, a timelike geodesic gets more energetic as one of its endpoints gets closer to the conformal boundary. This is exactly the regime in which we expect the WKB approximation to be reliable. 2 In addition, at leading order in an inverse energy expansion e −1 1, such highly energetic geodesics follow light rays and their gravitational backreaction can be exactly computed in the form of shock waves. It will hence be useful to equip the considered background spacetime with coordinates X = (u, v, x), where u, v are advanced and retarded lightlike coordinates and x labels the transverse direction.
We want to find a formula applying to arbitrary asymptotically AdS spacetimes. In the general case, it is not clear whether the method of [7] may be directly applied, especially in time-dependent backgrounds for which asymptotic in and out Hilbert spaces are not always available. However, the intuition that the OTOC is dominated by the scattering of highly energetic particles [21] should still hold. 3 In this section, treating field excitations in the WKB (geodesic) approximation, we obtain a general position-space formula for the OTOC (2.1) in terms of the scattering of highly energetic geodesics or point-like particles. We will proceed as follows, viewing the OTOC as the overlap Ψ |Ψ of the two states and First, we give a simple representation of |Ψ at earlier times and of |Ψ at later times, in the spirit presented in the Introduction. For this, we use Green's second identity to represent φ V (X 1 ) on a future null slice Σ u 0 and φ W (X 4 ) on a past null slice Σ v 0 . As a result, the overlap is purely expressed in terms of time-ordered four-point correlators (Fig. 4). Second, we write down the gravitational path integral representation of these time-ordered four-point correlators. Treating field propagators in the WKB approximation and taking the classical gravity limit, we show that their evaluation boils down to the computation of classical quantities: highly energetic geodesics and their associated shock wave gravitational fields. The method developed here shares many similarities with [22][23][24][25].
From out-of-time-order to time-ordered correlators We first consider the state |Ψ . We want to express the field operator φ V inserted at X 1 in terms of the same operator inserted at points X u 0 = (u 0 , v, x) belonging to a null slice Σ u 0 of constant u = u 0 such that any X u 0 ∈ Σ u 0 is either in the future of X 2 , X 4 or spacelikeseparated from them. For this we first choose a volume V containing all insertions, which is enclosed by the union of Σ u 0 , parts of the conformal boundary and of the past horizon, as well as transverse surfaces of constant x = x ±∞ being infinitely far away from all operator insertions. The projection of ∂V onto the (u, v)-plane is indicated in light green in Fig. 4a. We then apply Green's second identity, to the following solutions of the Klein-Gordon equation, The retarded Green function G R V (X, Y ) vanishes for (X − Y ) 2 > 0 and for X 0 < Y 0 , and satisfies We consider Dirichlet boundary conditions at the conformal boundary 4 and no incoming initial conditions at the past horizon for the scalar field φ V such that G R V (X; X 1 ) vanishes there. On the enclosing surface ∂V, this Green function has therefore support only on Σ u 0 ⊂ ∂V and we getφ giving a representation ofφ V (X 1 ) in terms of fields on the null surface Σ u 0 . We can now write the state |Ψ by simply multiplying the right hand side of (2.7) byφ W (X 2 ): Figure 4: (a) The operator φ V (X 1 ) is represented on a future null slice Σ u 0 in terms of the operator φ V (X u 0 ) integrated over the point X u 0 ∈ Σ u 0 . The enclosing surface ∂V is indicated in light green. (b) The operator φ W (X 4 ) is similarly represented on a past null slice Σ v 0 in terms of the operator φ W (X v 0 ) integrated over the point Their evaluation boils down to the computation of highly energetic geodesics indicated by the straight colored arrows, and their associated gravitational fields. As shown in the next sections, we anticipate that the overlap Ψ |Ψ is dominated by the contribution of a single pair of such geodesics.
Choosing a null surface Σ v 0 of constant v = v 0 such that any X v 0 ∈ Σ v 0 is either in the past of X 1 , X 3 or spacelike-separated from them, and applying a similar construction to the operator φ W (X 4 ), we end up with Taking the overlap of (2.8) and (2.9) then leads to an expression where only time-ordered four-point correlators appear:

WKB and classical gravity limits
The time-ordered correlation functions appearing in (2.10) may be represented as Feynman path integrals. For example, one can write where T indicates time-ordering. The correlators of φ V and φ W have factorized in this expression because we are taking the action to be the Einstein-Hilbert action plus the Klein-Gordon actions of the scalar fields, so that the fields φ V and φ W only interact through gravity. Note that we can safely drop the time-ordering operator T in (2.11) due to the respective location of the points X 2 , X v 0 , X u 0 , X 3 . Under the assumption that the insertions points X 1 , ..., X 4 are sufficiently close to the conformal boundary, one can now evaluate these propagators in the WKB approximation. As reviewed in Appendix B, the WKB approximation for a propagator in a background spacetime with metric g µν is, for X in the future of Y , (2.14) with S V being the proper time along the timelike geodesic going from Y to X. Thus, is the future-directed momentum of this geodesic at the point X. This phase is also the action for a point-like particle following the geodesic trajectory, and varying it with respect to the background metric g µν gives the associated stress-tensor which acts as a source for perturbations of the background metric. Plugging (2.13) back into the path integral (2.11), and taking the saddle point approximation in the classical gravity limit, we hence get The action to be evaluated on-shell is the one of gravity sourced by two highly energetic (lightlike) geodesics. At leading order in an inverse energy expansion, the on-shell metric is µν is the original background metric satisfying prescribed boundary conditions and h µν ≡ h V µν + h W µν describes the gravitational shock waves associated to the two lightlike geodesics. It should be noted that (2.17) would be an exact solution of Einstein's equations if only one geodesic were considered [26]. Here however, we neglect nonlinearities associated to the dynamical interaction of the shock waves with each other because these are higher order in the Newton constant. The on-shell action (2.16) can be expanded to linear order in h µν , giving [7,24] and where T µν ≡ T µν V + T µν W is the stress-tensor associated to the two geodesics. The linear piece is also known as the eikonal phase shift δ, and will play a central role in the next sections. Hence, we have Note that the phase factor e iS EH g (0) µν is common to all observables computed on this background geometry, and is therefore irrelevant. The overlap integral (2.10) is then It is left to apply the WKB approximation to the retarded propagators. For X in the future of Y , these are This final formula is the WKB version of the overlap integral presented in [7]. It has the important advantage of being expressed in position-space and is therefore easily applicable to any background spacetime. Using (2.7) and (2.22), one can write similar expressions for the scalar propagators, Hence, the crucial ingredient that prevents the OTOC (2.23) from factorizing into the product of the two propagators (2.24-2.25) is the eikonal phase shift e iδ which encodes the gravitational interaction of highly energetic field excitations. Again, as the radial location of the insertion points of the bulk operators (X i ) are taken to the boundary of AdS space, a radial rescaling of (2.1) computes the OTOC in the dual conformal field theory on the boundary [16][17][18][19]. We will perform all calculations in a nearboundary limit of the insertion points. So, the dependence of the correlation functions on time and the transverse space also compute these dependences in the dual field theory.

OTOC in BTZ-Vaidya
We want to apply the general OTOC formula (2.23) to the three-dimensional BTZ-Vaidya background spacetime. For this, we need to compute geodesics whose endpoints are close to the conformal boundary together with the eikonal phase shift δ associated to the gravitational interactions of these geodesics.
The three-dimensional BTZ-Vaidya spacetime can be described in retarded coordinates (v, r, x) whose associated metric is The surface v = v s corresponds to the location of a shell of dust falling in from the asymptotic boundary. This spacetime is also simply understood as a gluing of two planar BTZ black holes along the surface v = v s , of respective horizon radii r − and r + . This is made manifest by performing the change of coordinate v = t + 1 2r h ln r−r h r+r h , For simplicity, we refer to these two planar black hole regions as BTZ − and BTZ + , respectively. See Fig. 2. Note that we work in units where the AdS radius has been set to unity. Note also that unlike the retarded coordinate v, the time coordinate t is not continuous across the shell except in the limit r → ∞ corresponding to the location of the conformal boundary. We will also often make use of the advanced-retarded coordinates in which the BTZ metric takes the form In this gauge, r = r(u, v) is a function of the advanced and retarded coordinates through (3.3). We specify the four insertion points of the OTOC in Schwarzschild coordinates (t, r, x): Following Section 2, we choose two null surfaces Σ u 0 and Σ v 0 and study the timelike geodesics connecting these surfaces to the insertion points X 1 , X 2 , X 3 , X 4 . The requirement that no point X v 0 ∈ Σ v 0 can be in the causal future of the insertion points X 1 , X 3 means that this surface necessarily lies in BTZ − , and we can choose it to be any constant v = v 0 < t − surface. Similarly, Σ u 0 should be chosen such that no point X u 0 ∈ Σ u 0 lies in the causal past of X 2 , X 4 . This is achieved by choosing a surface of some constant u = u 0 in BTZ − and continuing it across the shell in BTZ + . The precise value of u 0 does not actually matter, only the existence of such a surface. These surfaces are illustrated in Fig. 4b-4a.
When working in the limit of small UV regulators + , − 1, the WKB phases in the overlap formula (2.23) simply differ by such that As reviewed in Appendix A, the momenta k t (X 1 ) and k t (X 2 ) of the ingoing and outgoing geodesics are directly related to their conserved energies, 5 and depend implicitly on the geodesic endpoints X 1 , X u 0 and X 2 , X v 0 . Hence, the overlap (2.23) and the coincident two-point functions (2.24-2.25) simplify to These integrals over the points In the spirit of [7], we perform these by stationary phase approximation in the regime of large field masses 6 m V − ε and m W + ε (in units 5 If X 2 lies in BTZ + , the outgoing geodesic crosses the shell and a jump in energy occurs at the crossing point. In that case, we need to distinguish between the energy e − out of the BTZ − segment of the geodesic and the energy e + out = −k t (X 2 ) of the BTZ + segment of the geodesic. 6 It is natural to consider ε and − , + on the same footing since they all act as UV regulators. In that case, the stationary phase approximation can be performed within the m V , m W 1 WKB regime.
14 where the AdS radius is set to unity). We will show that the functional dependence of e in and e out on the integration variables x, x , u, v are such that the saddle point of the integrands corresponds to a single pair of highly energetic radial geodesics. In that case, the normalized OTOC takes a particularly simple form, where it is implicit that all insertions happen close to the conformal boundary, r i = ε −1 . Hence, the normalized OTOC equals the eikonal phase shift associated to the gravitational interaction of this single pair of highly energetic radial geodesics. In order to determine which pair of geodesics, we have to find the saddle point of e in and e out . This is done in Sections 3.1-3.2. Having determined the dominant pair of geodesics, we compute the eikonal phase shift associated to their gravitational interaction in Section 3.3.
We start by stating some useful general results regarding timelike geodesics in a planar BTZ black hole, whose derivation may be found in Appendix A. It is known that any locally AdS space can be viewed as a patch of a hyperboloid embedded in a higher-dimensional flat space, and that timelike geodesics are circles on this hyperboloid. However, not any pair of timelike-separated points (t a , r a , x a ) and (t b , r b , x b ) with t b > t a in Schwarzschild coordinates, may be joined by circles. This is only possible for those satisfying If the constraint (3.14) is satisfied, L ab is the length of the geodesic connecting these points. In this case, the conserved energy and transverse momenta associated to the Killing vectors ∂ t and ∂ x are where the dot refers to the derivative with respect to the geodesic proper time. In addition, the v component of the momentum of an ingoing geodesic and the u component of the momentum of an outgoing geodesic are both equal to where one should keep in mind that the value of r varies along the geodesic.

Ingoing geodesics
In order to find the saddle point of the energy e in along the null surface Σ u 0 , we consider those geodesics that connect the insertion point Assuming that such a geodesic exists, its energy and transverse momentum are found from (3.15-3.17), with As previously explained, the integrals over x and v in (3.10-3.11) are performed by stationary phase approximation. Since the x integral localizes at values of x satisfying cos L in = 0 or ∂ cos L in ∂x = 0. Subsequently performing the v integration, we similarly look for saddles along the v direction. If cos L in = 0, there is no saddle point along the v direction provided that u 0 = t 1 . Hence, we consider the remaining case ∂ cos L in ∂x = 0, which is equivalent to x = x 1 and characterizes radial geodesics. Importantly, not every value of v corresponds to a point X u 0 = (u 0 , v, x 1 ) that can be connected by a timelike geodesic to the insertion point X 1 . From (3.14), the existence of such a geodesic requires In this limit, the timelike geodesic connecting X 1 to the unique saddle point X u 0 = (u 0 , t 1 , x 1 ) follows a radial null ray of constant retarded coordinate given by (3.24), and the geodesic length can be shown to satisfy From (3.16-3.20), we also deduce its tangent momentum, Figure 5: The region accessible to timelike geodesics connected to the insertion point X 1 is shown in yellow. It puts a restriction on the set of points X u 0 = (u 0 , v, x 1 ) ∈ Σ u 0 that may be connected by timelike geodesics to the insertion point X 1 = (t 1 , ε −1 , x 1 ) lying close to the conformal boundary.

Outgoing geodesics
We similarly want to find the saddle point of the energy e out along the null surface Σ v 0 . For this, we need to distinguish between two types of outgoing geodesics, depending whether the insertion point X 2 = (t 2 , ε −1 , x 2 ) lies in BTZ − or in BTZ + . In the former case, the treatment is parallel to the one of ingoing geodesics, such that the associated saddle lies at x = x 2 , u = t 2 + O(ε). We therefore turn to geodesics that connect X 2 in BTZ + to points Such a timelike geodesic necessarily crosses the shell at a point X s = (v s , r s , x s ). For that reason, we study segments of that geodesic on either side of the shell separately, and subsequently impose appropriate junction conditions at the shell that will fix the coordinates of the intermediate point X s .

BTZ + segment
Assuming that it exists, the energy and transverse momentum of the timelike geodesic connecting a point X s on the shell to the insertion point X 2 are found from (3.15-3.17), with (3.29) The superscript '+' makes explicit that these quantities are associated to the BTZ + segment of the considered outgoing geodesic. Similarly, we will use the superscript '−' to denote quantities associated to the BTZ − segment of the outgoing geodesic.

BTZ − segment
We turn to the timelike geodesic segment in BTZ − connecting a point X v 0 ∈ Σ v 0 to X s . The energy and transverse momentum of the timelike geodesic are found from (3.15-3.17), (3.32) Continuity of the retarded coordinates (v, r, x) across the shell and a careful analysis of the geodesic equation in BTZ-Vaidya [27] imply that both segments of a BTZ-Vaidya geodesic have to satisfy the following junction conditions: 34) or equivalently [28] j These junction conditions may be shown to coincide with those originally derived in [29] by extremizing the total geodesic length L out = L + out + L − out . They can be used to fix the coordinates r s , x s of the intermediate point X s in terms of the endpoints X v 0 and X 2 , but are very difficult to solve in general. In Appendix A, we show that the velocity of the geodesic at the shell is given bẏ .

Saddle point
In order to perform the x integration in (3.10) and (3.12) by stationary phase approximation in the regime of large field masses, we need to find the associated saddle point of e + out . Here, we argue that a saddle point lies at x = x 2 for any fixed value of u. From (3.27), one sees that e + out depends on x only through r s (x ) and cosh r + (x 2 − x s (x )). Taking the square of the first junction condition (3.35), one can express cosh r + (x 2 − x s (x )) in terms of r s (x ) and cosh r − (x s (x ) − x ), The second junction condition (3.36) allows to express r s (x ) in terms of cosh r − (x s (x )−x ), Hence, the energy e + out associated to the BTZ + segment of the outgoing geodesic is shown to depend on the integration variable x only through the combination cosh r − (x s (x ) − x ), and a saddle point therefore exists for x s (x ) = x . In turn, this implies that the associated outgoing geodesics have zero transverse momentum j − out = j + out = 0 and that the saddle point lies at x = x 2 .
After having set x = x 2 , we look for saddles of e + out along the u direction. Only those intermediate points X s = (v s , r s (u), x 2 ) satisfying the constraint (3.14) may be connected by timelike geodesics to the insertion point X 2 . To first order in ε 1 and using (3.29), this constraint reads − εr 2 Hence, only values of u satisfying the above constraint should be considered when looking for saddle points of e + out . Taking the derivative of (3.27), we get Hence, a saddle occurs when the term in square brackets vanishes. The corresponding value of r s (u) satisfies the constraint (3.41), and is given to first order in ε 1 by One can use the junction condition (3.36) to determine the corresponding value of u. In the limit ε 1, this problem considerably simplifies by noticing that the timelike geodesic connecting X 2 to the saddle point described above necessarily satisfies From (3.32) and (3.43-3.44), one can therefore infer the associated value of u, In deriving the location of the saddle point, it has been implicitly assumed that the functions r s (u, x ) and x s (u, x ) were differentiable. A posteriori, we have explicitly checked that this was indeed the case, at least at the saddle point. Moreover, we numerically studied the dependence of e + out on x and u, which suggested that the saddle we found is the only saddle.
As we pointed out, the timelike geodesic connecting X 2 to the saddle point described above follows a radial null ray in the small ε limit. In particular, its BTZ − segment has constant advanced coordinate given by (3.45). Using (3.43), we may simplify the expressions of the conserved energy (3.27) and initial velocity (3.37b) of the BTZ + geodesic segment, Then, the junction condition (3.36) allows us to evaluate the momentum (3.18) of the BTZ − segment of the outgoing geodesic,

Eikonal phase shift
Having found the pair of highly energetic geodesics corresponding to the saddle point of (3.10-3.12), we are now in position to give the formula of the eikonal phase shift δ which encodes the interaction of one geodesic with the gravitational shock wave sourced by the other geodesic, and vice versa. When those interactions happen in BTZ − , we show in Appendix C that the eikonal phase shift is given by 20 where x − , x + ,ū,v label the position of the saddle point. If all insertions take place before the shell (t − , t + < v s ), this interaction automatically happens in BTZ − and is nonzero only Provided that this causality condition is satisfied, the eikonal phase shift is found from (3.49) and (3.24-3.26), If the second insertion takes place above the shell (t + > v s ), the eikonal gravitational interaction may happen either in BTZ − or BTZ + . In Appendix C.3, we show that the interaction happens in BTZ − if the transverse separation between the operators is sufficiently small: (3.52) The eikonal phase shift is then obtained from (3.49), (3.26) and (3.48), When (3.52) is not satisfied, one of the two gravitational interactions takes place in BTZ + and its computation goes beyond the scope of this paper. In particular, it requires to evolve the shock wave produced by the ingoing geodesic across the shell, i.e. in the BTZ + region. For this one would have to solve the linearized Einstein's equations around the BTZ + background spacetime and subject to appropriate junction conditions at the shell.
As an interesting special case, we take the limit r − → 0 which corresponds to considering an AdS 3 -Vaidya background instead of BTZ-Vaidya. The associated eikonal phase shift is worked out in Appendix C.2 for the case where gravitational interaction occurs entirely in Poincaré AdS 3 , and is given by (3.54) Similarly to the BTZ-Vaidya case described above, two situations should be distinguished. If all insertions take place before the shell (t − , t + < v s ), this interaction automatically happens in Poincaré AdS 3 and is nonzero only for |∆x| < ∆t. (3.55) Provided that this causality condition is satisfied, the eikonal phase shift is found from (3.54) and (3.26), If the second insertion takes place after the shell (t + > v s ), the eikonal interaction does not necessarily happens in Poincaré AdS 3 . In Appendix C.3, we show that it is only the case if where trajectories of the interacting radial geodesics follow constant advanced and retarded coordinates on the AdS 3 side, respectively, The eikonal phase shift is then obtained from (3.54), (3.26) and (3.48), When (3.57) is not satisfied, one of the two eikonal interactions happens in BTZ + and we do not have the appropriate formula at our disposal, as in the more general BTZ-Vaidya case.

Results
We now discuss the results obtained in the WKB and stationary phase approximations, for the various cases studied in the previous sections. The expectation value of the commutator squared D(t + , x + , t − , x − ), which is the object of interest in the diagnosis of quantum chaos, is entirely controlled by the OTOC in the regime of time separation larger than the dissipation time scale, Two inverse 'local' temperatures may be associated to a quenched state with BTZ-Vaidya dual, The inverse temperature β − is the one of the unperturbed thermal state before the quench is applied at t = v s , while β + is the inverse temperature associated to the quenched state after it relaxed to equilibrium. The dissipation time scale is set by one of these two inverse temperatures in the following way. If all insertions happen before the quench (t − , t + < v s ), then we are in the configuration described in [7] and the approximation (4.1) requires t + − t − β − . If the first and second pairs of operators are inserted before and after the quench (t − < v s < t + ), respectively, then the above approximation applies in the regime t + − v s β + . In the previous sections, we derived the simple formula (3.13) for the OTOC in the limit of large scalar field masses m V , m W 1 and to first order in a small G N expansion, such that However, the WKB approximation has limitations that should be pointed out. Indeed, it corresponds to a high-frequency perturbative expansion for which propagating bulk fields only feel the local geometry of the background spacetime. Any physical effect associated to low-frequency or long-wavelength modes are nonperturbative from this perspective and are therefore missing from our description. In particular, this means that we cannot account for dissipation effects. On the one hand, these should manifest themselves for small time separations that coincide with the region of invalidity of (4.1). As explained in the Introduction, dissipation is also responsible for the saturation of the CS at very late times. The saturation behavior controlled by the quasinormal decay rate has been explicitly obtained in [6,7] in case of a thermal state, but is missing from (4.3) due to the use of the WKB approximation. Hence, we should not trust the result for time separations such that δ 1.
The precise expression of the eikonal phase δ has been given in (3.51), (3.53) for the BTZ-Vaidya background, and in (3.56), (3.59) for the AdS 3 -Vaidya background as a particular limiting case. We use it to produce a series of contour plots in Figs. 6-7 depicting the behavior of δ = D(t + , x + , t − , x − ) for different values of the shell insertion time v s , indicated by the position of the black horizontal line. We set t = t + , x = x + and t − = x − = 0 without loss of generality.

BTZ-Vaidya
This regime is therefore characterized by the maximal Lyapunov exponent allowed by a thermal state at inverse temperature β − , namely λ − = r − = 2πβ −1 − [2]. To first order in a small G N expansion and in the regime δ 1, equation (4.4) is in agreement with the results presented in [6,7]. For t > v s corresponding to the upper parts of Figs. 6a-6c, new features appear. On the upper left side, one observes contour lines with approximately constant slope. This behavior is recovered from (3.53) together with (3.24) and (3.45) Hence, in this regime we observe the maximal Lyapunov exponent associated to the second equilibrium temperature, λ + = r + = 2πβ −1 + . From this equation, we see that a Lyapunov growth e λ − vs has been accumulated before the quench, and that the subsequent Lyapunov the eikonal gravitational interaction happening in BTZ + , in order to derive an expression appropriate to that region of the plot. Next to this hatched region, contour lines are bending upwards and we further comment on this behavior in Appendix C.4. Finally, it can be noted that the transition between the regimes t < v s and t > v s is continuous. Indeed, we have such that (3.53) smoothly transitions to (3.51).
As discussed in the Introduction, significant values of the commutator squared characterize regions where quantum information has been efficiently scrambled. The corresponding contour lines (δ = 1) have been highlighted in red in Figs. 6a-6c. In this regard, two effects should be distinguished. The first one is the fast scrambling among local degrees of freedom, while the second one is the spatial spreading of information. Depending on the values of the parameters, local scrambling is achieved either before or after the quench is applied. In the former case, it happens at the scrambling time Moreover, the boundary of the spatial region where information has propagated and has been significantly scrambled first grows ballistically at the speed of light (v B = 1) until the quench is applied, and subsequently faster than the speed of light Hence, the butterfly effect cone "opens up"; in other words chaos is superluminal. Eventually, the butterfly velocity v B decreases as contour lines are bending upwards, next to the hatched region. This situation is illustrated in Fig. 6a. As a result, there is no conflict with causality since the commutator squared vanishes outside of the lightcone. When local scrambling is achieved after the quench is applied, the scrambling time is given by In this case, the boundary of the spatial region where information has propagated and has been significantly scrambled grows faster than the speed of light from the onset (v B = r + /r − > 1) before eventually slowing down close to the hatched region. This situation is illustrated in Figs. 6b-6c. Again, there is no conflict with causality. In fact, this result agrees to first order in a small G N expansion with the one derived in [6] by CFT methods. In particular, the commutator growth is quadratic in t − |x| rather than exponential, and has been qualified as slow scrambling characteristic of vacuum in [6]. It should however be emphasized that the dissipation time is formally infinite in Poincaré AdS 3 , such that the approximation (4.1) is not reliable, and no direct relation between the eikonal phase δ and the commutator squared is expected. For t > v s corresponding to the upper parts of Figs. 7a-7b, the eikonal phase quickly behaves as Hence, the commutator growth is controlled by the maximal Lyapunov exponent allowed by the late-time equilibrium temperature, λ + = 2πβ −1 + , and starts from t = v s onwards. As in BTZ-Vaidya, contour lines are bending upwards. This latter behavior extends all the way from the left of the figures to the hatched regions, leaving no linear/ballistic regime. Indeed, the time and space dependences in (4.10) do not occur through the combination t − |x|, such that there is no regime in which the spatial spreading is ballistic (v B ≈ const), in contrast to the BTZ-Vaidya case. Neglecting the contribution from the initial period of slow scrambling, the scrambling time characterizing the successful scrambling of information among local degrees of freedom is simply given by The boundary of the spatial region where information has propagated and has been significantly scrambled grows faster than the speed of light from the onset and gradually decreases without ever violating causality.

Conclusions and outlook
In this paper, we computed the OTOC holographically in a state far from equilibrium, by investigating the chaotic behavior of a thermal state undergoing a sudden global injection of energy. The OTOC was constructed by considering the scattering of geodesics propagating in the BTZ-Vaidya spacetime, expressed in position space and using the WKB approximation. With this method, we recover the standard Lyapunov behavior in the case where the two pairs of operators are inserted before the quench. However, when one pair of operators is inserted before and the other pair after the quench, and the insertion locations are sufficiently close in space, the Lyapunov behavior saturates the bound set by the local temperature while the spatial suppression factor is set by the smallest temperature, causing the butterfly effect cones to become superluminal, albeit in a way that preserves causality. There are several further directions that could be investigated. One next step would be to fill the unexplored (hatched) regions of Figs. 6 and 7. In the gravitational description, this would require solving Einstein's equations for the evolution of the shock wave created by an ingoing geodesic travelling in BTZ − past the infalling shell of matter. Another interesting question to ask is whether distinct quantities characterizing the spreading of information, such as entanglement and information velocities, also feature superluminal behavior. These velocities have been studied in a non-equilibrium context in [30].
One shortcoming of our method is that it does not describe the saturation behavior of the commutator squared, as the WKB approximation only accounts for high frequency modes. In [7] the states |Ψ and |Ψ were decomposed in a momentum basis and the knowledge of the bulk-to-boundary propagators on the horizons made it possible to perform at least one of the four necessary integrals exactly. In this way, the Lyapunov behavior of the commutator squared as well as its saturation, where the lowest quasinormal mode dominates the integrals, were captured in a single formula. Inclusion of the saturation behavior in our setting would involve the bulk-to-boundary propagators in Vaidya spacetimes, which have been studied numerically in [31][32][33] for the case of AdS-Vaidya.
It would also be interesting to understand how different aspects of our discussion, and especially the superluminal chaos, arise in the dual field theory picture. In [34,35] a framework that allows for the computation of n-point functions in a quenched state starting from the vacuum was proposed, assuming domination of the Virasoro vacuum block exchange at large central charge. Note that the field theory calculation of the OTOC in a thermal CFT ensemble [6], performed using large central charge Virasoro block techniques, was sensitive to both the Lyapunov and the saturation regions. Therefore, formulating the computation of the OTOC in a quenched state in these terms might provide another way to tackle several of the questions raised by this paper.
leading to the induced metric In planar BTZ, timelike geodesics have conserved energy e and transverse momentum j associated to the Killing vectors ∂ t and ∂ x , such that the geodesic equation simply readṡ We can express the conserved quantities in terms of any two points (t a , r a , x a ) and (t b , r b , x b ) with t b > t a along a geodesic, by using (A.2) together with the coordinate transformation (A.5), Note that the M summation runs over any subset of three independent coordinates. In the main sections, we use the momentum of the geodesic in the u − v directions defined in (3.4). For an ingoing geodesic,ṙ < 0 such that For an outgoing geodesic,ṙ > 0 such that (A.14) We also need the velocity of an outgoing geodesic at its initial endpoint,

B WKB approximation
In this appendix, we describe the WKB approximation [37,38] to Klein-Gordon propagators satisfying The leading (quadratic) order in m of (B.1) leads to the so-called eikonal equation which may be written alternatively as Taking the covariant derivative of this equation and using ∇ µ k ν = ∇ ν k µ , it may be shown that k µ satisfies the geodesic equation and we conclude that the phase S(X) is simply related to the action of a timelike geodesic ending at the point X, However, there is still freedom in choosing whether the tangent momentum k µ is future-or past-directed. For the time-ordered Wightman propagator one has to choose k 0 < 0 or equivalently k 0 > 0, which corresponds to a geodesic traveling forward in time from Y to X. As we use it in the main text, we also write the expression of the retarded propagator, At next-to-leading (linear) order in m, (B.1) leads to The solution for A in (d+1) spacetime dimensions is given in terms of the Van Vleck-Morette determinant ∆ [39,40], Here, L(X) = −S(X) is the length of the geodesic. In pure AdS d+1 , the Van Vleck-Morette determinant is given by [41] ∆(X) = L(X) sin L(X) d . (B.12)

Amplitude in BTZ-Vaidya and validity of the WKB approximation
We now discuss the WKB amplitude of the Wightman propagator (B.7) in the case of a BTZ-Vaidya background. Let us first consider X, Y as being both on the BTZ − side of the shell. As we have seen, the amplitude may be expressed entirely in terms of the geodesic distance L − (X) between X and some other reference point along the BTZ − segment of the geodesic, As the reference point can be chosen arbitrarily, the most general form of the solution is Then, the fact that the source is located at Y fixes the integration constant L 0 , For example, we can conveniently choose The other integration constant C 0 should be fixed by normalization of the amplitude. We turn to the case where X is lying in BTZ + . The amplitude still assumes the same form, where now L + (X) is the geodesic distance between X and any other point along the BTZ + segment of the geodesic. One then uses continuity of the wave function at the shell in order to fix L 1 and C 1 . This leads to In the bulk of the paper, we do not explicitly need these formula. However, they are useful in assessing the validity of the WKB approximation. Indeed, it is assumed that variations of A(X) are adiabatic with respect to the time scale m −1 , Looking at (B.10) and (B.12), we see that this breaks down for a short geodesic distance between X and Y , i.e. when The geodesics corresponding to the saddle points of Section 3 have geodesic length π/2 up to O(ε) corrections, such that the WKB approximation for the associated propagators is valid. In addition, we avoid the evaluation of propagators whose endpoints are connected by a geodesic traveling very close to the shell. In such a case, moving X on the other side of the shell might completely change the nature of the connecting geodesic and therefore lead to a discontinuity in the amplitude (B.10). This would of course lead to a breakdown of (B.19).

C Shock waves and eikonal phases in AdS 3
Shock wave solutions in pure AdS have been described in various dimensions and coordinate systems in [42][43][44][45]. We provide an alternative derivation of shock waves in any locally AdS 3 spacetime by solving Einstein's equations sourced by a null particle, starting with its embedding in 4-dimensional flat space where isometries are most straightforwardly realized. This allows to make direct contact between null geodesics and their gravitational field. We then apply this general construction to the case of radial null geodesics in planar BTZ. The shock wave associated to radial null geodesics in Poincaré AdS 3 is obtained as a limiting case.
First, we recall how to obtain the stress-tensor of a point-like particle. The action of a massive particle following a trajectory X µ (s) is given by where s is any parameter along the geodesic. The particle stress-tensor is then obtained by functional differentiation with respect to the metric, We further choose a set of coordinates x µ = (u, v, y) where u, v are lightlike, such that X µ (s) = (U (s), V (s), Y (s)). Changing of integration variable s → U directly leads to where we specialized to the proper time s = τ in the last line. In a limit where the particle follows a null ray at V (τ ) = v 0 , Y (τ ) = y 0 , the stress-tensor becomes We now consider a locally AdS 3 spacetime, and work in flat embedding coordinates (A.1). We eliminate the coordinate X 3 such that the AdS metric becomes We then define lightcone coordinates, such that We consider a particle following a null trajectory along V = 0 and Z = Z 0 , such that its stress-tensor is given by Note that the momentum k V is constant over the trajectory. We then look for the associated shock wave geometry from the following ansatz, The only nontrivial component of Einstein's equations which needs to be solved is, using Vδ(V) = −δ(V), This can be recognized as the equation for a scalar propagator of unit mass in the transverse hyperbolic space spanned by X = (X 2 , X 3 ), This equation appeared in [46,47], and its solution is well-known to be 7 where we have used equation 9.121.24 of [49] for the last step, and where the two-point invariant in transverse hyperbolic space is The shock wave solution (C.9) describes the gravitational field of a highly energetic particle in any locally AdS 3 spacetime, assuming that one knows the coordinate transformation that maps it to the flat embedding coordinates. In addition, one can make use of AdS isometries in order to map other null trajectories to the one described above.

C.1 BTZ
We give the shock wave geometry associated to a radial null geodesic in planar BTZ whose trajectory is characterized by V =V , x = x in in Kruskal coordinates. We write the transformation between Kruskal coordinates and some embedding coordinate systemX i , The following AdS isometry maps the particle trajectory to the null ray V = Z = 0, where a ≡ ln 1−V 1+V . Hence, we have From this coordinate transformation and the general solution (C.9), we find the V V component of the shock wave, where we identified the momentum of the geodesic as Note that without the step function θ(V −V ), (C.17) would also contain the shock wave associated to the reflected geodesic at the conformal boundary, namely the outgoing geodesic with trajectory U = −V −1 , x = x in . The other components of h µν are not required for what follows. ForV = 0, we recover the shock wave solution discussed in [7], In addition, the shock wave produced by an outgoing particle following the U =Ū , x = x out trajectory is found by interchanging U and V in (C.17).

Eikonal phase shift
The eikonal interaction of an ingoing particle at V =V , x = x in and an outgoing particle at U =Ū , x = x out may now be computed. This interaction will be nontrivial only if the outgoing particle crosses the shock wave of the ingoing particle, and vice versa. For this, there should exist a point (Ū , V, x out ) along the outgoing geodesic such that the argument of the delta function in (C.17) is zero, i.e. such that Since V ∈ 0, −Ū −1 for such an outgoing geodesic, we immediately getŪV > −1, and a nontrivial interaction requires V <V cosh r h ∆x +ŪV 1 +ŪV cosh r h ∆x < − 1 U , (C. 21) or equivalently, cosh r h ∆x < −ŪV + (ŪV ) −1 2 .
(C. 22) In this case, using (C.4) and (C.17), one has δ = S on−shell 1 +ŪV cosh r h ∆x 1 +ŪVŪV + cosh r h ∆x 1 +ŪV cosh r h ∆x Π(ξ 2 ) (C.23b) = 2πG N m in m out k in V k out U e −r h |∆x| 1 +ŪV e r h |∆x| 2 . (C.23c) For eitherŪ = 0 orV = 0, the eikonal phase shift discussed in [7] is recovered, The shock wave geometry associated to an ingoing radial null geodesic in Poincaré AdS 3 with trajectory characterized by v =v, x = x in , can be obtained from (C.17) in the limit r h → 0. Before taking this limit, one first needs to change to advanced-retarded Schwarzschild coordinates (u, v, x) defined in (C. 25), such that the vv component of the shock wave is The limit r h → 0 may now be taken, resulting in where evaluation of the transverse propagator gives Π(ξ 2 ) = 1 2 Putting these together, we find the following shock wave component This result agrees with the shock wave solution associated to a radial geodesic in Poincaré AdS 3 presented in [45].

Eikonal phase
The eikonal interaction of an ingoing particle at v =v, x = x in and an outgoing particle at u =ū, x = x out may now be computed. This interaction will be nontrivial only if the outgoing particle crosses the shock wave of the ingoing particle, and vice versa. For this, there should exist a point (ū, v, x out ) along the outgoing geodesic such that the argument of the delta function in (C.33) is zero, i.e. such that Since v ∈ (−∞,ū) for such an outgoing geodesic, a nontrivial interaction again imposes the condition (C.26). If it is fulfilled, one has This could have been obtained directly from (C.27) in the limit r h → 0.

C.3 BTZ-Vaidya and AdS 3 -Vaidya
The shock wave geometry associated to an ingoing radial null geodesic in BTZ − with trajectory characterized by v =v < v s , x = x in , is given below the shell by the solution (C.17). However, part of the shock wave reaches the shell at v = v s , at which point the solution (C.17) is not valid anymore. Instead, one would have to evolve the shock wave solution across the shell through Einstein's equations, subject to appropriate junction conditions at the shell. This goes beyond the scope of this paper. Hence, with the shock wave solution in BTZ − , we are only able to compute the eikonal phase shift associated to the gravitational interaction of this ingoing geodesic with an outgoing radial geodesic following the u =ū, x = x out trajectory in case that the latter interacts with the shock wave before crossing the shell. The condition for this to happen requires a slight modification of (C.21), namelyV <V cosh r − ∆x +ŪV 1 +ŪV cosh r − ∆x < e r − vs , (C. 36) or equivalently, cosh r − ∆x < 1 + e −r − (ū+vs−2v) e −r − (vs−v) + e −r − (ū−v) .
(C. 37) This condition is of course stronger than (C.21) or (C.26). If it is satisfied, the eikonal phase shift is still given by (C.27). If not, we do not have the appropriate formula at our disposal. The same comments apply to the computation of the eikonal phase shift in AdS 3 -Vaidya.
In particular, it agrees with the Poincaré AdS 3 result (C.35) if the following condition is satisfied,

C.4 Upward bending of contour lines
In Figs. 6-7 of Section 4, it has been noticed that contour lines are bending upwards close to the hatched regions, where we do not have the appropriate formula describing the interaction of the outgoing geodesic with the gravitational shock wave produced by the 38 ingoing geodesic. However, we have a formula for the other contribution to the eikonal phase shift (2.19), namely the one corresponding to the interaction of the ingoing geodesic with the gravitational shock wave produced by the outgoing geodesic, given by Since this interaction always occurs in BTZ − , it is given by half of (3.49). In Fig. 8, we plot contour lines associated to this known contribution. Strikingly, there is a transverse distance |x max | on which all contour lines accumulate. Indeed, (3.45) implies that the positionū of the outgoing geodesic is bounded, such that (C.26) implies the existence of a maximal transverse separation |x max | =ū max beyond which this contribution is trivial. Although x max lies in the hatched regions of Figs. 6-7, part of the accumulation behavior is retained in the form of upward bending of contour lines. As a final comment, although both contributions to the eikonal phase are identical outside of the hatched regions as they both correspond to gravitational interactions occurring in BTZ − , we may expect them to largely differ inside the hatched regions.