Dynamics of the Area Law of Entanglement Entropy

We study the evolution of the universal area law of entanglement entropy when the Hamiltonian of the system undergoes a time dependent perturbation. In particular, we derive a general formula for the time dependent first order correction to the area law under the assumption that the field theory resides in a vacuum state when a small time-dependent perturbation of a relevant coupling constant is turned on. Using this formula, we carry out explicit calculations in free field theories deformed by a time dependent mass, whereas for a generic QFT we show that the time dependent first order correction is governed by the spectral function defining the two-point correlation function of the trace of the energy-momentum tensor. We also carry out holographic calculations based on the HRT proposal and find qualitative and, in certain cases, quantitative agreement with the field theory calculations.


Introduction
In recent years, the study of the time-dependence of entanglement entropy has received a considerable attention [1][2][3][4][5][6][7][8][9][10] . In the simplest scenario one studies the evolution of the spatial entanglement entropy under the assumption that the system begins in the vacuum state in the asymptotic past while a relevant (i.e., positive mass dimension) coupling constant in the Hamiltonian is given a non-trivial time dependence which asymptotes to a constant value in the far future.
Although the assumption that the system starts in the vacuum state might seem to be an oversimplification, certain questions about the entanglement may still be addressed without loss of generality. For example, in quantum field theory the entropy of a region is always UV divergent [11,12], whereas the structure of these divergences is state independent, see e.g., [13]. Hence, assuming that the system resides in the vacuum state is enough to determine the pattern of divergences in any state. Indeed, this pattern is fixed by the entanglement of degrees of freedom near the UV cutoff scale, δ, and therefore the characteristic scale of an excited state plays no role in this regime. Furthermore, the dynamics of the divergent part of the entanglement entropy in a time dependent setup is not trivial, since the asymptotic vacua, in general, are not identical.
Even though the time dependence of the full entanglement entropy in a generic setup is very complicated, we expect that evolution of the UV divergences is quasi-adiabatic, and therefore tractable. Our expectation rests on the locality of UV divergences when all length scales, ℓ i , including the characteristic time scale of the change in couplings, are much longer than the UV cut off, ℓ i ≫ δ, i.e., when the system is found in the physical regime. 1 In particular, the divergent part of the entropy at time t will only depend on the instantaneous value of the coupling and its time derivatives.
There is an additional outcome that stems directly from the locality of divergences. It is now a well-established fact that entanglement entropy is sensitive to the geometry of the background and shape of the entangling surface [14][15][16][17][18][19][20][21][22]. The full answer can be expressed as a sum of certain geometric structures with coefficients fixed by QFT data. Not all of these structures are necessarily local, however those which are multiplied by divergent coefficients as δ → 0 should be local.
The ubiquitous 'area law' [11,12] is probably the most prominent and simple example of a local geometric structure. In a static d-dimensional CFT it is proportional to δ 2−d A since δ in this case is the only admissible scale at hand, while entropy is dimensionless and area of the surface, A, has length dimension d − 2. Therefore the 'area law' term is not universal in a CFT -it is sensitive to the choice of regularization scheme. However, the presence of an additional scale in a non-conformal field theory may overturn this conclusion and result in a universal 'area law' of entanglement, e.g., a logarithmically divergent term proportional to the area multiplied by an appropriate power of some built-in scale.
In this work we focus on the evolution of the 'area law' terms when a QFT undergoes a quantum quench, and calculate the time-dependence of the universal coefficients in a variety of examples. To make the field theory analysis analytically tractable, we expand the entanglement entropy to linear order in the time dependent coupling λ(t). This is not necessary to recover the exact pattern of the universal divergences in the holographic approach, since it is completely fixed by the asymptotic region of the dual geometry. Even though our focus is on the 'area law' only, we expect qualitatively similar behavior from the subleading divergences as well.
If the mass dimension of λ(t) is α, then for α = (d − 2)/k we may have a universal term proportional to λ k (t)A. Similarly, if α = (d − 2 − l)/k with l being a non-negative integer, then a term proportional to λ k−m (t)∂ l t λ m (t)A with 1 ≤ m ≤ k is allowed. It may happen that for a given α there is more than one choice of l and k, in which case there is a proliferation of possible universal structures. Generically, one would expect all allowed terms to appear, and with independent coefficients. Note, however, that for a given α there are only finitely many allowed terms.
In sections 2 and 3 we carry out explicit calculations for a set of special values of k and l in holographic and free field theories. In both cases we observe the quasi-adiabatic behavior that we had argued for. Moreover, in some cases our holographic and free field calculations exactly match. Such an agreement between the time-dependent terms of entanglement entropy in the weakly and strongly coupled regimes resembles similar match in the static case [23,24] and can be attributed to the universality of certain two-and three-point correlation functions in a CFT. In particular, we observe explicit match when α = 1, 2 and k = 2 in various dimensions. 2 These values of α correspond to the Dirac and scalar mass operators, respectively. Recalling now that α = (d − 2 − l)/k results in the following cases explicitly elaborated in the text: • α = 1, d = 4, l = 0. This case corresponds to a 4D massive Dirac field. Up to a numerical factor the universal entanglement entropy is given by ∼ λ 2 (t)A, where λ(t) is the time-dependent mass of the Dirac field.
• α = 2, d = 6, l = 0. This case of a 6D massive scalar field is analogous to the previous case, i.e., the entanglement entropy is given by ∼ λ 2 (t)A, except that λ(t) is mass squared this time.
• α = 1, d = 6, l = 2. In six dimensions, there is a universal term involving two time derivatives of the Dirac mass ∼ λ(t)λ(t)A. We calculate the numerical coefficient for the first time both holographically and in the free field theory, and find exact agreement.
• α = 2, d = 8, l = 2. The behavior of the universal term for the massive scalar field in eight dimensions ∼ λ(t)λ(t)A resembles the massive Dirac field in six dimensions, and again we calculate the precise numerical coefficient holographically and in free field theory, finding agreement.
The remainder of this paper is organized as follows. In Section 2 we build on the holographic proposal [26,27] to compute the time-dependent universal entanglement entropy for holographic field theories. Since we are only interested in divergent terms, the gravity calculation is localized in the asymptotic region of the dual geometry; this feature allows us to calculate the universal coefficients exactly even though the geometry is only determined perturbatively. We calculate the coefficients of the terms for which agreement with the free fields is expected, as well as for the terms where agreement is not anticipated. In the latter case the answer on the field theory side depends on more than just the universal correlation functions of the underlying CFT.
In Section 3 we derive a general expression for the time-dependent first order correction to the 'area law', and show that for a generic QFT it can be expressed in terms of the spectral function defining the two-point correlation function of the trace of the energymomentum tensor. We use these formulas to calculate the universal 'area law' for massive free fields with time dependent masses. We end with conclusions in Section 4, and a number of Appendices detailing some of the more technical aspects of the calculation.

Setup
We consider a CFT in d spacetime dimensions on the boundary, dual to Einstein gravity in the d + 1-dimensional bulk. The computation of the entanglement entropy of some boundary region, B, can be done using HRT prescription [26,27]. The prescription tells us to find the codimension-2 surface in the bulk with extremal area, Σ B , subject to the condition that B and Σ B are homologous. If there are multiple extremal surfaces, we must choose the one with the minimum area. We will call this minimal extremal surface the HRT surface from now on. The entanglement entropy of B and the area of the HRT surface are related by Beginning with the vacuum state of the CFT, we turn on a relevant operator O of dimension ∆ < d at t = 0 by giving it a nonzero, time-dependent, spatially homogenous coupling λ(t). It will be useful to express our results in terms of the mass dimension of the coupling, α ≡ d − ∆. The operator O is dual to a scalar field Φ in the bulk of mass m 2 = ∆(∆ − d) [28,29]. A nonzero coupling λ(t) corresponds to nontrivial boundary conditions for Φ. Through the Einstein-scalar field equations, this will change the bulk metric and therefore affect the calculation of the HRT surface and its area.
The Einstein-scalar equations are given by the variation of the bulk action [13] where we have set L AdS = 1. We have chosen to normalize the bulk scalar action so that there is an overall factor of 1/16πG d+1 in front; this amounts to a normalization condition on the operator O. In particular, with this normalization, we can read off Newton's constant from the two-point function of O [30]: This will allow us to make contact with the free field calculations in Section 3. 3 The bulk interaction parameters κ, ω, etc., are similarly determined by the three, four, and higherpoint functions of O. The most interesting terms in the entanglement entropy are the ones which depend only on the two-point function, since these may be the same across theories once the normalization has been fixed using (2.4). Thus we should pay special attention below to those terms which depend only on G d+1 . The coupled Einstein-scalar equations are Following [31], we choose our bulk 'radial' coordinate z to be orthogonal to all of the boundary coordinates. This corresponds to the gauge fixing where x i for i = 1, ..., d − 1 are the spatial coordinates of the boundary. We also rescale the bulk coordinate to fix g zz = 1 z 2 . Since our quench is homogenous and isotropic, the bulk geometry must be homogenous and isotropic as well. Then we can write the bulk metric as with the boundary conditions that f (t, z = 0) = h(t, z = 0) = 1, and f (t < 0, z) = h(t < 0, z) = 1.
With this choice of metric, (2.5) and (2.6) take the following form 4 : The UV divergences in the entanglement entropy arise from the divergences of the area functional of the HRT surface, Σ B , near the asymptotic boundary (z = 0). To compute these divergences, we need to find the area functional, or the induced metric on Σ B , near z = 0. In [13], the authors studied the solutions of (2.9)-(2.11) as a power series in z. They have shown that the solutions can be written as the sum of two separate power series. One of them is independent of the state of the boundary theory, whereas the other series solution carries information about the state of the boundary theory. They also showed that the latter solution does not contribute to the divergences in the entanglement entropy, which implies that the divergences in the entanglement entropy are state-independent (as expected).
Having learned that the state-dependent power series solution does not contribute to the divergences, we will only be writing the state-independent solution in the following analysis. 5 In [13], all the powers of z that arise in the asymptotic solution have been identified. Using their result, we make the following ansatz for the metric and the scalar field: where φ 0,0 (t) is equal to the time dependent coupling, λ(t). That is,

HRT Surface
In this work we will only consider terms which are independent of the curvature of the entangling surface. Then it will be enough to consider the entanglements of a half-space in the field theory, which means the entangling surface is a flat plane. That is, the region B in the boundary is given by The natural coordinates on the codimension-2 HRT surface in the bulk are {z, x i } for i = 2, ..., d − 1, with t and x 1 left as functions determining the position of the surface in the bulk. The HRT surface must be invariant under translation in x i direction, for i = 2, ..., d − 1. This means t = t(z) and x 1 = x 1 (z) are functions of z alone. The state is also translation-invariant in the x 1 direction, and so we can fix x 1 (z) = 0. To summarize, the HRT surface, Σ B , is with the boundary condition, t(z = 0) = t 0 , the time at which we are computing the entropy.
The area functional of the HRT surface is given by where A is the area of the entangling surface in the boundary, and prime denotes the derivative with respect to z. Also note that we have introduced a cut-off surface near the boundary at z = δ.
Applying the variational principle to the above functional gives us an equation for t(z). As mentioned above, we are only interested in the asymptotic behavior of the solution near z = 0. As we will show in Appendix A.1, the lowest nontrivial power of z that appears in the solution of t(z) is z 2+2α . That is, at the lowest order, the solution of the stationary surface is where t 2+2α is a constant.
By expanding the integrand of (2.18) in z, we see that we can see that t 2+2α in (2.19) enters as the coefficient of z 4α+3−d . This will not affect the divergent terms in the entropy as long as In the following, we will perform explicit calculations for α = 1, 2 in d = 4, 5, 6, and for α = 2 in d = 8. These values of α corresponds to the mass term of a Dirac fermion and scalar boson, respectively. For these values of α and d, the above inequality is true. This allows us to forget about the details of the HRT surface, and simply set t(z) = t 0 in (2.18). With this simplification, (2.18) becomes The entanglement entropy of the boundary region B is then given by where we can give G d+1 a field-theory definition by relating it to the two-point function of the operator O as in (2.4), reproduced here: (2.23)

Scalar mass operator: α = 2
In this subsection, we will specialize to the case where α = 2. In the free field case, this corresponds to a time-dependent mass: The entanglement entropy is given by (2.22). The ansatz for h(z, t 0 ) given in By counting powers of z in the integral (2.22), we deduce that there is no time-dependent universal log-divergence in the entanglement entropy for d = 4 and d = 5. This result is consistent with the time-independent holographic calculations [13]. Now we turn to d = 6 and d = 8. To facilitate comparison to the free field case (wherever applicable), we will use (2.4) and the result [23] which says that the effective value of Newton's constant is Using (2.22), we get After inserting (2.26) in the above expression, and extracting the log-divergence, we get Using (A.7), we get This is an example of a result which depends only on Newton's constant and not any of the other parameters in V (Φ). Then it depends only on the two-point function of O, through (2.4). Using (2.28), we get G 7 = 12π 2 . Then the holographic result, extrapolated to the free scalar, is If m 2 (t) is constant, this agrees with the well-known answer for the static case [13]. We will see below in Section 3 that we have full agreement also in the time-dependent case.

d = 8
Using (2.22), we get After inserting (2.26) in the above expression, and extracting the log-divergence, we get Using (A.10) and using (2.15), we get where G 9 = 120π 3 is given by (2.28). Note that the coefficent of the λ 3 (t) term is proportional to κ and hence it's not expected to be universal for all unperturbed CFTs.

Fermionic mass operator: α = 1
In this subsection, we turn to the case α = 1. This corresponds in the free field limit to a Dirac mass for a fermion: The ansatz for h(z, t) given in (2.12) for α = 1 is By a power counting argument, we find that there is a log divergence for all d > 3. The absence of the log term in d = 3 is consistent with the known result [13,25] that the change in the entanglement entropy to first order in a perturbation away from a CFT is zero: the coefficient of the log term, if present, would have been proportional to m.
To determine the effective value of G d+1 using (2.4), we will use the following two The log divergence in (2.22) where h 2,0 (t 0 ) is given in (A.12): Using (2.4) and (2.41), we get G 5 = π 2 . Then we have In this case, the log divergence in the entanglement entropy is equal to We copy the result for h 3,0 (t 0 ) from (A.15): Using (2.4) and (2.41), we get G 6 = 8. Therefore, the entanglement entropy becomes In this case the answer depends on the bulk interaction strength κ, and hence is not expected to be universal across all theories.

d = 6
In this case, the entanglement entropy in (2.22) simplifies to By expanding (2.40) in the above expression and extracting the log term, we get Using (2.4) and (2.41), we get G 7 = 3 2 π 2 . Combining these results gives us the log term divergence which not only depends on the instantaneous value of the coupling, but also on the instantaneous values of its derivatives. Notice that the derivative term does not depend on the higher-point interactions. Putting together the above results gives

Field Theory Calculation
In this section, we use the Hamiltonian formalism to find the time evolution of the log term in the entanglement entropy for a general field theory to first order in perturbation theory. We will identify the linear term in the time-dependent entropy with a certain correlation function in the original theory. Using the free scalar field as a guiding example, we will go on to evaluate the divergent parts of this correlation function in terms of the spectral function c (0) (µ) of the theory, defined in terms of the two-point function of the trace of the energy-momentum tensor [32]. Finally, we use this general prescription to calculate the answer for a free Dirac field.

Entropy Perturbation as a Correlation Function
The time dependent Hamiltonian of our theory is where H 0 is the Hamiltonian of the unperturbed theory, and λ(t ≤ 0) = 0. We start with the system in the vacuum state, |0 , of the unperturbed Hamiltonian. The global state at any time t > 0 is given by the time-dependent density matrix,ρ(t) = |χ(t) χ(t)|, where |χ(t) solves the Schrödinger equation generated by H(t) with initial condition |0 . The time-dependent reduced density matrix of some spatial region B is given by the tracing over the Hilbert space of complementary region,B. That is We define the time dependent modular Hamiltonian K(t) of region B as where K(t) has only support in the region B. In other words, this operator acts as the identity operator on the Hilbert space on the regionB. The entanglement entropy of B can be written using the modular Hamiltonian: To extract the term linear in λ(t), we define the state in the interaction picture. That is where the interaction picture state solves the equation whereÕ(t) = e iH 0 t Oe −iH 0 t . We can solve the above equation perturbatively where the first order correction is We expand the modular Hamiltonian in the same fashion where K 0 is the modular Hamiltonian of the vacuum state of the CFT, and hence is time independent, and K (n) involves n powers of λ.
Combining these expansions we get the first order change in the entanglement entropy: The normalization of the reduced density matrix provides a constraint on the modular Hamiltonian: Expanding this constraint to leading order yields which simplifies (3.10) to This is one of our main results. The content of this formula is summarized in Fig. 1. The formula (3.13) is applicable for any spatial region B. However, the vacuum modular Hamiltonian K 0 is not known for a general region. We will therefore again restrict ourselves to the case where entangling surface is an infinite plane, so that the region B is a half-space. In that case, the modular Hamiltonian of the vacuum state is well known [33][34][35][36] (3.14) Here we have adopted the following notation for spatial vectors, v.
where v ⊥ is the component of v in the direction orthogonal to the entangling surface and v is the projection of v parallel to the entangling surface.
Putting this together, we see that the first-order change in the entropy is related to the commutator of the energy density and the perturbation: Notice that the interaction picture operatorÕ(t) is identical to the Heisenberg operator O(t) in the original λ(t) = 0 theory, i.e., when the total Hamiltonian is equal to H 0 . In the following, O itself will be the integral of a local operator, so we effectively reduce the problem to a commutator of local operators (or equivalently the imaginary part of a two-point function) in the original theory.

Example: free scalar field
In the static perturbation theory the first order change in the entanglement entropy after perturbing a CFT by a relevant operator vanishes [25]. It follows from (3.15) that similar conclusion holds in the time dependent case. To resolve this difficulty, we will start with the free scalar with non-zero mass, m, and then introduce a time dependent perturbation to the mass, m 2 (t) = m 2 + 2λ(t). 6 The action of our theory is 7 where λ(t ≤ 0) = 0. The perturbing operator is O(t ′ ) = dx ′ φ 2 (t ′ , x ′ ), where the time dependence of these operators is governed by the free Hamiltonian with λ = 0. From (3.15) we see that the entropy is determined by correlation functions of T 00 and φ 2 . The energy-momentum tensor of the free scalar field is given by where the last term is a possible improvement term. The energy momentum tensor is traceless if it is massless and ξ = ξ c ≡ d−2 4(d−1) , in which case the scalar field is said to be conformally coupled. The other noteworthy case is the minimally coupled scalar, where ξ = 0. The required correlation functions are easily computed using Wick's theorem: is the (unordered Lorentzian) two-point function of the scalar field. The expectation value of the commutator T 00 (x), φ(x ′ ) 2 is simply twice the imaginary part of the correlation function, up to a factor of i. Before taking the imaginary part, we first perform the required integrations, the details of which we leave to Appendix B. From (B.5) we find Using this result in (3.15) yields where Ω d = 2π d/2 Γ(d/2) . To isolate UV divergences in the above expression, we expand the integrand around Next we note that for any given d, the integral over x in the above expression results in a UV-divergent S (1) if and only if d − 4 − 2n ≥ 0. There are only a finite number of such terms, and therefore we regularize them by assuming that n is sufficiently large, carry out the integral over x and treat the special values of n by analytic continuation. The final result depends on the parity of d. For odd d there are no divergences, and S (1) is finite. 8 On the other hand, for even d we obtain + finite terms . 8 The non universal divergences vanish for our choice of regularization scheme.
Introducing a UV cut off surface t ′ = t − δ to ensures that the integral over t ′ is finite, and expanding λ(t ′ ) in the vicinity of the cutoff yields S (1) There are two points to make about this result: (i) It gives us the expected adiabatic expansion of the time-dependence of the divergent terms in the entanglement entropy, where in accord with local nature of the divergences only the instantaneous value of the coupling and its derivatives matter. (ii) Only even numbers of derivatives appear in the divergent terms, which is consistent with the holographic result.
For the time-dependent λ(t) in d = 4, 6, 8, the explicit form of the first order correction is Setting ξ = ξ c we find that terms proportional to m 2 λ and m 2λ in d = 6 and d = 8 respectively agree with their holographic counterparts. 9 Furthermore, linear terms in λ(t) and its derivatives vanish in all dimensions which is also in full agreement with holography. Such a match between the calculations for the free field theory and strongly coupled N = 4 SYM might seem surprising at first sight. However, as we argue in the next section, similarly to the time independent case [13,43], this match can be attributed to the universality of T µν O , OO and T µν OO in a conformal field theory [23,24]. We note that tuning ξ to a conformal value ξ c was important to ensure tracelessness of the stress tensor in a CFT.

Spectral representation
The scalar field theory example of the previous section suggests that the time dependent universal terms of EE can be easily identified if we rewrite (3.15) as an integral over energy states. The main goal of this section is to perform such spectral analysis for an interacting QFT which is driven by the action with only one operator, O(x), of scaling dimension ∆ < d. The corresponding relevant coupling will be denoted by λ. We assume that λ is equal to some fixed value plus a time-dependent perturbation: λ(t) = λ 0 + δλ(t).
We start from noting that our main result (3.13) (or equivalently, (3.15) in the case of planar entangling surface) can be written as follows where ... λ 0 denotes the correlator in the unperturbed theory with coupling λ 0 . Recalling now that by definition K λ 0 is sitting at the t = 0 slice, we deduce that the above correlator can be easily obtained from its Euclidean counterpart by analytic continuation. One simply evaluates the Euclidean correlator K λ 0 O(t ′ E ) λ 0 and then substitutes Euclidean time, t ′ E , with i(t ′ + iǫ). Moreover, if we implement conformal perturbation theory and expand this correlator around λ 0 = 0, we will get [23] where we have used the fact that the correlator K 0 O 0 vanishes for the CFT 10 . This observation implies that the result for S (1) (t) to linear order in λ 0 will be independent of the underlying CFT. 11 This follows from the universality of conformal correlators OO and T µν OO [44,45]. Thus we conclude that all terms of S (1) (t) which are proportional to λ 0 δλ(t) (or λ 0 ∂ 2 t δλ(t) and so on) are identical for any perturbed CFT. In particular, this explains the perfect matching, mentioned in the previous section, between certain terms in the holographic and free field theory calculations. This can be seen in the following way: Note that the terms proportional to λ 0 δλ(t) in S (1) (t) complete into terms which are quadratic in λ(t) or its derivatives in the full time dependent entropy, S(t). Recall from the section-2 that the terms that are quadratic in λ(t) or its derivatives (eg. λ 2 (t) or λ(t)λ(t)) in the holographic entanglement entropy only multiply the effective Newton constant, G d+1 , whereas λ 3 (t) contributions also multiply κ and λ 4 (t) also multiply ω, where κ and ω are coupling constants in (2.2). Since the effective Newton constant is determined by the boundary two point function according to (2.4), the coefficients of terms quadratic in λ(t) or its derivatives in the holographic calculation are therefore universal. See [23,24] where this universality was observed in a static setup.
By assumption there is only one relevant coupling, g, in the action. Hence, it can be shown that [40,46,47] where c (0) (µ) is the spectral function defined by the 2-point function of the energy-momentum trace [32], K 0 on the right-hand side is the modified Bessel function of the second kind 12 , β is the anomalous dimension of g and Ω d = 2π d/2 /Γ(d/2).
Like in the free scalar case we first integrate this correlator over x and x ⊥ , (3.30) Next, we analytically continue this result back to Lorentzian signature t E → i(t ′ + iǫ), substitute it into (3.27) and get the following spectral representation of the general formula (3.15) where we have absorbed the time dependence in the function Λ(t, µ) defined as The divergences of (3.31) depend on the µ → ∞ behavior of c (0) (µ). Expanding c (0) (µ) near µ = ∞ results in a finite number of terms of the form µ s for s > 0, and upon substituting these terms into (3.31) we find integrals of the form Note that the oscillatory nature of the integrand means that integrating over µ alone is not enough to give a divergence: integrating over µ would only produce a factor ∝ |t ′ | −s−2 , but the subsequent integral over t ′ could have a divergence near t ′ = 0, depending on the behavior of λ(t + t ′ ).
To isolate the divergent terms, we begin by integrating (3.32) by parts with respect to t ′ a total of 2N + 1 times. This results in the identity If 2N > s + 1, then the last line produces no UV divergences upon substitution back into (3.33). The oscillatory terms in the middle two lines likewise produce no divergences in (3.33). All that remains is the first line. There we see clearly that a logarithmic divergences is produced if and only if s is an odd integer, with the coefficient of the divergence proportional to λ (s+1) (t). The non-universal coefficients can also be determined this way. Thus we learn that the quasi-adiabatic expansion holds for the divergences in the entanglement entropy in a general interacting theory to this order in perturbations, and that the coefficients are entirely determined by the large-µ behavior of the spectral function c (0) (µ).

Example: free Dirac field
Let us implement (3.31) in the case of free massive Dirac field with a time-dependent mass deformation. In this case ∆ = d − 1, β = 0, and the spectral function is given by [32] c (0) Or equivalently, (3.36) This representation is useful for analyzing the UV divergences, which are built from the modes satisfying µ/(2m) ≫ 1.
Thus, for the free Dirac fermion we get S (1) where we introduced a dimensionless integration variable x = µ/(2m). From the above expression, it is obvious that only terms satisfying d − 4 − 2n ≥ 0 result in UV divergences, whereas for d − 4 − 2n < 0 integrals over t ′ and x are finite. In fact, the universal divergences can be readily isolated. We assume that n is sufficiently large and carry out integral over x. Small values of n are treated by analytic continuation. Hence, In particular, there are no universal divergences for odd d, 13 while for even d the above formula reads S (1) This result is in full agreement with known calculations in the literature [37][38][39]. For the time-dependent case, we find Note that if we take the holographic result in d = 4, (2.44), and set m(t) = m + λ(t), then we find agreement with (3.41). Similarly, the double derivative term in d = 6 in (2.51) matches with that in (3.42). The precise match between the calculations for the free field theory and strongly coupled holographic CFT has to do with universality of certain correlators that we discussed in the previous section, see also [23,24].

Conclusions
In this work we studied the evolution of the 'area law' of the spatial entanglement when the relevant coupling of the field theory undergoes a quantum quench. Our main results, (3.13) and (3.15), for the time dependent first order correction to the area law were derived under the assumption that the relevant coupling satisfies δλ(t)/λ 0 ≪ 1, where λ 0 is the unperturbed value and δλ(t) is a time-dependent perturbation which vanishes in the past and asymptotes to a constant value in the future. We argued that if λ 0 is the only relevant parameter in the system, then further simplification takes place and the final expression (3.31) for the time dependent correction to the area law can be expressed in terms of the spectral function, c (0) (µ), which is defined in terms of the two-point correlation function of the trace of the energy-momentum tensor in the theory.
As an application of these general formulas, we scrutinized the time dependence of the area law in a generic QFT. In particular, we illustrated that 'area law' divergences reveal a quasi-adiabatic behavior, i.e., they are local and depend on the instantaneous value of the relevant coupling and finitely many derivatives. We showed that this qualitative behavior is in full agreement with the holographic prediction based on the HRT proposal [27]. 13 Recall that non-universal divergences vanish within our choice of regularization scheme.
We also carried out explicit calculations of the universal terms in the case of free fields perturbed by a time-dependent mass, and found quantitative agreement with the HRT predictions of the coefficient of λ 2 (t) in the entanglement entropy. 14 From the point of view of the holographic calculation, we observed in section-2 that the coefficients of the terms which are quadratic in λ(t) or its derivatives ( eg. λ 2 (t) or λ(t)λ(t)) in the entanglement entropy only multiply the effective Newton's constant, which is fixed by the boundary two point functions. Since CFT two point functions have a universal form, this ensures that coefficients of these terms are independent of the exact nature of the unperturbed CFT. From the purely field-theoretic point of view, we observed in equation (3.28) that the coefficient of λ 0 δλ(t) in S (1) (t) is fixed by the universal CFT correlation functions T µν OO and OO . As the term proportional to λ 0 δλ(t) in S (1) (t) completes into the terms quadratic in λ(t) and its derivative, this explains the matching of the coefficients of these quadratic terms between free field theory and holographic calculations.
Up to a set of mild assumptions, the RT proposal [26] for holographic entanglement entropy in the static case was derived in [48]; see also [49,50]. The time-dependent case is more complicated: we view the match between the time-dependent terms in holography and QFT as non-trivial evidence in favor of the time-dependent HRT proposal for holographic entanglement entropy. 15 It would be interesting to extend our analysis to the next-to-leading order correction in δλ(t)/λ 0 ≪ 1. Such an extension might shed light on the fundamentals of the linear growth of entanglement entropy sufficiently far from the instant of quench [1][2][3][4][5][6][7][8][9][10], and help to uncover the underlying reason for the universality of this behaviour. To make a step in this direction, it is necessary to establish a consistent perturbative expansion for the modular Hamiltonian. In fact, this problem is interesting by itself since it has many valuable practical applications beyond analysis of the time-dependent program studied in this work, and we plan to address it in the future.

A.1 HRT surface calculation
In this appendix, we will solve for the position of the HRT surface and verify the claim that we made in (2.19) that at the lowest order, the solution of the stationary surface is As mentioned in the Section 2.2, the HRT surface is given by t = t(z) and x 1 = 0. The equation of motion for t = t(z) is obtained by varying We find where dot represents derivative with respect to t and prime represents derivative with respect to z.
The ansatz (2.12)-(2.13) for the metric functions, h(z, t) and f (z, t), at leading order gives Lets make the ansatz t(z) = t 0 + z β t β . The first term of the equation of motion (A.3) is proportional to z 2α+1 . Out of all the other terms, the term with the smallest power of z is proportional to z β−1 . Thus we learn β = 2 + 2α.
However, the coefficient of the z β−1 term is equal to 2β(d − β)t β . Therefore, the above result is only true if d = β. For d = β, we would have to modify the ansatz to be t(z) = t 0 + z β log(z)t β . We do not consider this case in the main argument of the paper, so it can be ignored.

A.2 Perturbative Solution of Einstein-scalar equations
In this section, we will find the time-dependent coefficients h m,n (t), f m,n (t) and φ m,n (t) in terms of φ 0,0 (t) = λ(t) by substituting the ansatz (2.12)-(2.14) in the Einstein-scalar equations (2.9)-(2.11). We use this result in the holographic calculation of the entanglement entropy in Sections 2.3-2.4.
(B.2) Introducing q = p 1 + p 2 , ∆t = t − t ′ , ∆x ⊥ = x ⊥ − x ′ ⊥ and integrating (B.1) with respect to x and x ′ , yields It is important that we perform the integrals over x ⊥ and x ′ ⊥ in the proper order so as to avoid an ambiguous "0 × ∞" result. First we integrate over x ⊥ , using the identity 16 and then over x ′ ⊥ which results in a 2πδ(q ⊥ ). The result is