Local quenches, bulk entanglement entropy and a unitary Page curve

Quantum corrections to the entanglement entropy of matter fields interacting with dynamical gravity have proven to be very important in the study of the black hole information problem. We consider a one-particle excited state of a massive scalar field infalling in a pure AdS3 geometry and compute these corrections for bulk subregions anchored on the AdS boundary. In the dual CFT2, the state is given by the insertion of a local primary operator and its evolution thereafter. We calculate the area and bulk entanglement entropy corrections at order ON0,\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$ \mathcal{O}\left({N}^0\right), $$\end{document} both in AdS and its CFT dual. The two calculations match, thus providing a non-trivial check of the FLM formula in a dynamical setting. Further, we observe that the bulk entanglement entropy follows a Page curve. We explain the precise sense in which our setup can be interpreted as a simple model of black hole evaporation and comment on the implications for the information problem.


The big picture
The presence of entanglement is an essential and ubiquitous feature of quantum systems. It is often quantified using entanglement entropy S A . Given a quantum system in a state |ψ , one starts by partitioning its space into subsystems A and A c . The entanglement entropy between the subsystems is then defined as the von Neumann entropy S A ≡ −tr[ρ A log ρ A ], where ρ A ≡ tr A c [|ψ ψ|] is the reduced density matrix associated with the region A. It is interesting to study this quantity for matter fields interacting with gravity, particularly, in situations where the state is fully dynamical, e.g., for evaporating black holes. Indeed, one of the observations that led to the black hole information paradox is that the semiclassical result for the entanglement entropy of matter fields outside an evaporating black hole initially in a pure state, keeps increasing with time [1], a behaviour that contradicts the expectations from unitarity [2]. This result could be interpreted in two ways: either, one should accept that the information is lost in black holes, challenging the rules of quantum mechanics; or one should understand the breakdown of unitarity as an artifact of the semiclassical approximation.
Recently, a resolution of this tension was put forward by suggesting a new prescription for the computation of the entanglement entropy of a quantum gravitational system in the semiclassical approximation. It amounts to adding a new term dubbed the island contribution which accounts for the possibility of new saddle points in the semiclassical path integral [3] S A (t) = min Here I denotes the so-called island, ∂I is its boundary, G is the Newton's constant of the gravity theory and S eff A is the effective, or coarse-grained entropy, which is calculated using the semiclassical density matrix, i.e., the ρ eff A that describes the quantum dynamics of matter fields living on a classical geometry. On the other hand, the entropy in the lefthand side of (1.1) uses the fine-grained density matrix of quantum fields, including gravity. The problem is that there is no simple way to calculate the latter from a first principle calculation since this task would require starting from a complete and consistent quantum theory of gravity.
There is some evidence, owing to studies in low-dimensional models of conformal matter interacting with gravity, that such a formula may be correct [4]. These models study gravity in AdS spacetime and make use of holography. One of the precise entries in holographic dictionary is the statement that the entanglement entropy of a subregion A in the dual CFT is equal to the area of an extremal codimension-two surface γ ext A anchored at the boundary of AdS and homologous to the subregion A [5] S A = min This is known as the RT/HRT formula and the extremal surface γ ext A is often referred to as the RT surface. However, this formula gives only the leading O(1/G) contribution to JHEP08(2020)152 the entanglement entropy, and ignores corrections coming from quantum fields in the bulk. To obtain a more accurate answer one must naturally include the entanglement entropy of these fields, an idea first advocated in [6] and now called the FLM formula Here Σ A is the codimension-one region in the bulk between the RT surface and the boundary region A in a time-slice that contains the two, and S bulk is the von Neumann entropy of the semiclassical bulk density matrix of quantum fields in Σ A . There is a further proposal [7] that extends the FLM formula and is analogous to the island formula in equation (1.1) [8,9]. According to this proposal, the entanglement entropy is given in terms of a different extremal surfaceγ ext A called quantum extremal surface so that This new prescription agrees with FLM at order O(1), but generalizes it to all orders in G. However, the implementation of (1.4) would require understanding of the bulk entanglement entropy for general regions, a task that is very difficult to accomplish in practice, in higher than two bulk dimensions. These proposals highlight the utility of holography in the study of non-perturbative aspects of the black hole information paradox. In this paper we will use this framework to start a systematic study of entanglement entropy in cases where the bulk state is fully time dependent and completely under control, and interpret the results in terms of the microscopic CFT description. As a first step, we will focus on the leading O(1) corrections to the RT/HRT prescription which are given by the FLM formula (1.3). In particular, our goal is to discuss 1. Evidence that the FLM formula works in dynamical settings. 3. Properties of bulk entanglement entropy S bulk .
We will do so in a simple model of a quantum field theory interacting with dynamical gravity in AdS 3 . We will also check all our calculations in the dual CFT 2 and explain how the different terms in the left-hand side and right-hand side of equation (1.3) match. As a byproduct, we will show that we can interpret our setup as a simple toy model for black hole evaporation, thus, providing fresh insight on the black hole information problem.

Road map and summary
Dynamical gravity backgrounds can be obtained by starting with a time-independent geometry and then perturbing it globally or locally. These backgrounds are holographically dual to time-dependent states in the dual CFT where these perturbations are also known as quantum quenches. The two categories of global or local perturbations are dubbed JHEP08(2020)152 global and local quenches, respectively. Both types of quenches have been widely studied in the context of holography [10][11][12][13][14][15][16][17][18][19]. They have been used to study the problem of thermalization of closed quantum systems [20][21][22][23][24][25][26] and to describe the process of black hole formation [27][28][29][30][31][32][33][34], an important problem in gravity. In this paper, we will however focus only on local quantum quenches since they provide more realistic models of perturbations as well as observers.
Local quenches have been used to describe the spread of local perturbations in QFTs and many-body systems [35,36]. For holographic field theories, the holographic description of local quenches is useful to study quantum chaos [37,38]. They also serve as models for generation of quantum entanglement [39][40][41], a resource for quantum computations. In fact, local quenches are approximately equivalent to a quantum gate [42] and sequence of quantum gates are ubiquitous in quantum computations [43]. The study of local quenches is exciting also because there is hope that local quenches could be simulated in condensedmatter systems like cold atoms [44][45][46]. From a holographic perspective, local quantum quenches provide models for the dynamics of localized perturbations in the dual gravity theory [15]. One can use these states to study the real time dynamics of bulk fields coupled to gravity. For example, if the unperturbed background is that of a black hole, one can compute quantities such as quasi-normal modes and other time-dependent bulk observables of interest. This is the setup that we will mostly work with. In particular, we will be interested in understanding entanglement structure of matter fields interacting with dynamical gravity with local perturbations.
In section 2, we start by describing local quench states in 2d CFTs. We calculate the time-dependent entanglement entropy for a single interval of length 2R centered at x = x c , using the replica trick. The general formula for the change of entanglement entropy in the excited state is given by equation (2.40), and can be written as a sum of two contributions, defined by equations (2.41) and (2.42), respectively. The final result for the two terms are: given by (2.69) and (2.73), respectively. Here, ∆ is the conformal dimension of the perturbing operator, which is inserted at x = x 0 and spatially smeared over a length-scale of order α. Later in the paper, we set x 0 = 0 which can be done without any loss of generality, given the translation invariance of the vacuum state. The δS uni A term is fixed by kinematics and is a universal, i.e., valid for any CFT. The δS dyn A term depends on dynamical data, in particular, on a specific correlation function that is not fixed by symmetries. For the evaluation of this term we use an OPE expansion and keep only the first non-vanishing contribution, therefore the final result is only valid in the small R limit. Further, we use OPE data for large-c CFTs and include all corrections from single-and multiple-trace operators at order O(c 0 ).

JHEP08(2020)152
In section 3, we discuss the holographic time-dependent geometry dual to the local quench state. We argue that the time-dependent metric is given by the action of a specific large diffeomorphism, given by equations (3.28)- (3.30), on the backreacted metric of a one-particle state in global AdS, given in (3.12). We also compute the expectation value of CFT operators dual to light bulk fields, with particular emphasis on the stress-energy tensor and its time-evolution following the quench. The result of this calculation is given in equation (3.39).
In section 4, we discuss the holographic computation of entanglement entropy at order O(G 0 ) in the excited state. The calculation naturally splits in two contributions, corresponding to the two terms in the FLM formula (1.3). Firstly, there is the correction to the area term due to change in geometry, dubbed geometric correction δS geom A , which amounts to compute (4.5). For intervals centered at x c = x 0 = 0, we show in equation (4.20) that it is equal to with a(t) given in (4.18). We observe that the first term above matches exactly the universal correction δS uni A obtained from CFT calculation. For non-centered intervals, we show that the geometric correction can be written as a simple integral expression, given in equation (4.32), which can be explicitly integrated for integer values of ∆. For arbitrary values of ∆, we show that it can generically be expressed as a sum of two series expansions, Once again, the first series matches the universal correction δS uni A obtained from CFT calculation. We further show explicitly that this series follows from an application of first law of entanglement in the CFT, equation (4.57), which gives the result at linear order in the change of the density matrix δρ. The appearance of the second series, with terms ∝ R 2∆+2i not seen in the CFT expressions, raises a puzzle which we resolve in subsection 4.2.
In subsection 4.2, we compute the corrections of order O(G 0 ) coming from bulk entanglement entropy δS bulk in the one-particle excited state. We show that the change in the bulk entanglement entropy up to linear order in the bulk density matrix δρ bulk is as shown in equation (4.68). This clearly cancels the second series in δS geom A , thus removing the terms that do not appear in the CFT calculation. We explain that this cancelation is, in fact, expected from the precise relation between the CFT and bulk modular Hamiltonians, given in equation (4.60). Then, we use the bulk replica trick to calculate the JHEP08(2020)152 bulk entanglement entropy at quadratic order in δρ bulk . This calculation involves the Bogoliubov coefficients that bring the global state to the Rindler state associated with the entanglement wedge of A, which we only obtain in the small R limit. The final result, given in equation (4.79), reads This term is in perfect agreement with the dynamical correction δS dyn A obtained in the CFT, thus providing a explicit check of the FLM formula in our dynamical setting.
In subsection 4.2.3, we argue that our setup can be interpreted as a simple toy model for black hole evaporation. In order to reach this conclusion, we use the CHM map [47] and specialize to a Rindler observer adapted to the entanglement wedge of the region A, as shown in figure 11. Indeed, for such an observer the bulk geometry coincides with that of a planar BTZ black hole, with metric given in (4.83), so that δS bulk computes the entropy of matter fields outside the black hole. We further observe that, as seen from the global perspective, the RT surface and hence the black hole horizon generically vary in time. For the particular case of centered intervals, the entanglement wedge corresponds to that of an interval with δτ ≡ |τ 2 − τ 1 | = 0 and opening angle δθ ≡ |θ 2 − θ 1 | given by equation (4.86), Thus the horizon decreases monotonically in time as time goes from t = 0 to t = ∞. In this way, we can interpret δS bulk as the entropy of Hawking radiation in an eternally evaporating black hole geometry. In fact, we show that δS bulk in this case follows the expected behavior for a unitary Page curve, as seen in figure 12, i.e., increasing from t = 0 up to t = t Page and then decreasing as the black hole fully evaporates. We explain that this behavior follows directly from the purity of the quantum state in combination with the Araki-Lieb inequality (4.90). This shows that unitarity can be preserved in a semiclassical analysis, and implies that in more realistic models of black hole evaporation, the information loss or the lack thereof should not be interpreted as an artifact of the approximation.
In section 5, we conclude with some important remarks and a list of open questions.

CFT entanglement entropy
In this section we will discuss the CFT calculation of the entanglement entropy of a single interval after a local quench.

JHEP08(2020)152
2.1 Reduced density matrix for the local quench Regularized quenched state. Consider the vacuum state of a two-dimensional, largec CFT: |0 . We are interested in a special class of excited states produced by locally quenching the CFT at time t = 0 by inserting an operator O(0, x 0 ) at a point x = x 0 . Due to the quench, the state |0 changes to the excited state Exactly localized states of the form (2.1) contain modes of unbounded frequency, which is problematic if we want quantities like energy (density) and entanglement entropy in the state |ψ to be finite. Hence, we need to regularize the state. A convenient way to do this is to give the time coordinate of O a small imaginary part α. We will thus take the state to be |ψ = Notice that in this state α acts as a UV regulator, so high energy modes are effectively suppressed by the factor e −αH . The constant √ N is an appropriate normalization that ensures that the state |ψ as written above has norm one. We will omit this factor from here on. We can recover the proper normalization of various physical quantities whenever necessary.
In the Schrödinger picture, this excited state under consideration evolves under the CFT Hamiltonian H for t > 0, leading to the spread of energy and entanglement. Specifically, at time t > 0 the state can be written as Consequently, the density matrix of the CFT is given by where the operator O(τ, x) is now in the (Euclidean) Heisenberg picture, i.e., and we have defined the Euclidean times τ 1 = −τ + α and τ 2 = −τ − α. This object has a well defined Euclidean path integral representation, namely its matrix elements are given in terms of an Euclidean path integral over R 2 with open cuts on the τ = 0 + and τ = 0 − surfaces respectively. Schematically this is

JHEP08(2020)152
Reduced density matrix and its moments. Let us now consider a subsystem Starting from (2.8), we can arrive to a formula for the reduced density matrix on A, defined as the partial trace In practice, this is implemented by identifying the open cuts at τ = 0 + with the one at τ = 0 − on A c and summing over all field configurations there. The result is a path integral over R 2 with open cuts along A at τ = 0, where the field configurations φ and boundary conditions ψ A and ψ A are defined on A.
The path integral representation is quite convenient since it relates the calculation of the moments of the reduced density matrix trρ n A with a path integral computation over a multi-sheeted Riemann surface Σ A n , obtained by taking n copies of the QFT on R 2 ⊗n and sewing them together along the region A. We denote that path integral on this state by Z n (A), so that the above relation can be written as For the normalized density matrixρ A = ρ A /trρ A one has the formula where Z 1 is the Euclidean path integral evaluated on a single copy of R 2 , with the operator insertions that create the state. Specifically, taking the trace of (2.10), we arrive at It is easy to see that Z 1 can be further related to the two-point function of Euclidean operators, (2.14) Identifying the denominator in (2.14) as the trace of the density operator in the ground state, tr |0 0| = Z gs 1 , we obtain that A similar analysis leads to an expression on the multi-sheeted surface that generalizes (2.15): where the index k represents the fact that the operator is inserted in the k th sheet of Σ A n and Z gs n (A) represents the ground state partition function on the replicated manifold Σ A n , i.e. one without any operator insertions. JHEP08(2020)152

The replica trick
We are now interested in calculating the von Neumann entropy associated with ρ A , (2.17) However, in QFT the Hilbert space is infinite dimensional and it is hard to calculate the logarithm ofρ A . Instead, one often uses the so-called replica trick, where one first computes the Rényi entropies, defined as 18) and then extract the entanglement entropy as the formal limit In practice, most methods allow us to calculate S n for n ∈ Z. To take the limit in (2.19) one analytically continues n to real numbers and defines S n in the neighbourhood of n = 1. Such a continuation exists and is unique provided S n has proper asymptotics. 1 Using the path integral representation of the reduced density operator and its moments, one can write down the following formula for the Rényi entropies where Z 1 and Z n (A) are given by (2.15) and (2.16), respectively. This expression can be massaged into the following form (2.21) Identifying the second term in (2.21) as the ground state Rényi entropy S gs n , one can then write down a simple formula for the regularized Rényi entropy, δS n ≡ S n − S gs n , Provided one can compute this quantity and find its analytic continuation for n ≈ 1, the replica trick gives us the regularized von Neumann entropy δS A , We will now proceed to evaluate this quantity and study it in detail. 1 The exact conditions are stated in what is known as Carlson's theorem. JHEP08(2020)152

Conformal mapping
In order to exploit the full power of conformal invariance in two dimensions, it is convenient to use complex coordinates w andw to label our space-time points, where In these coordinates the full density operator becomes and with the τ i given below (2.7). Expanding it out, these coordinates are explicitly given by Next, we need to express (2.23) in complex coordinates. For that purpose we adopt the following prescription. We label the coordinates on the k th sheet of Σ A n as (w (k) ,w (k) ) and leave the operators O and O † unlabeled. We further define the insertion points of operators on the k th sheet to be With these changes in mind, the entanglement entropy (2.23) becomes (2.28) The uniformization map. In general, computing a 2n-point correlation function on a non-trivial Riemann surface such as Σ A n is a very complicated task. Fortunately, it is simplified in two dimensions due to large conformal symmetry. The group of conformal transformations coincides with the analytic coordinate transformations: (2.29) By using a suitable function f such that f : Σ A n → C, one can map the 2n-point correlation function on Σ A n in (2.28) to a 2n-point correlation function on C, which is simpler by the virtue of the analytic transformations mentioned above. Such a map exists and is known as the uniformization map where w ∈ Σ A n , and z ∈ C. Under its action, the left end-point x L of the open cut on Σ A n is mapped to the origin of the complex plane C, whereas the right end-point x R is mapped to complex infinity along a direction that differs from one sheet to another. In particular, on the k th sheet, x R is mapped to complex infinity along the angle θ k = 2π(k − 1)/n with JHEP08(2020)152 respect to the real axis. This is easy to prove by studying the w (k) → x L and w (k) → x R limits of (2.30). Near the left end-point on the k th sheet, w (k) is obtained by circling counter-clockwise (k − 1) times around x L , at a fixed but infinitesimal distance away from it. That is, w (k) = x L + e 2πi(k−1) . Similarly, near the right end-point x R (k − 1) on the k th sheet, w (k) is given by circling clockwise around x R at a fixed distance , i.e. w (k) = x R − e −2πi(k−1) . In the limit, this gives us where ≡ x R − x L is the length of the interval. Moreover, under this map the k th sheet of Σ A n is mapped to a sector of the complex plane C defined within angles 2π(k − 1)/n ≤ θ ≤ 2πk/n, as shown in figure 1. Finally, it is easy to see that the insertion points w 2k+1 , w 2k+1 , w 2k+2 andw 2k+2 in Σ A n are mapped to z 2k+1 = e i2πk/n z 1 , z 2k+2 = e i2πk/n z 2 , on C, respectively, where z 1 ,z 1 , z 2 , andz 2 are the principal roots in (2.38) JHEP08(2020)152 Figure 1. Schematic representation of the uniformization map (2.30). After this transformation is implemented, each sheet of the original Riemann surface Σ A n maps to a wedge in the complex plane C with 2π(k − 1)/n ≤ θ ≤ 2πk/n. The insertion points are mapped according to (2.33) and (2.34) which means that i) only two of these points are inserted in each wedge and ii) points with k ≥ 1 differ only by a phase to those with k = 0. For this plot we have set n = 12 as an illustration, but for our particular calculation we are interested in smaller values of n, specifically, in the vicinity of n ≈ 1.
It is also convenient to define the factor J 1 With these notations in mind, we can now rewrite the expression for the regularized entanglement entropy (2.28) as 4 which is naturally separated into two contributions: one which we call the universal part, and a dynamical part, defined as (2.42)

JHEP08(2020)152
Notice that the former one depends only on the dimension of the quench operator O and the underlying geometry of Σ A n , but is otherwise independent of the coupling constants in the CFT under consideration. The latter one depends on a higher-point function, which is not fixed by conformal symmetry, and it is therefore not universal.

Regularized entanglement entropy
In the previous subsection we showed that the regularized entanglement entropy after the local quench can be naturally separated into two contributions, with universal and dynamical parts, δS uni A and δS dyn A , given in (2.41) and (2.42), respectively. We will now compute each of these contributions.

Universal contribution
The calculation of the universal part (2.41) can be split in two steps: first, we need to compute the Jacobian factors J n O † , J n O , J 1 O † and J 1 O , and then we need to analytically continue the result for n ∈ R and take the limit n → 1. For the computation of the Jacobians we need the derivatives of z andz, with ≡ x R −x L , evaluated at the points z 2k+1 , z 2k+2 ,z 2k+1 andz 2k+2 respectively. Notice that we have expressed the answers directly in terms of the z-coordinates so the evaluation is now straightforward. To get a closed expression for J n O † J n O notice that the first factor of (2.37) can be expressed as (2.45) where we have used (2.33) and the fact that e 2πik = 1 for k ∈ Z. The second factor of (2.37) is just the complex conjugate of (2.45) so the product of the two cancels their phases. Repeating the same analysis for (2.38) and combining the two results, we obtain We also need the factor J 1 Writing down the explicit form of the two-point correlators in (2.39), it follows that which, in combination with (2.46), leads to a closed expression for δS uni A solely in terms of the insertion points z 1 and z 2 , (2.48) It only remains to carry out the analytic continuation for n ∈ Z and take the explicit limit n → 1. We will merely transcribe the outcome of the calculation and relegate the details of this analysis to appendix A.1. The final result yields: where ∆ = 2h is the dimension of the quench operator O. Here, the insertion points z 1 , z 1 , z 2 andz 2 are given by setting n = 1 in the general expressions (2.34), i.e., and similarly for z 2 andz 2 . Finally, we notice that by using the identity we can rewrite (2.49) in a slightly more useful form: which is conveniently written in terms of (z 1 − z 2 )/(z 1 + z 2 ) and its complex conjugate. It also has the advantage that its small interval expansion converges much faster than the one obtained from (2.49), a property that will prove useful later.

Dynamical contribution
Next, let us consider the dynamical contribution (2.42). To calculate the 2n-point correlator in the numerator we will work in the special class of holographic CFTs. We will further simplify the calculation by taking the limit of small intervals. First, notice that in the small interval limit x R → x L or → 0, all the insertion points localize around the unit circle |z| = 1, as shown in figure 2. This behavior can be traced back to the uniformization map (2.30), and can be seen more directly from the expressions (2.33) and (2.34). In this limit we obtain where we have used x c and R to label the center and the half-length of the interval, respectively Similarly, it is easy to see that in this limit |z 2k+1 − z 2k+2 | ∼ O(R) for all other k since the insertion points with k ≥ 1 differ only by a phase to those with k = 0. This implies that all pairs of operators with equal k in the 2n-point correlator in (2.42) become arbitrarily close in the limit R → 0. See figure 2 for an illustrative representation of this limit. Assuming JHEP08(2020)152 Figure 2. Small interval limit of the 2n-correlator in the complex plane C. Assuming R → 0 all pairs of insertion points with equal k (i.e., in the same sheet of the original Riemann surface Σ A n ) approach to each other and one can carry out an OPE expansion between the two. This is possible provided that n is not too large, so points with different k are still a finite distance away from each other as R → 0. For this plot we have set n = 12 as an illustration, but for our particular calculation we are interested in smaller values of n, specifically, in the vicinity of n ≈ 1.
that n is not too large, points with different k will be a finite distance apart in this limit, and we can replace each product by an OPE expansion of the form (2.55) In such an expansion, the sum runs over all primary operators Φ while the dots represent contributions from their descendents. Notice, however, that we have isolated the contribution from the identity operator Φ = 1 1 (the "1" outside the sum), whose coefficient is fixed by the normalization of the two-point correlator. Finally, the operators Φ inside the sum must be evaluated at an arbitrary point z c k in C within the radius of convergence of the expansion, which we can choose to place the fused operators. A simple and convenient choice would be to pick the "center point" (rotated to the appropriate sector corresponding to a given sheet), i.e., At this point, and in order to proceed with the calculation, we need to specify more data about the CFT of interest. Since we want to compare our results with a bulk calculation using the AdS/CFT correspondence, we will focus on CFTs with holographic duals or, in other words, holographic CFTs. Generically, these are theories with large central charge c and a sparse spectrum of low-dimension operators [48][49][50][51]. One of the important features that follow from these properties is large-c factorization [52]. Large-c factorization is the statement that for single-trace operators Φ the OPE coefficients C OOΦ are all JHEP08(2020)152 suppressed, i.e., Since we are interested in the contribution to the entanglement entropy at order O(1), we can safely ignore these contributions. On the other hand, multi-trace operators can indeed have OPE coefficients at O(1) so they must be considered in the sum of (2.55). In particular, in the limit of small R, the leading contribution coming from the sum in (2.55) is given by the lightest multi-trace operator that can appear in the OPE, while all others are suppressed by a higher power of R. This operator is Φ = :O 2 :, a double-trace, and has conformal dimension ∆ Φ = 2∆. In holographic CFTs, its OPE coefficient is given by [53] Next, we insert the OPE expansions in equation (2.55) for all pairs of operators with the same index k in equation (2.42). At leading order, the result is given by all pair of operators replaced by the identity contribution. This gives rise to (2.59) If we plug this into (2.42) we obtain exactly zero, but it is easy to understand why. The reason is that we have already factored out the Jacobians of the correlators in the numerator and denominator of (2.23), in what we have called the universal contribution (2.41). This means that we can in fact interpret the universal term as the contribution coming from the identity operator. Similarly, we can analyze the sub-leading corrections by imagining the case when only a few of the pairs are replaced by the operator O 2 , while the remaining ones are replaced by the identity. Now if only one pair is replaced by O 2 , this sub-leading correction to (2.42) vanishes. This is because all one-point functions of local operators are exactly zero in the vacuum due to conformal invariance. Therefore, the first non-trivial contribution to the sub-leading term appears when two pairs are replaced by the operator O 2 . Considering all possible Wick contractions, and using (2.58), this yields , with z c k andz c k given in (2.56). Using the replica symmetry one can fix the location of one of the operators appearing in the double sum of (2.60), say to the location of the operator with = 0, and multiply the result by n. This leads to:

JHEP08(2020)152
where the extra factor of 1/2 is introduced to avoid double counting. The sum in (2.61) can be carried out for n ≈ 1 as was done in [54] and further generalized to generic thermal green functions in [55]. We will merely write the answer here, and relegate the details of this analysis to appendix A.2. At the end, the outcome of the calculation yields: . (2.62) Combining (2.60)-(2.62) leads to an expression for the leading term in the n → 1 expansion of the logarithm in (2.42), i.e., .
(2.63) Finally, taking the n → 1 limit in (2.42) we obtain We emphasize that this is only the first term in the small R expansion of δS dyn A . On the other hand, our result for δS uni A given by equation (2.52), is valid for any R.

Analytic continuation and real-time interpretation
Given our final expressions for δS uni A and δS dyn A given in (2.52) and (2.64), the final step is to analytically continue these results to real time. We will do so in this subsection. Along the way, we will uncover a clear physical picture for the spread of entanglement after local quenches that will allow us interpret our results in a transparent way.
Universal contribution. Consider our result for δS uni A given by equation (2.52). As explained before, this is conveniently written in terms of (z 1 − z 2 )/(z 1 + z 2 ) and its complex conjugate so we will start by writing this combination in terms of the parameters {τ, x 0 , x c , R, α} using (2.26) and (2.50): Performing the analytic continuation to Lorentzian time, τ → it, leads to and Notice that after the analytic continuation, both parameters become purely imaginary and are no longer the complex conjugate of each other. Defining the parameters Figure 3. Pictorial representation of the universal contribution to the spread of entanglement entropy after a local quantum quench in 2d CFTs. The state at t = 0 is the vacuum state perturbed by a local operator smeared over a region of compact support ∼ α. The state evolves under the CFT Hamiltonian H for t > 0, generating an entangled pair of wave packets that move in opposite directions at the speed of light. The wave packets eventually increase the entanglement entropy of , and then decrease it as they disperse to infinity. and using the identity arctanh(ix) = i arctan(x), we can rewrite δS uni A as Equation (2.69) allows an interesting interpretation. The entanglement generated by the perturbation can be interpreted as coming from two independent sources (corresponding to the two independent components η ± moving in opposite directions and at the speed of light. Each of these sources contributes by a term −(∆/η) arctan η with the specific η ± depending on the motion of the source, see figure 3 for a pictorial representation.
Physically, the two sources appear because the operator inserted at (t = 0, x = x 0 ) creates an entangled pair of wave packets moving in opposite directions (as required by momentum conservation) each with a characteristic size of order α. If we further assume that quantum correlations between the two are negligible compared to the quantum correlations of each component with itself whenever the distance between the two is much larger than α, we can predict that the entanglement entropy of a region A will be the sum of two contributions, consistent with the separation found in (2.69).
Furthermore, this physical picture predicts that, as a function of time, entanglement entropy is maximal when the location of the two endpoints of ∂A, {x L , x R } symmetrically divide each of the traveling wave packets in two. We can calculate the time when this happens.
• Operator inserted inside the interval. In this case we expect that the centers of both right and left moving components η + , and η − would reach the endpoints of the entangling region, x R = x c + R and x L = x c − R, at different times t ± > 0. In order to determine these times we need to equate x(t) = x 0 ± t = x c ± R and solve for t. This gives rise to the times t ± = R ± (x c − x 0 ).

JHEP08(2020)152
• Operator inserted outside the interval. Let us assume for concreteness that the operator has been inserted to the left of the interval. In this case, the left moving wave packet will not cross the interval and so its contribution will not reach a maximum. On the other hand the right moving wave packet will reach local maximum at two different times t ± given by the solution of Both situations are consistent with our assumptions provided that α is sufficiently small. If the above picture is correct then one should be able to reproduce the local times derived above in both scenarios using the formula (2.69). To do so, notice that the function −(∆/η) arctan η is an even function of η and a monotonically increasing function of |η|, therefore its local maximum values occur when |η| is globally maximal. From (2.68), and for R − |x 0 − x c | α, it is clear that this happens at the times Since our setup only considers positive times, then depending on whether the operator is inserted inside or outside the interval, we will have the following values for t ± : • Operator inserted inside the interval. In this case which means that out of the four times given in (2.70) the only positive ones are These results are in perfect agreement with the expectations from the interpretation in terms of entangled wave packets traveling at the speed of light.
which, after the analytic continuation to Lorentzian time, τ → it, becomes This expression does not have a similar separation as the one found for the universal piece δS uni A . This fact is manifest if we rewrite (2.72) as (2.73)

JHEP08(2020)152
Indeed, this product factorization shows that the full quantum state has non-vanishing entanglement between the two wave packets.
Let us briefly discuss δS dyn A as a function of time. It is maximal when the denominator in (2.72) is minimal. This function cannot be arbitrarily small as it comes from a sub-leading term in the OPE expansion (2.60). The smallest it can be is zero, which happens at the roots of (2.74) which are These roots have real, positive values for R > α. Hence, consistency with the OPE expansion requires R < α. Further, local extrema of (2.74) occur when For t ≥ 0, it is easy to see that the global minima occur at In the next sections, we will recover terms of this form from the bulk calculation and show that they generically arise from contributions due to bulk entanglement entropy.

Aspects of the gravity duals
We will now describe the geometry of quench states considered in the previous section via holography. We will follow the method described in [15], wherein it was developed for operators of large conformal dimension ∆ c. It involves finding the backreacted metric for a localized perturbation which is initially localized near the boundary and falls into the interior of AdS. In the following, we will generalize this calculation to a one-particle excited state of a light scalar field coupled to gravity. The lightness of the scalar field implies that the dual operator have small conformal dimension, ∆ c.

Light operator excited states and bulk backreaction
We begin by discussing the backreaction of a quantum scalar field on the metric of a pure AdS 3 spacetime. We start with the action where G is the Newton's constant in 3 dimensions, R denotes the Ricci scalar and Λ is the cosmological constant which is fixed in terms of the AdS radius Λ = − 1 L 2 . This theory is dual to a CFT 2 with central charge

JHEP08(2020)152
which is large provided that G/L is small. Moreover, the mass of the scalar field m is related to the conformal dimension of the dual operator O ∆ , through Since we are interested in the case of light operators, ∆ c, we require that mG 1. To find the backreacted metric for the state dual to the local quench, following [15], we start in global coordinates with the vacuum AdS 3 solution The scalar field can be expanded in terms of modes on this background, each mode labeled by two quantum numbers, corresponding to an expansion in either of the two space coordinates (r, θ). The lowest energy mode is an S-wave, as described in [56]. The wavefunction for this mode is spherically symmetric and is annihilated by the isometries L 1 andL 1 of AdS 3 , i.e., L 1 |ψ =L 1 |ψ = 0. This mode defines a one-particle excited state for the scalar field on the pure AdS 3 background where a † 0,0 denotes the creation operator. Solving for the wavefunction of this mode one can show that The normal-ordered stress-energy tensor of the scalar field is whose one-point function ψ|T µν |ψ can be evaluated in the one-particle excited state as defined above. Using spherical symmetry, the off-diagonal components of this one-point function can be shown to vanish. The diagonal components can be evaluated using (3.6) and they are given by [57] ψ|T τ τ |ψ = ∆(∆ − 1) This stress-energy tensor backreacts on the AdS 3 vacuum and the backreacted geometry can be obtained by solving semi-classical Einstein's equations

JHEP08(2020)152
where R µν denotes the Ricci tensor. Since the source ψ|T µν |ψ is diagonal, one immediately concludes that the backreacted metric has to be diagonal. So one can propose the following ansatz The functions F i (r) were first determined in [57], and at order O(G∆/L) are given by (3.11) Hence, the backreacted metric at this order is (3.12) The stress-energy tensor, although smooth, becomes sharply peaked around the origin as one increases ∆. In the limit of heavy operators ∆ c or, equivalently, when the mass of the scalar field is large mG 1, the backreacted background can be shown to approximate to that of a conical defect. This was indeed the case considered in [15]. However, we will work in the opposite regime of light operators and a non-trivial wavefunction for the scalar field.

Local quenches via large diffeomorphisms
As discussed in the previous section, our local quenches are defined on a plane R × R. We will take this to be the boundary of Poincaré AdS 3 . The local quench can then be described by the motion of a localized perturbation, which is localized near the boundary at time t = 0 and then falls into the interior of AdS. We will use the backreacted global metric in equation (3.12) to compute the perturbed Poincare geometry. As discussed in [15], the idea is to perform a large diffeomorphism on the circle R × S to obtain the plane R × R. The diffeomorphism extends naturally inside the AdS spacetime.
The transformation has two parts. First let us recall that AdS can be thought of as an hyperboloid embedded in a higher dimensional Minkowski spacetime with two timelike coordinates. The definition of the global patch and Poincaré patch of AdS 3 in terms of the R 2 × R 2 coordinates is the following, and they satisfy the constraint Moving between the two patches, one finds that a stationary point at r = 0 in the global patch maps into a non-trivial geodesic in the Poincaré patch, with x = 0 and z 2 − t 2 = L 2 . This is problematic, since the particle cannot reach arbitrarily close to the boundary. An easy way to fix this is to consider boost along the T and X directions. This leaves the pure AdS space invariant, but modifies the geodesic to where β is a boost parameter. This trajectory captures the desired behavior in the Poincaré patch: for small enough α it gets arbitrarily close to the boundary at t = 0. The backreaction of the perturbation following this geodesic gives us our local quench. For finite α, the geodesic does not reach the boundary. This is equivalent to preparing the state at t = 0 by smearing the operator over a region with finite support of order α, consistent with the standard notion of UV/IR connection [58][59][60]. See figure 4 for graphical representation.
As an intermediate step, and for future reference, we give here the explicit transformation between the original global frame (τ, r, θ) and a boosted global frame (τ , r , θ ): In the left figure we have plotted the original global frame, represented with coordinates (τ, r, θ), where the particle lies at the origin r = 0. Specializing to a Poincaré patch maps this geodesic to z 2 = L 2 + t 2 , so the minimum approach to the boundary is z = L at t = 0. This surface is depicted in brown. In the right figure we have plotted the boosted global frame (τ , r , θ ) where the particle oscillates according to (3.23).
Specializing to a Poincaré patch of this boosted frame maps the particle's trajectory to z 2 = α 2 +t 2 , with α ≡ L e β , so it can now get arbitrarily close to the boundary. The minimum approach is now z = α at t = 0. This surface is also depicted in brown in the corresponding cylinder. In this boosted frame, the r = 0 geodesic maps to which is periodic in τ . The geometry in this frame is dual to a CFT state with perpetual collective oscillations, of the kind studied in [61]. The final transformation consists of specializing to a Poincaré patch of the boosted global frame. The form of this transformation is the same as given in (3.13)-(3.16) but with (τ, r, θ) replaced by (τ , r , θ ). In figure 5 we show pictorially the effects of these transformations. We can also obtain the explicit combined transformation from the original global frame (τ, r, θ) to the latter Poincaré patch. This map is given by:

25)
JHEP08(2020)152 Equivalently, inverting these relations, and writing them explicitly in terms of α we get: The full metric after the combined coordinate transformation is straightforward to obtain but is very lengthy and not particularly illuminating. Hence, we will not transcribe it here. We want to calculate all O(G 0 ) corrections to the holographic entanglement entropy due to the presence of scalar field in the bulk. As we will see in section 4, they include a geometric correction, due to the backreaction of the quantum fields on the geometry and a quantum correction due to the bulk entanglement entropy of the quantum fields. Altogether, these two corrections will give all leading O(c 0 ) corrections to the local quench in the CFT. One important thing to note is that the bulk entanglement entropy piece is sensitive to the scalar field profile in the bulk. For finite ∆, this is not a problem, since φ(t, x, z) is smooth everywhere. However, as ∆ increases the profile becomes sharply peaked around the origin, rendering this problem largely degenerate. In fact, the requirement ∆ c is precisely what will enable us to carry out such a calculation without problem.

One-point function of local operators
Before proceeding to the computation of the entanglement entropy, we study here other observables of interest, i.e., the one-point function of local operators dual to light bulk fields. In our theory (3.1) we have two of these fields: the metric g µν and the scalar field φ. In Fefferman-Graham (FG) coordinates, we can write the following near-boundary expansions: In terms of these expansions, and given the normalization of the action (3.1), the one-point function of the CFT stress-energy tensor and scalar operator are given by [62]: However, a brief comment is in order. Since we are working in the semiclassical approximation, φ must be treated as an operator. Hence, in order to use the above formulas, JHEP08(2020)152 we first need to find the expectation value of φ in the quantum state |ψ = a † 0,0 |0 . This amounts to compute a 3-point function in the bulk, which should vanish for a free theory. Indeed, a quick calculation shows that Therefore, it immediately follows that Next, for the computation of the stress-energy tensor of the boundary CFT we must write the bulk metric in the Fefferman-Graham gauge (3.31). In order to do so we need a second coordinate transformation, which can be obtained perturbatively as (3.37) The function f (t, x) which satisfies the conditions of the FG gauge is found to be With this change of coordinates, we can make use of (3.33) to obtain the stress-energy tensor in the CFT: From these expressions one can obtain quantities of interest, such as the energy density E = T tt and momentum density P = T tx , which by symmetry equals the energy flux T xt . The pressure in this case also equals the energy density, since T xx = T tt , as expected for a conformal theory in 2-dimensions. We now make a couple of comments. First notice that both the traceless condition and stress-energy conservation are satisfied, Second, the total energy is constant, and in agreement with the expectation for a oneparticle state, for an insertion of a primary of dimension ∆ The extra term e −β here accounts for the boost factor. The total momentum vanishes, p = 0, because the excitations generated by the quench move both to the left and to the right. In order to better understand this time dependence, we plot in figure 6 the two non-trivial components of the stress-energy tensor, namely the energy density E and the momentum density P. As expected, we observe profiles that are peaked on the light-cone −t 2 + x 2 = 0, which can be understood as shock waves that move at the speed of light [63] due to the initial excitation at t = x = 0. JHEP08(2020)152 Figure 6. Profiles of the normalized energy densityẼ = E/E and momentum densityP = P/E, with E = ∆/α, as a function of t and x for a local quench with α = 1/2. The parameter α measures the initial spread of the excitations around x = 0, which is estimated to be of order O(α). The initial perturbation generates two shocks, moving to the left and to the right respectively, that move at the speed of light. The heights of these shocks remain constant due to energy conservation in the CFT.
In the CFT, these shock waves arise from the action of a primary operator O ∆ on the vacuum [64][65][66][67], as discussed near equation (2.2). In fact, we have deliberately identified the parameter α in the bulk trajectory (3.18) with the UV regulator in the quenched state (2.2). For finite values of α, the parameter gives the smearing of the operator O ∆ around x = 0. Only in the limit α → 0, the bulk trajectory touches the boundary and we get an exactly local quench. In this limits, the energy (3.41) post-quench diverges, and so will energy density as seen in figure 6. Analytically, the density in this limit looks like The equal constant E/2 upfront implies that the heights of the two peaks remain constant, a consequence of energy conservation in the CFT. 5 Our results for the stress-energy tensor are in agreement with the results of [15], even though they worked in the regime of heavy operators. This is because the value of the conformal dimension ∆ only appears as a multiplicative factor in the stress-energy tensor. In the bulk, this is manifested in the fact that the backreacted metric of the one-particle state looks identical to that of a conical defect for an observer at infinity [57].

Holographic entanglement entropy
On the gravity side, order 1/N corrections to entanglement entropy can be computed using the FLM prescription [6,7], which states that (4.1)

JHEP08(2020)152
In this formula, γ A is a codimension-2 bulk surface anchored at the boundary, with ∂γ A = ∂A, and S bulk is the entanglement entropy of bulk fields across the extremal surface γ ext A , in a Cauchy slice Σ ⊃ A ∪ γ A . For perturbative excited states over the vacuum, Σ can be taken to be a constant-t slice (or a boosted version of it, for intervals in generic time-slices).
The leading order term in (4.1) can be calculated from the standard RT/HRT prescription [5,68], using only the extremal area in the unperturbed geometry. At order O(1) there are two contributions: one due to the correction to the area term in the backreacted geometry, and another due to the entanglement entropy of bulk fields in the given quantum state. The latter contribution can be computed in the unperturbed geometry because the bulk entanglement in the perturbed geometry would be further suppressed in 1/N and would only appear at higher order in the expansion. The dots in the above equation represent such higher order contributions and could be computed in the framework of quantum extremal surfaces [7].
In this section, we will compute these O(1) contributions to entanglement entropy in the perturbed geometry constructed in section 3. We will study them separately in sections 4.1 and 4.2 respectively.

Geometric corrections to entanglement entropy
Here we will discuss some O(1) corrections to the holographic entanglement entropy of a single interval in an excited state following a local quench. These corrections come purely because the bulk metric that corresponds to the excited state is different from the one that corresponds to the pure state in the CFT. Hence we denote these corrections as "geometric corrections". We relegate the discussion of the O(1) corrections due to bulk entanglement entropy to section 4.2.
Consider an arbitrary perturbation over pure AdS due to matter fields, such that the metric takes the form 6 g µν = g (0) The corrections to entanglement entropy due to the change in the geometry can be computed from the area term in (4.1). The change in geometry (4.2) induces a linear variation in area as follows Area(γ A ) = Area (0) (γ A ) + δArea(γ A ) . The leading order correction in the metric is of order O(G) so it is clear that the δArea term contributes at order O(1) to the entanglement entropy. The leading term and first order variation of the area are found to be and is the induced metric on the extremal surface, with and ξ α are coordinates parametrizing the surface. An important point here is that, at this order in the perturbation, the embedding functions X µ (ξ) can be taken to be the same as in empty AdS. This means that we do not need to know the precise shape of γ A in the perturbed geometry to evaluate (4.5). This simple but useful observation can be nicely illustrated by making use of the variational principle [69,70], and holds true regardless of the expansion parameter.
In the presence of a 1-particle excited state of a light scalar field (3.5) in AdS 3 , the backreacted geometry in global coordinates is given by (3.12). 7 We can expand this metric as in (4.2), to obtain Next, we need to implement the bulk diffeomorphism (3.28)-(3.30) to arrive to a Poincaré wedge where (3.12) is viewed as a local quench. The resulting metric is lengthy so we will not transcribe it here. Next, we would like to compute (4.4) and (4.5) in this new frame. However, since areas are invariant under coordinate transformations, we can work directly in the original global coordinate system but taking care of properly transforming the embedding functions X µ (ξ). In the Poincaré patch of AdS 3 , the metric at order zero is given by In these coordinates, the extremal surfaces that we are interested in are given by (4.10) These are semicircles at a constant-t slice, with radius R and centered at x = x c . The two endpoints of these geodesics are x ± = x c ± R, so they naturally span boundary intervals of length = 2R . Importantly, note that the local quench considered here is due to the insertion of a primary operator at x = 0, whereas the interval of the entangling region is centered at x = x c . The 7 When the scalar field is heavy, the backreaction is instead given by a conical defect geometry. The geometric corrections to entanglement entropy in this case were calculated in [15]. We will focus here on the case where the scalar field has a small mass, or equivalently, is dual to a light operator in the CFT, i.e. ∆ c.

JHEP08(2020)152
special case x c = 0 naturally has more symmetry than the generic case with x c = 0. In the following we will specialize to the more symmetric case with x c = 0, since the calculations will be simpler, but at a later stage we will study the most general case.
Centered intervals. For intervals centered at the origin (x c = 0) we choose to parametrize the geodesic in terms of ξ = z ∈ (0, R), i.e., with X µ = {t(z), x(z), z}, where Now we use the transformations (3.28)-(3.30) to obtain the embeddings in the global coordinate system. In terms of the parameter ξ = z, X µ = {τ (z), r(z), θ(z)}, where and x(z) is given in (4.12). We can compute the leading term and first variation of entanglement entropy by evaluating the embedding functions (4.13)-(4.15) in (4.4) and (4.5), respectively, and using the RT/HRT formula. Since we are considering the x c = 0 case, the symmetry of the problem allows us to take one branch of x(z), say the positive one, and multiply the resulting integrals by a factor of two. After some manipulations we arrive at and where we have defined the function The first integral gives rise to the standard result for the entanglement entropy in the vacuum of a 2D CFT, The second integral is a bit more involved, but can be explicitly performed to obtain (4.20)

JHEP08(2020)152
This expression can be massaged into a more familiar form, by defining a function ϑ(t) according to R a(t) = sin ϑ(t) 2 , (4. 21) or equivalently, With this definition, and using the following hypergeometric identity we arrive at (4.24) Upon identifying ϑ(t) ↔ θ, this formula coincides exactly with the result of [57] for the geometric corrections to the entanglement entropy of a 1-particle excited state of a light scalar field in global AdS. Notice that this is expected, since the length of the geodesic is invariant under general coordinate transformations. Indeed, in the Poincaré coordinates, a constant time-slice interval with endpoints at x = ± /2 is mapped through (3.28)-(3.30) to an interval in the global coordinates at a constant-τ slice with opening angle δθ = ϑ(t) given by (4.22). Reference [57] considered an interval with opening angle δθ, so the agreement of the two results is not suprising. On the other hand, non-centered intervals in Poincaré coordinates (with x c = 0) map to intervals in global coordinates that are tilted in the time direction (with δτ = 0). These were not considered originally in [57]. Before moving to the most general case of non-centered intervals, let us briefly analyze our final result of geometric correction for centered intervals. First, notice that when ∆ is an integer (4.24) takes a much simpler form. For example, for a marginal operator dual to a massless scalar field (∆ = 2), we have (4.25) For other integer values of ∆, the hypergeometric function simplifies to a rational function as well, but the final expressions are longer as we increase the value of ∆. Second, note that the first part of (4.24) matches exactly with the universal geometric term (2.69) obtained in section 2. It turns out that this piece can be easily extracted from the modular Hamiltonian in the CFT, as we will show in section 4.1.1.
Lastly, for future reference, we discuss the explicit expansion of the final result (4.24) for small intervals. In this limit we get two distinct contributions The first series arises from expanding the first part of (4.24) while the second one from the hypergeometric function. The second series might seem puzzling at this point. Recall that the full CFT result at order O(1), obtained in section 2, does not contain terms proportional to R 2∆+2i with i ∈ N. Later in section 4.2 we will see that the second series precisely cancels out with the O(1) corrections coming from the bulk modular Hamiltonian. This remarkable matching is in fact expected and follows directly from the exact relation between the CFT and bulk modular Hamiltonians. We will come back to this point in section 4.2.
In figure 7 we show various plots of the geometric correction (4.24) for centered intervals in different physical scenarios. In the left panel, we plot entanglement entropy for several values of ∆, in a situation where α/ 1, so the quench is localized and the perturbation is sharply peaked. The entanglement itself peaks at t = R = /2, which shows that the excitations created by the quench carry entanglement at the speed of light v E = 1. This behavior resembles the seminal result for evolution of entanglement entropy after global quenches in 2d CFTs [71], which can be interpreted in terms of EPR pairs. In holographic CFTs, the entanglement still propagates at the speed of light, however, the quasi-particle picture is replaced by the entanglement tsunami interpretation [72][73][74][75]. The difference between global quenches and our setup is that the energy injected in a local quench eventually disperses to infinity, so the entanglement entropy drops to zero as the state relaxes to the vacuum state. This sharp transition has a clear interpretation from the bulk perspective: it happens exactly at the moment that the bulk particle crosses the entangling surface γ A . The entanglement entropy peaks at this moment because the particle has largest relative backreaction on the area of the entangling surface. At late times though, the particle falls into the deep infrared, eventually escaping the Poincaré patch at t → ∞.
In the right panel we show the opposite situation, where α/ 1. In this case the quench is smoothed over a region much bigger than the entangling interval and the bulk particle is always outside and far from the entanglement wedge. The evolution of entanglement entropy in this case is monotonically decreasing in time, since the particle is JHEP08(2020)152 always moving away from the surface γ A , hence reducing its effect over time. Remarkably, its full behavior strongly resembles the late-time evolution of entanglement entropy in an expanding boost-invariant plasma, studied in [76,77], specially for large ∆. This can be explained by the fact that, at leading order in hydrodynamics, the bulk dual to the latter system is given precisely by a boosted black brane that moves into the radial direction z. Indeed, in the limit /α 1, the backreacted metric near γ A (close to the boundary) enjoys also an approximate planar symmetry. As a consequence of these two facts, the entanglement entropy obeys the 'first law'-like relation δS A = δE A /T A [78,79] even though the state is time-dependent [70]. Moreover, such a relation can be easily derived from the first-law of entanglement δS A = δ H A provided that the time-evolution is sufficiently slow and the region is small. We will discuss this point more in detail in section 4.1.1.
Non-centered intervals. For generic intervals with center x c = 0, we parametrize the RT surface using a variable ξ ∈ (0, 1) such that As one varies ξ from 0 to 1, x(ξ) varies linearly between x − = x c − R and x + = x c + R.
We again use the transformations (3.28)-(3.30) to obtain the embeddings in the global coordinates, in terms of this parameter. It is X µ = {τ (ξ), r(ξ), θ(ξ)} and takes the form (4.28) 30) and x(ξ) and z(ξ) are given in (4.27). Next, we compute the leading term and first variation of entanglement entropy by evaluating the embedding functions (4.28)-(4.30) in (4.4) and (4.5), respectively, and using the RT/HRT formula. The leading order term yields the same result as for the centered interval (4.19), i.e., the value of entanglement entropy in the vacuum, as expected 8 (4.31) The first variation of entanglement entropy yields an integral of the form

JHEP08(2020)152
where we have defined (4.33) The full integrand as a function of ξ is lengthy so we will not transcribe it here. The next step is to perform the integral. This can be done analytically for integer values of ∆ but the final expressions are cumbersome and not very enlightening. To give a flavor for the kind of expressions one obtains, the explicit result for an operator dual to a massless scalar field (∆ = 2) is given by where F (t, R, x c , α) and G(t, R, x c , α) are rational functions given by respectively. This expression looks singular in the limit of centered intervals x c → 0. However, although each individual piece has terms that diverge in this limit, they cancel out amongst each other upon resummation. Therefore, expanding for small 9 x c we obtain where δS geom A | xc=0, ∆=2 is given in (4.25). Similarly, the geometric correction for the noncentred intervals in the limits t → 0 and R → 0 is given by respectively, which can be easily verified to be well-behaved.

JHEP08(2020)152
Now for arbitrary ∆, the geometric correction for the non-centered intervals can be obtained by expanding the integrand of (4.32) in powers of R and performing the individual integrals. This gives which agrees with leading term of the ∆ = 2 result (4.39) and generalizes the expansion found for centered intervals (4.26) to non-centered intervals. As we will show later in section 4.1.1, the complete series in the first term above agrees with the result obtained from the first law of entanglement entropy in the CFT, δS A = δ H A . Further, the first term of that series captures the 'first law'-like relation δS A = δE A /T A , valid only for small intervals. Lastly, resumming the series in the first term gives rise to the universal term (2.69) obtained from the CFT calculation in section 2.
The series in the second term in (4.40) turns out to cancel out with the result coming from the first law of entanglement entropy in the bulk, δS A = δ H bulk , as will be shown in section 4.2. This explains the absence of terms proportional to R 2∆+2i (i ∈ N), in the final result for the entanglement entropy at order O(1), as is the case for the CFT calculation of section 2.
Finally, it is worth noticing that the two series in (4.40) can be identified with specific contributions coming from the integral in (4.32). More specifically, we can split the numerator of (4.32) such that 42) and D(ξ) is given in (4.33). From these definitions, together with the explicit form of the embedding functions (4.28)-(4.30), we see that these terms have expansions denomiator of N Q . This non-trivial dependence is also manifested in the small interval approximation which gives terms like R 2∆+2i . As we will see, the quantities δS P A and δS Q A will prove useful in sections 4.1.1 and 4.2.1 respectively.
Lastly, we study the time evolution of the geometric correction for the non-centered intervals. Figure 8 shows plots for different physical scenarios. In the left panel, we consider the situation α/ 1, so the quench is localized and the perturbation is sharply peaked around x = 0. We see two distinct behaviors: for 0 ≤ x c / 1 2 the peak of the entanglement entropy is located at t ∼ R, but its value increases with x c / . In this situation the excitations are created within A, and they take on average t ∼ R to exit region. However, as x c / is increased, the perturbations are created closer and closer to the boundary of A, contributing more to the short-range entanglement across ∂A. For 1 2 x c / , the peak happens at t ∼ x c while its value slightly increases as x c /R is increased. The quench in this case is localized around x = 0, which happens to be outside A, and the excitations take on average t ∼ x c to localize within the entangling region. We can attribute the latter behavior to the fact that α/ is in fact finite, so for x c / just above 1 2 , a small portion of the quench profile is still supported within A. All these results are consistent with entanglement propagation at the speed of light v E = 1, be it via quasi-particles [71] or wave packets described in section 2.4.
In the right panel, we show plots for the case α/ 1, so the quench is supported over a region α much bigger than the entangling interval. From the bulk perspective, the behavior follows because the particle is always falling in the deep IR, far from the entanglement wedge. We observe two distinct behaviors: for 0 ≤ x c /α 1 2 the peak is still located at t = 0, while its value decreases as x c / is increased. In this range of x c /α, the A lies in the region where the quench is supported, however as x c / increases, less and less energy is initially contained within A. For 1 2 x c /α, most of the quench energy is initially supported outside A so the peak moves and is now located at t ∼ x c , while slowly decreasing in amplitude. Similar to the case of centered intervals, the late-time evolution JHEP08(2020)152 of the geometric correction for the non-centered intervals resemble that in an expanding boost-invariant plasma [76,77]. In particular, the non-monotonic behavior as we vary x c is reminiscent of the non-monotonic behavior of the correlators with respect to the rapidity [76,77]. The comments made in section 4.1.1 about emergent planar symmetry and 'first law'-like behavior also apply here.

First law of entanglement and CFT modular Hamiltonian
For small perturbations of a reference state, ρ = ρ (0) + δρ, entanglement entropy satisfies a first law relation, where H A is the so-called modular Hamiltonian. To see this, recall that by definition this operator is related to the reduced density matrix . (4.46) The small perturbation of the full state translates generically into a small perturbation of the reduced density matrix ρ A = ρ (0) A + δρ A . Hence, to linear order in this perturbation, the variation of entanglement entropy S A = −tr[−ρ A log ρ A ] is given by (4.47) The last term in (4.47) is identically zero, since the trace of the reduced density matrix is one by definition. Hence, the leading order variation of the entanglement entropy is given by (4.45), as advertised. However, there are very few cases for which H A is known explicitly. The most famous example is the case where A is half-space, say x 1 > 0, and ρ corresponds to the vacuum state. In this case [80,81] For generic CFTs, this setup can be conformally mapped to the case where A is a ball of radius R, centered at x c , in which case [47,82] Now, local quenches are not perturbatively close to the vacuum. Although the energy injected is small, some of the eigenvalues of the full density matrix ρ will drastically differ from those of the vacuum ρ (0) due to the sharply peaked perturbation. However, in some limiting cases, the reduced density matrix of a subsystem may still satisfy ρ A = ρ (0) A + δρ A . This happens when there is an additional small parameter to carry out the expansion, e.g., when the size of the subsystem A is small in comparison to other scales of the state. If this JHEP08(2020)152 is true, then the expectation value of the energy density operator would be approximately constant inside region A, T 00 (t, x) | A E(t), so (4.45) becomes Here, Ω d−2 = 2π is the surface area of a (d−2)-dimensional unit sphere. Defining δE A as the energy inside region A, where V A is the volume of region A, it follows that where T A is known as the entanglement temperature. This 'first law'-like relation holds true for arbitrary static states provided that E R d 1 [78,79]. For time-dependent states, equation (4.52) is still expected to be valid if additional conditions are satisfied: in addition to E R d 1, the size of the region must be smaller than all characteristic time scales of the state, i.e.Ė(t)R d+1 1,Ë(t)R d+2 1 and so on [70]. In the following, we will specialize to the case of 2d CFTs and take R = /2 as the half-length of the interval.
In our setup of a local quench in 2D CFT, we can distinguish between three scenarios: (a) R α and ∀ {t, x c }, (b) R t and ∀ {α, x c } or (c) R x c and ∀ {t, α}. These three cases are depicted in figure 9. The common feature among them is that the bulk region where the backrection is large is parametrically far away from the entanglement wedge of region A. Since the gravitational potential of a point particle decays with the distance, the small parameter in each case guarantees that the metric around the entanglement surface (i.e., near the boundary) is i) perturbatively small in such parameter and ii) approximately spatially homogeneous. The former condition implies that E R 2 1, while the latter one ensures that T 00 (t, x) | A E(t), as required for the derivation of (4.52). Now, for the state that we are considering, all the time variations are smooth provided that α is finite. It is only in the strict limit α → 0 that E(t) approaches the sum of two delta functions (3.42). However, in this limit the two conditions above are not satisfied so (4.52) is not expected to hold anyway. For any finite value of α, we can chose a sufficiently small R such that the extra conditions are satisfied:Ė(t)R 3 1,Ë(t)R 4 1, and so on, so that we can expect (4.52) to be valid in a corner of the space of parameters.
Let us check in detail the above claims. At this point it will be useful to recall the analytic expressions for the stress-energy tensor in the CFT (3.39) and entanglement entropy for centered and non-centered intervals, (4.24) and (4.34), respectively. As explained above, either if we are in situation (a), (b) or (c), we can assume that the energy density T 00 (t, x) is approximately constant in the region A, namely for x ∈ (x c − R, x c + R). For simplicity, we will evaluate it at the center of the interval, Performing the integral of the kernel over the region A, the first law (4.45) then implies that where . (4.55) Observe that this expression matches the leading order term of the entanglement entropy that we obtained from the direct computation in (4.26) and (4.40), for centered and generic intervals, respectively. This concludes our check.
Going beyond the leading order in the size of the region, we still expect the first law (4.45) to give the full result at linear order in the density matrix, O(δρ), but receive corrections at higher orders in δρ. Since we know the exact result for expectation value of the energy density operator (3.39), it is easy to compute this contribution. The easiest way to perform the integral is to change variables, x → ζ = x − x c + t so that the integrand becomes: (4.56) These two integrals can be written in terms of the arctan(z) function. Using the sum identity arctan(z 1 ) ± arctan(z 2 ) = arctan( z 1 ±z 2 1∓z 1 z 2 ) we can massage the full result into the following form: (4.57) This expression matches exactly the universal term (2.69) obtained from the CFT calculation in section 2. Expanding the above for small R, we obtain a complete match with the full first series of (4.26) and (4.40), which can be checked order by order in the small R expansion. These series can alternatively be written in a compact form by considering the integral expression for δS P A presented in (4.43), thus also matching the universal term JHEP08(2020)152

Corrections due to bulk entanglement entropy
The second class of O(1) corrections to entanglement entropy are due to the entanglement of bulk fields in the given quantum state. The bulk is state specified by the density matrix ρ bulk , defined on a Cauchy slice Σ ⊃ A ∪ γ A . Defining Σ A as the codimension-one region bounded by γ A and A, and Σ A c = Σ\Σ A , one then performs the trace over Σ A c and compute the von Neumann entropy associated to the reduced density matrix ρ Σ A = Tr Σ A c (ρ bulk ), S bulk = −Tr(ρ Σ A log ρ Σ A ). In practise, this can be accomplished by implementing the replica trick in the bulk, as was done in [57]. Inspired by the CFT discussion around equation (2.21), we will want compute the difference between the entanglement in the one-particle excited state |ψ , defined in (3.5), and the vacuum state |0 : bulk . This quantity is free of UV divergences. Given the crucial success of the first law of entanglement in the CFT, not only in reproducing the leading behavior of the geometric correction to entanglement entropy in the limit of small intervals but also in recovering the full universal piece of the entanglement entropy, it will prove useful to directly apply the first law of entanglement in the bulk to isolate the contribution coming from the entanglement of bulk fields. This first law stipulates that where H bulk is the bulk modular Hamiltonian. We emphasize that equation (4.59) is not expected to hold generally for the states we consider here, but only at leading order in bulk . The arguments follow closely to those for the first law in the CFT, as explained in section 4.1.1. Nevertheless, this should suffice to capture the difference between the universal term (2.69) expected from CFT considerations and the full result for the geometric corrections, obtained in section 4.1. In particular, we should be able to extract an infinite series of terms proportional to R 2∆+2i (i ∈ N) that are needed to cancel the second series of (4.26) and (4.40). This expectation is based on the exact relation between the bulk modular Hamiltonian H bulk and the CFT modular Hamiltonian H A [83] In the following section we will use the first law (4.59) to compute the contribution of the bulk entanglement entropy at linear order in δρ bulk , which should enable us to verify the above claims. We relegate the study of the next order corrections to section 4.2.2 and the interpretation in the context of the black hole information problem to section 4.2.3.

Linear order corrections and bulk modular Hamiltonian
To calculate the contribution at linear order in the density matrix δρ bulk , we will need an explicit expression for the bulk modular Hamiltonian in Poincaré coordinates: 61) where N µ is the unit normal vector associated to Σ A , ξ ν is the Killing vector that generates the entanglement wedge (normalized such that the surface gravity at the bifurcate horizon JHEP08(2020)152 Figure 10. Calculation of the quantum correction to the entanglement entropy associated to a region A = {x c − R ≤ x ≤ x c + R, t = t c }. Given a quantum state on Σ = {M : t = t c }, specified by its density matrix ρ Σ , we define Σ A ⊂ Σ as the codimension-one bulk region within γ A and A, and Σ A c = Σ \ Σ A . The entanglement entropy in the bulk is computed by tracing over the degrees of freedom in Σ A c and then finding the von Neumann entropy associated to the reduced state ρ Σ A . γ A is κ = 2π), and T µν is the bulk stress energy tensor. For simplicity, we will take Σ to be constant-t slice, as depicted in figure 10. This implies that The Killing vector ξ A takes the following form in the Poincaré patch [84] which is already correctly normalized since the surface gravity at the horizon yields κ = 2π. The expression (4.63) simplifies further for a constant t = t c slice and we get Defining the coordinates (ρ, φ) such that x = x c + ρ cos φ, z = ρ sin φ (0 ≤ ρ ≤ R, 0 ≤ φ ≤ π), we arrive to the following expression for the change of entanglement entropy in the bulk, to linear order in the density matrix T 00 (t, ρ, φ) .

(4.65)
This expression resembles the change of entanglement in flat space, but with an extra factor of L/z that arises from the volume form in AdS. We now need the expectation value of the energy density operator T 00 in the Poincaré frame to evaluate (4.65). In the original global coordinates, we have explicit expressions JHEP08(2020)152 for the components of the stress-energy tensor in (3.8), while the coordinate transformation from global to Poincaré coordinates x µ = (τ, r, θ) → x µ = (t, x, z) is given by equations (3.28)- (3.30). The standard transformation rule gives The explicit expression is lengthy and not particularly illuminating, so we will not display it here. For our purposes it will suffice to analyze the small interval limit of (4.65). Expanding the evaluated form of (4.66) for small ρ, we obtain a leading order term proportional to T 00 ∼ O(ρ 2∆−2 ) and further corrections suppressed by higher powers of ρ. Plugging this expansion into the equation (4.65) leads to and cancels the whole series O(R 2∆+2i ), i ∈ N, in (4.40). This is expected from field theory considerations, and the exact relation between bulk and CFT modular Hamiltonians (4.60). For x c = 0 this series can be resummed as with a(t) given in (4.18). For x c = 0 we do not know the most general form but we can formally write it in a compact integral form as where δS Q A is given in (4.44). We can obtain explicit expressions for integer ∆, by considering the difference between the universal term (2.69) and the result coming from the area term (4.32) (e.g. (4.34) for ∆ = 2), however the expressions are lengthy and we do not show them. JHEP08(2020)152

Quadratic corrections and bulk replica trick
We now go further and compute the vacuum subtracted bulk entanglement entropy to higher orders in δρ bulk . One way to do this is by implementing the replica trick in the bulk where ρ ψ = tr Σ A c |ψ ψ| and ρ 0 = tr Σ A c |0 0| are the reduced density matrices in the oneparticle excited state and the vacuum respectively. We start by expanding the operator 10 trρ n ψ = tr(ρ 0 + δρ) n = trρ n 0 + n trρ n−1  where δρ 1 ≡ ∂ n (δρ n ) | n=1 and we have used some formal properties of δρ n . 11 The calculation of this term is more naturally performed in the Hilbert space of the Rindler observer, which leads to a relatively straightforward answer in terms of the Bogoliubov coefficients that relates the global and Rindler modes [57]. We relegate the computation of the Bogoliubov coefficients for our combined transformation (global → Poincaré 10 The second order term has an interesting similarity with the sum (2.61) if one relates each O 2 insertion on a given sheet with a δρ factor and the index m with the distance in number of sheets between such operators. In particular, the factor of n/2 is explained in the exact same way as in (2.61): the replica symmetry ensures that for fixed m one has n equal contributions related by an overall translation, while the extra 1/2 is introduced to avoid the double counting of equivalent configurations. 11 We assume that δρn has a proper analytic continuation for n ≈ 1. In that case one can show that trδρn = (n − 1) trδρ = 0, and δρ1 ≡ lim n→1 δρn = 0 .
The first property follows from the cyclicity of the trace and the normalization of both, ρ and ρ0. The second is a consistency condition of the n → 1 limit of (4.72). Namely, since every order contribution to the left hand side of (4.72) must vanish in the n → 1 limit, this requires that This must hold for arbitrary δρ and ρ0 and therefore δρ1 = 0 must hold as an operator equation.

JHEP08(2020)152
→ Rindler), given by (3.28)-(3.30) and (4.80)-(4.82), to appendix B. Here we will merely point out that, in the small R limit, the Bogoliubov coefficients in our case differ only to those in [57] by a constant phase and hence, the resulting integral coincides with the one in [57], up to the identification of the corresponding small parameters (B.14). With this in mind, we obtain (4.77) and As in [57], we could not get a closed expression for the integral but by numerical evaluation we can confirm that it is consistent with the expected result, which matches with our formula for the dynamical contribution δS dyn A obtained from the CFT analysis (2.72). This completes our check of the FLM formula in the dynamical setting of local quenches. We emphasize, again, that our result for δS (δρ 2 ) bulk (and similarly for δS dyn A ) is only valid at the leading order in the small R expansion, while our result for δS (δρ) bulk is valid for any R.

Interpretation as the entropy of Hawking radiation
We now study the behavior of δS bulk and interpret the results. In our interpretation, two transformations play an important role: the large gauge transformation in equations (3.28)-(3.30) and the CHM map that takes a Poincare patch to a Rindler patch. The CHM [47] map is given by where r ∈ (1, ∞), t ∈ (−∞, ∞) and x ∈ (−∞, ∞). In this coordinate system, the bulk metric takes the form of a planar BTZ geometry It is useful to study the effect of the composite map created from transformations (3.28)-(3.30) and the CHM map (4.80)-(4.82). When we start in global coordinates with a static, single-particle wavefunction for the bulk excitation, this map takes it to a multi-particle time-dependent wavefunction. This is because of a Bogoliubov transformation implicit in it, see appendix B. We thus get a thermally-populated state of the scalar field in the planar BTZ geometry (4.83). We interpret these excitations as Hawking radiation outside the BTZ black brane. Consequently, when viewed from the Rindler frame, S bulk computes the entanglement entropy of bulk fields with the interior of the brane traced over. In black hole context, this is often referred to as the entropy of the Hawking radiation.
The dependence of the multi-particle wavefunction on Rindler time t has an interesting consequence. As a result, the boundary condition for the scalar field in the CHM frame is time-dependent. This allows a leakage of the wavefunction from the wedge described by (4.83) to its complement. We interpret the complement as a bath that absorbs the Hawking radiation from outside the black brane (4.83). This should be contrasted with the setup in [3] (and related works) where a non-gravitating region needs to be coupled to the black hole to absorb the Hawking radiation. The global frame point of view discussed later further supports these claims.
In fact, if we go to the global patch, we can make precise the sense in which our model is a toy model of black hole evaporation and also see the absorption of Hawking radiation in a different way. First notice that the black brane (4.83) is actually time dependent. There are different Rindler wedges associated to the region A at different time slices. In the Poincaré frame, these Rindler wedges are adapted to a region of fixed length . However, they have different sizes from the point of view of the global frame. For centered intervals, x c = 0, a careful calculation shows that otherwise . (4.84) δτ ≡ |τ 2 − τ 1 | = 0 , (4.85) in the boosted global frame, or (4.86) in the original global frame. The reason that δθ and δθ have two branches is due to the periodicity of the two coordinate systems, θ ∼ θ + 2π and θ ∼ θ + 2π. In both cases we show that the size of the RT surface in fact decreases in time as t goes from t = 0 to t → ∞, as depicted in figure 11. In the left figure, we show an example as seen in the boosted global JHEP08(2020)152 frame with R < L so we only have the second branch of δθ for t ∈ [0, ∞]. In this case we have δθ | t=0 < π and decreases monotonically as t evolves, with δθ | t→∞ → 0. In the right figure, we show an example as seen in the original global frame. In this case, one has to always include both branches of δθ. In fact, the branches are switched at exactly the Page time since then δθ| t=t Page = π. This implies that δθ| t<t Page > π and δθ| t>t Page < π, as can be confirmed from (4.86). With the inclusion of both branches, δθ decreases monotonically for t ∈ [0, ∞] also in this frame. Hence in either boosted or original global frame, the size of the black brane always decreases in time for centered intervals. We interpret this as a toy model of an eternally evaporating black hole. We can also see the absorption of the Hawking radiation mentioned earlier while working in the original global frame, where the wavefunction for the bulk excitation is static and spherically-symmetric. It is useful to identify the region which maps to the exterior of the BTZ black brane under the aforementioned composite map. In the global frame, the overlap of the wavefunction with this region is proportional to the probability of the particle being in the exterior of the black brane. However, the size of that region decreases monotonically as plotted in figure 11. This decreasing probability fits in with our interpretation of the leakage of Hawking radiation. Conversely, as the size of the complementary region increases, the probability of the particle to be in the complement increases with time. This further justifies our interpretation that the interior of the black brane serves the purpose of a bath.
In our toy model of black hole evaporation, we can ask how the entropy of the radiation evolves with time. Let us consider the case where the multi-particle wavefunction is initially inside the entanglement wedge (i.e. exterior of the brane). This is given by the condition In this case, the state as seen by the Rindler observer is a thermal perturbation of the black brane geometry, which is a quench of a thermal state in the dual CFT. The entropy of Hawking radiation is given by the bulk entanglement entropy δS bulk . For concreteness, we will only consider the expressions for δS (δρ) bulk , which are valid for any R (recall that our result for δS (δρ 2 ) bulk requires R α, so it is not valid in the regime where (4.88) holds true). In figure 12 we plot δS (δρ) bulk for different scenarios when (4.88) is satisfied. In the left panel, we plot δS (δρ) bulk for the case of centered intervals, x c = 0, and different values of ∆. Strikingly, we observe that the bulk entanglement at this order follows the expected behavior for a Page curve, i.e., increasing up to a time t = t Page and then decreasing as t → ∞. The Page time t Page in our setup can be identified as the time at which the particle exits the wedge, We can perform an analysis similar to that in equation (4.86) for non-centered intervals. It is straightforward so we will omit the specifics here and directly state our results. In this case, the size of the black brane as viewed from the global perspective changes nonmonotonically as we vary the Poincaré time from t = 0 to t → ∞. We take this to mean that for x c = 0 our setup cannot be interpreted as a toy model of black hole evaporation. JHEP08(2020)152 δθ| t=t Page δθ| t>t Page Figure 11. Evolution of the size of the RT surface (i.e. black brane horizon from the perspective of the Rindler observer) for centered intervals from the point of view of the global frame. Left: we plot the situation as seen from the boosted global coordinates. Right: we plot the situation in the original global frame. In both cases we observe that the size of γA (hence the size of the black brane) decreases monotonically as the Poincaré time evolves from t = 0 to t → ∞. This implies that, for centered intervals, our setup can indeed be interpreted as a toy model of black hole evaporation. We remark that δS (δρ) bulk indeed behaves as a Page curve for t ∈ [0, ∞], reaching a peak at t = tPage and then decreasing as the state evolves in time. In our setup, the Page time tPage can be identified as the time at which the particle crosses the RT surface γA.
The plot in the right panel of figure 12 confirms our justification. In this case, the evolution of S (δρ) bulk does not follow a standard Page curve, but instead, develop two peaks at different times t 1,2 > 0.
Before closing this section, we comment on the implications of the above analysis for the information problem. Since our state is pure, the von Neumann entropy of a full time slice vanishes, S Σ = 0 and the Araki-Lieb inequality implies For excited states, we can subtract the vacuum contributions to obtain When applied to our discussion of decreasing entanglement wedges in the original global frame, we can easily explain the Page curve in the left panel of figure 12 in terms of the purity of the state. Recall that δθ| t<t Page > π while δθ| t>t Page < π. If we now replace δS (δρ) bulk for t < t Page by the entanglement of the complementary region Σ c A using (4.91), the increasing and decreasing parts of the Page curve immediately follow from the monotonicity of δS (δρ) bulk for δθ ∈ [0, π] (see red dashed line in the right plot of figure 11). Although this result is only manifest in the unboosted global frame, it must be true in general, even in the Poincaré frame. We thus conclude that, at least in our toy model for black hole JHEP08(2020)152 bulk , for cases that satisfy the condition (4.88). The latter condition is imposed such that the particle is initially inside the entanglement wedge so for a Rindler observer the state can be viewed as a local perturbation of a black brane geometry. In the left panel we have plotted the case of centered intervals, so we have set x c = 0, and we have varied ∆ = 3/2, 2, 5/2, 3, 7/2, corresponding to the red, orange, green, blue and purple curves, respectively. In the right panel we have plotted the case of non-centered intervals, specializing to ∆ = 2 and varying x c / = 0, 10 −1 , 2 × 10 −1 , 3 × 10 −1 , 4 × 10 −1 , corresponding to the red, orange, green, blue and purple curves, respectively. For concreteness we have fixed α/ = 10 −1 in all the plots and, accordingly, chosen values for x c / such that (4.88) is satisfied. evaporation, the behavior of the Page curve can be directly associated with the purity of the global state.
Interestingly, the above result follows from a semiclassical calculation in the bulk, without the requirement of the so-called island contribution. This is due to the fact that our quantum state is simple enough so that there are no competing saddles as time evolves. Nevertheless, the mere fact that we are finding a Page curve consistent with unitarity in the framework of semiclassical gravity implies that, in more realistic models, information loss or the lack thereof should not be interpreted as an artifact of the semiclassical approximation.

Conclusions and outlook
In this paper, we studied the evolution of entanglement entropy after a local quench in twodimensional conformal field theories with holographic duals both from the CFT perspective and using its gravitational description. Our CFT computation was carried out using the replica trick. This entailed the calculation of a 2n-point correlator on a topologically nontrivial multi-sheeted Riemann surface. The computation was considerably simplified by virtue of the enhanced conformal symmetry special to two dimensions which allowed a map that relates the 2n-point correlator on the aforementioned Riemann surface to a 2npoint correlator on the complex plane. As a consequence, a separation of the entanglement entropy into two pieces (universal and dynamical) emerged naturally. We provided insights into their nature from the CFT perspective. In the bulk, the aforementioned calculation amounted to compute two separate contributions, the area of an extremal surface in the backreacted geometry and the bulk entanglement entropy associated to the bulk fields in JHEP08(2020)152 the excited state. We carried out explicit computations of both contributions and showed that the sum matched the CFT result, thus, providing a non-trivial check of the FLM proposal in a fully dynamical situation.
The calculation of the bulk entanglement entropy piece was particularly interesting when interpreted in the context of the black hole information paradox. Due to the purity of the bulk quantum state, the bulk entanglement entropy of a subregion in a compact space must follow a Page curve as a function of size. Interestingly, under various coordinate transformations we could relate the bulk entanglement entropy of a fixed region in the quench state with the bulk entanglement entropy of a region that decreases with time from the point of view of the global state. For fixed time, this could be further associated to the entropy of bulk fields on a planar BTZ black hole geometry. Thus, one could relate the Page curve that follows from purity of the bulk quantum state with the Page curve associated to the entropy of Hawking radiation in an eternally evaporating black hole background. Is in this sense that our setup could be interpreted as a simple toy model for unitary black hole evaporation. We emphasize that, even though our bulk calculation was semiclassical, we did not require the so-called island contribution to restore unitarity. This was the case because our particular quantum state was simple, leading to only one saddle point dominating the path integral throughout the evolution. By simple, we mean states with an order O(1) number of excitations in the bulk which change the total energy and entropy by an order O(1) amount. For more generic complex states this is not necessarily the case, and one must consider multiple saddle points that eventually exchange dominance. It is nevertheless interesting to note that unitarity can be preserved within the semiclassical approximation, at least in simple setups as the one considered here. This implies that information loss or the lack thereof is, in fact, not an artifact of the approximation.
There are many open problems and generalizations, some of which we discuss now: • Higher dimensions. Many of the steps used in the bulk computations can be generalized for boundary spherical regions or bulk hemispheres on higher dimensional AdS spaces. However, our CFT computation relies heavily on the enhanced conformal symmetry special to two dimensions. To generalize our analysis to higher dimensions, one could imagine instead using higher dimensional twists operators in the CFT, as the ones introduced in [85] to mimic the non-trivial geometry of the replica manifold. Further approximations might be required to gain analytic control. A natural would be the small radius limit of the entangling sphere which was successfully used in the analysis of [85] and which nevertheless could provide a non-trivial match with the bulk computations.
• Relative entropy and mutual information. Relative entropy and mutual information are interesting quantities from quantum information theory perspective and have been widely studied in recent works. There exist replica like techniques for the computation of these quantities in QFT as well as concrete proposals for their holographic duals [83]. Therefore, it would be very interesting to extend our framework to the computations of these quantities in our set up. For instance, the mutual information for disjoint regions in excited states was computed recently in [86], and relative en-JHEP08(2020)152 tropy for similar states was discussed in [87], for general CFTs. See [88] for a recent dynamical analysis of mutual information.
• More general local quenches. One can consider local quantum quenches of a different type compared to the ones studied in this paper. For example, to consider operators with spin. Recently the backreaction of operators with spin was considered in [89], which generalized the work of [57] to more general one-particle states. One can imagine performing the same set of bulk transformations that we considered here to obtain quench states in these scenarios. Another option would be to consider a reference state that is not the vacuum. For instance, if the unperturbed state is taken to be a thermal state, then the gravity dual would be a localized perturbation of a black hole geometry [16]. This setup is interesting because from the entanglement entropy calculation one could obtain the leading quantum corrections to the butterfly velocity v B , e.g. by extending the work of [90] to order O(1). Finally, one could consider bilocal quenches as in [17,18]. These are obtained by inserting a local operator at two different points in space and evolving the state with the Hamiltonian. In particular, computing quantum corrections in these scenarios could shed light on the nature of gravitational interactions beyond the classical regime.
• Quantum extremal surfaces. Generalizing beyond the leading O(1) corrections to holographic entanglement entropy, there is a proposal for an all-order resummed result in terms of the so-called quantum extremal surfaces [7]. The aforementioned surfaces are obtained as the extremal surface with respect to the generalized entropy as opposed to the area term, as in the FLM prescription. Elevating our bulk calculation to this next level of complexity requires the understanding of the bulk entanglement entropy for general regions or arbitrary small perturbations of the semicircle. One would also need to find the various maps of interest for slightly deformed regions. An interesting work along these lines is the perturbative analysis in [91] for the modular Hamiltonian of deformed spheres and its generalizations for the bulk modular Hamiltonian following [92]. It would be interesting to work out this example in detail.
We hope to come back to some of these problems in the near future.

JHEP08(2020)152
γ n Im(s) = 0 Im(s) = 2πn Figure 13. Integration contour (depicted in solid gray) used to evaluate the complex integral in (A.7). Assuming that integrand vanishes for Re(s) → ±∞, we can then deform the contour γ n to just be the dashed gray lines at Im(s) = δ and at Im(s) = 2πn − δ with δ > 0. For illustrative purposes, we have picked the value n = 4 to make this figure.
JHEP08(2020)152 this coefficient, we divide (A.11) by (n − 1) and take the n → 1 limit on both sides, wherek (s) ≡ − s 4 sinh 2 (s/2) . (A.13) Finally, we deform the integration contour to Im(s) = π in the first term and to Im(s) = −π in the second term, so that the δ → 0 limit is non-singular. By doing so, we obtain an expression that can be readily integrated . (A.14) This leads to the result reported in (2.62).

B Bogoliubov coefficients and bulk entanglement
The scalar field has two different mode expansions depending on the choice of coordinates: in global coordinates we can write The creation and annihilation operators in the two frames are related by the Bogoliubov coefficients α, β, e.g., a m,n = I,ω,k α m,n;ω,k,I b ω,k,I + α * m,n;ω,k,I b † ω,k,I , (B .3) and similarly for a † m,n . We will compute the asymptotic form of these coefficients, focussing on the large-r limit. As we will see below, this will indeed suffice to extract the leading result for δS (δρ 2 ) bulk in the small R limit. Before going further, we need to determine how does the state |ψ = a † 0,0 |0 look in Rindler coordinates. This was worked out in detail in [57], so we will only quote the result here |ψ = a † 0,0 |0 = ω,k (1 − e −2πω )α * ω,k,R b † ω,k,R + (1 − e 2πω )β ω,k,R b ω,k,R . |0 (B.4) Notice that this relation only involves Rindler modes associated to the right wedge and so it is very convenient for the computation of the reduced density matrix obtained after JHEP08(2020)152