Dynamics of logarithmic negativity and mutual information in smooth quenches

In this paper, we study how quantum correlation between subsystems changes in time by investigating time evolution of mutual information and logarithmic negativity in two protocols of mass quench. Hamiltonian in both protocols is for 2-dimensional free scalar theory with time-dependent mass: the mass in one case decreases monotonically and vanishes asymptotically (ECP), and that in the other decreases monotonically before t = 0, but increases monotonically afterward, and becomes constant asymptotically (CCP). We study the time evolution of the quantum correlations under those protocols in two different limits of the mass quench; fast limit and slow limit depending on the speed with which the mass is changed. We obtain the following two results: (1) For the ECP, we find that the time evolution of logarithmic negativity is, when the distance between the two subsystems is large enough, well-interpreted in terms of the propagation of relativistic particles created at a time determined by the limit of the quench we take. On the other hand, the evolution of mutual information in the ECP depends not only on the relativistic particles but also on slowly-moving particles. (2) For the CCP, both logarithmic negativity and mutual information oscillate in time after the quench. When the subsystems are well-separated, the oscillation of the quantum correlations in the fast limit is suppressed, and the time evolution looks similar to that under the ECP in the fast limit.


Abstract
In this paper, we study how quantum correlation between subsystems changes in time by investigating time evolution of mutual information and logarithmic negativity in two protocols of mass quench. Hamiltonian in both protocols is for 2-dimensional free scalar theory with time-dependent mass: the mass in one case decreases monotonically and vanishes asymptotically (ECP), and that in the other decreases monotonically before t = 0, but increases monotonically afterward, and becomes constant asymptotically (CCP). We study the time evolution of the quantum correlations under those protocols in two different limits of the mass quench; fast limit and slow limit depending on the speed with which the mass is changed. We obtain the following two results: (1) For the ECP, we find that the time evolution of logarithmic negativity is, when the distance between the two subsystems is large enough, well-interpreted in terms of the propagation of relativistic particles created at a time determined by the limit of the quench we take. On the other hand, the evolution of mutual information in the ECP depends not only on the relativistic particles but also on slowly-moving particles. (2) For the CCP, both logarithmic negativity and mutual information oscillate in time after the quench. When the subsystems are well-separated, the oscillation of the quantum correlations in the fast limit is suppressed, and the time evolution looks similar to that under the ECP in the fast limit.

Introduction and summary
Understanding how quantum states thermalize is one of research themes on the frontier of theoretical physics. Thermodynamics well-describes properties of systems with macroscopic degrees of freedom, and thermodynamic systems time-evolve irreversibly. On the other hand, quantum field theories (QFTs) well-describes properties of systems with microscopic degrees of freedom. The states in QFTs time-evolve unitarily, and a pure state never thermalizes. Therefore, many researchers have been studying how thermalization is interpreted in terms of microscopic degrees of freedom. Since the area of black hole is subject to thermodynamic laws, the black hole is a thermal object [1,2,3,4,5,6]. In terms of holography, thermalization is closely related to the dynamics of black hole such as its formation and evaporation [7,8]. Therefore, we expect a microscopic description of thermalization to help us to study (quantum) dynamics of black holes. Thus, understanding thermalization in terms of microscopic degrees of freedom is a challenging but worthwhile problem.
Thermalization of a whole system of QFTs is unknown, but that of a subsystem can be studied by looking at its dynamics. When we focus on the local dynamics in the subsystem A owing to quantum entanglement of a state, the degrees of freedom in the outside of A play as a thermal bath. Therefore, late-time expectation values of local observables, n-point functions and entanglement entropy, will be approximated by those of a thermal state. Authors in [9] have studied time evolution of entanglement entropy in a 2-dimensional dynamical system, where a gapped theory changes suddenly to a CFT (a global quench). They have found that the lattice-independent term in a certain scaling limit increases linearly in time and is proportional to the subsystem size at late time. Authors in [10] have numerically studied the time evolution of entanglement entropy for a global quenched state in holographic theories, which have gravity duals. The authors in [11,12,13,14] have proposed holographic models, which analytically describe its linear growth in time in holographic theories. In those articles, researchers have been studying the evolution of entanglement entropy when Hamiltonian changes globally and suddenly.
Also, authors in [15] have studied the time evolution of entanglement entropy for excited states in local quenches, where a Hamiltonian changes locally and suddenly. Since the late-time entropy for locally-quenched states is not approximated by thermal entropy, the subsystem in the local quenches does not thermalize. Various holographic models describing local quenched states are proposed [16,17,18].
The time evolution of entanglement entropy for the excited states by the quenches shows thermalization depends on detail of quenches. Since three of authors in this paper has been interested in the time evolution of entanglement entropy in more general setups, we have studied the dynamics of quantum entanglement when a time-dependent Hamiltonian changes smoothly in time [19,20]. The time-dependent Hamiltonian is for a 2-dimensional free scalar theory with two kinds of time-dependent mass [21] (See Figure 1). If the time-dependent mass is initially constant, it starts to decrease around t = 0, and it asymptotically vanishes, then the quench is called End-Critical-Protocol (ECP). If the time-dependent mass is initially constant, it starts to decrease before t = 0, it vanishes at t = 0, it starts to increase after t = 0, and it asymptotically becomes constant, then the protocol is Cis-Critical-Protocol (CCP). In those smooth protocols, scaling behaviors of observables in QFTs are studied, and their time-evolution is compared to the one in the sudden quenches [22,23]. The two protocols in this paper have two parameters, which are an initial correlation lengthξ and a quench rate δt. We have studied the time evolution of entanglement entropy in the protocols in two limits, δt ξ 1 and δt ξ 1, which are first considered in [24]. The time evolution of entanglement entropy in both limits for ECP shows that the late-time entropy is proportional to the subsystem size because quasi-particles are created at a characteristic time, which is determined by the limit. On the other hand, the entropy in both limits for CCP oscillates in time, and we have found the entropy for the quenched state in the fast limit δt ξ 1 is proportional to the subsystem size, but we could not find whether the entropy in the slow limit δt ξ 1 is proportional to the size. In this paper, we study the time evolution of mutual information I A,B and logarithmic negativity E in order to take a more closely look at how the dynamics of quantum entanglement depends on these protocols. These quantum quantities, I A,B and E, measure quantum correlation between subsystems A and B with several distances between them, and those in 2D CFTs depend on the detail of CFTs because they depend on n-point functions where n > 2.

Summary
We study the dynamics of quantum entanglement by measuring the time evolution of logarithmic negativity E and mutual information I A,B in two protocols ECP and CCP. We define two dimensionless parameters ξ and δt by ξ =ξ a and δt = δt a . Here, a is a lattice spacing. The summary of main results is: ECP 1. The evolution of E and I A,B is interpreted in terms of the relativistic propagation of quasi-particles. In the fast limit where δt/ξ 1, the quasi-particles are created around t = 0. In the slow limit where δt/ξ 1, the particles are created at t ∼ t kz = δt, when adiabaticity breaks down.
2. As we expect from quasi-particle picture, slopes of E and I A,B for l a , l b ξ at early time are independent of the subsystem sizes of A and B when the distance between them is much larger than typical scales, which are determined by these limits.
3. The late-time I A,B is affected by modes propagating slowly for the disjoint interval where the distance between A and B is much larger than the characteristic length. However, the quasi-particles which slowly propagates does not affect E for the disjoint interval.

CCP
4. If the distance d between A and B is smaller than ξ, then I A,B and E in the fast CCP oscillate in time. Their late-time period is given by πξ. The amplitude of oscillation for the large distance d ξ is smaller than that for d ξ. The time evolution of I A,B and E in this limit is similar to that in the fast ECP. The evolution of I A,B and E for the subsystems with d ξ is expected to be independent of the late-time protocol qualitatively.
5. The measures I A,B and E in the slow CCP oscillate with the period, which approaches πξ at late time. As the distance d between the subsystems increases, the time, when I A,B and E in the slow CCP start to increase, becomes later.

Organization
In section 1, we have explained why we study the time evolution of logarithmic negativity and mutual information, and we have summarized what we find in this project. In section 2, we will define E and I A,B , and we will explain how to compute them in quantum field theories. In section 3, we will explain the details of two protocols of quench where we study the time evolution of E and I A,B . In the following section, we will explain what we obtain in this project. Finally, we would like to discuss our results and comments on a few of future directions.

Logarithmic negativity and mutual information
In this section, we define quantum measures for states, logarithmic negativity and mutual information, and explain their properties briefly. After that, we explain how to compute them in path-integral formalism and in terms of correlation function, replica trick and correlator method. Roughly speaking, the quantum measures in the replica trick are given by a kind of free energy in a partition function on a replicated geometry [25,26]. The correlator method is applicable to computation of quantum measures for a gaussian state in free field theories. It is a powerful tool to compute them numerically [27,28,29,30,31,32].

Definition of logarithmic negativity and mutual information
When the total space is A ∪ B ∪ C, logarithmic negativity and mutual information E and I A,B can measure non-local correlation which comes from quantum entanglement between subsystems, A and B, as defined in Figure 2. Let us explain the definition of E and I A,B , for which Hilbert space on a time slice is H = H A + H B + H C . Here, subsystem sizes of A and B are l a and l b , and the distance between A and B is d. For simplicity, we consider (1 + 1) dimensional quantum field theories.

Definition of mutual information between A and B
Mutual information I A,B is defined by a combination of von Neumann entropies for reduced density matrix ρ α , where ρ α is the reduced density matrix which is defined by tracing out the degrees of freedom outside α: Entanglement entropy S α is given by Rényi entanglement entropy S (n) α in von Neumann limit, Mutual information is non-negative quantity and has a upper bound due to Araki-Lieb inequality S A∪B ≥ |S A − S B | [33], Tripartite information, a linear combination of mutual information, is expected to be a useful tool to characterize conformal field theories (CFTs), which has gravity dual [34]. The tripartite information in holographic CFTs is a non-positive quantity, and the information in the strongly chaotic systems captures the scrambling properties of holographic CFT [35,36].

Definition of logarithmic negativity
Both I A,B and E can capture the non-local correlation between A and B. However, only E can measure "pure" quantum entanglement [37,38]. The quantum measure, which can measure pure quantum entanglement, should be an "entanglement monotone". That is, it should decease monotonically under a local manipulation which is called local measurement and classical communication (LOCC) [39,40,41]. Entanglement entropy for a pure state is an entanglement monotone because the entropy decreases monotonically under LOCC. However, the entropy for a mixed state is not. Logarithmic negativity is an entanglement monotone even for a mixed state. Let us explain how to define E. The definition of logarithmic negativity E is E ≡ log ||ρ T B A∪B || = log tr|ρ T B A∪B |, (2.5) where ||ρ T B A∪B || ≡ i |λ i | is the trace norm of ρ T B A∪B , and λ i are eigenvalues of ρ T B A∪B . The operation · T B is called by partial transpose. This operation transposes the matrix component as follows, Thus, the partially-transposed reduced density matrix ρ T B A∪B is given by (2.7)

Replica trick
We review the replica trick for logarithmic negativity E based on Refs. [42,43,30]. Instead of tr|ρ T B A∪B |, consider tr(ρ T B A∪B ) n for an integer n. The λ i -dependence of tr(ρ T B A∪B ) n for even n e is different from that for odd n o as follows, Since the trace norm tr|ρ T B A∪B | is given by taking the analytic continuation of the even sequence at n e → 1, the negativity E in the replica trick is given by We are able to compute E in 2D CFT by using twist fields. For example, consider the entire space as A ∪ B and A = [u 1 , u 2 ]. In this system, E can be expressed as where T 2 n is the twist field that connects the j-th to (j + 2)-th replica fields. For this configuration, E corresponds to Rényi entanglement entropy with n = 1/2 [44]. As another example, consider two disjoint intervals A = [u 1 , u 2 ] and B = [u 3 , u 4 ]. Logarithmic negativity E for this configuration is where T n is the twist field that connects the j-th sheet to (j + 1)-th one. It is well known that the conformal dimensions of T n andT n in 2D CFT with central charge c are [25,26] ∆ n = c 12 n − 1 n . (2.13) Because of the partial transpose T B , the connection between the replica sheets at B is different from that at A as shown in Figure 3. This partial transpose T B has the effect to exchange the twist operators at u 3 and u 4 compared to the ones in the correlation function We consider the time evolution of E between the two disjoint intervals in 2D CFT with central charge c in a sudden quench, where a mass gap suddenly vanishes. After the gap vanishes, the ground state for an initial Hamiltonian time-evolves under a CFT Hamiltonian. When we compute E in the replica trick, E in Euclidean space is given by the 4-point function of the twist fields on the boundaries of two strips. Real time-evolution of E is given by performing an analytic continuation τ = τ 0 + it, where τ is the imaginary time of the twist fields, and t is the real time.
We introduce a parameter τ 0 as width of the two strips and interpret 1/τ 0 as the energy scale of the initial ground state of the massive theory. In the limit of t τ 0 and |u j −u i | τ 0 , (2.14) where we assumed |u 2 − u 1 | < |u 4 − u 3 |. We note that mutual information I A,B shows the same behavior up to normalization (for example, see [30]).

Setup and correlator method
Here, we explain the detail of the dynamical system, a 2D free scalar theory with time dependent masses. We also explain how to compute I A,B and E by correlation functions.

Setup
The excited state in our setup time-evolves with a time-dependent Hamiltonian from the ground state for the initial Hamiltonian. In order to study the time evolution of quantum correlation between subsystem A and B, we define the change of mutual information and logarithmic negativity ∆I A,B and ∆E by where I A,B (t) and E(t) are measured at t, and I A,B (t = −∞) and E(t = −∞) are measured at t = −∞. We compute ∆E and ∆I A,B with the time-dependent Hamiltonian, which is for (1 + 1) dimensional free boson theory with time-dependent massm(t), where the time-dependent massm(t) is dimensionful and changes continuously in time. The time-dependent massm(t) in early time limit t → −∞ approaches a constant mass asymptotically. Therefore, we approximate I A,B (t = −∞) and E(t = −∞) by mutual information and logarithmic negativity for the ground state of the free scalar theory with a static mass m.
We compute I A,B (t) and E(t) by a method, correlator method, where these quantities for gaussian states are given by the summations of eigenvalue-dependent entropies and logarithmic negativities s A (γ) and (χ) [27,28,29,30,31,32]. The eigenvalues γ and χ are for matrices whose components are correlators in the subregions. Since the dimension of these matrices depends on the degrees of freedom in the subsystems, it becomes finite in a discretized system. After discretizing the system, we are able to diagonalize the correlator matrices, and we are able to compute ∆I A,B and ∆E. We discretize the system by putting it on a circle with the circumference L = N and the lattice spacing . After replacing dx → N −1 n=0 and φ → q n , π → p n / ,m(t) → m(t)/ as well as ∇φ → (q n+1 − q n )/ , the discrete Hamiltonian on the circle is, where we assume that N = 2M + 1 (M ∈ Z). We choose as a unit and define the dimensionless parameters t, l, δt and ξ wheret,l, δt andξ are dimensionful. The variables q n and p n are subject to a periodic boundary condition: Also, they have a canonical commutation relation: Discrete Fourier transforms are given by, whereq † l =q −l andp † l =p −l because p κ and q κ are real. The Hamiltonian in the momentum space is where [a k , a l ] = a † k , a † l = 0 and a k , a † l = δ k,l . The function f κ (t) has a property f κ (t) = f −κ (t). Then, f κ (t) satisfies the following Wronskian condition: After taking a thermodynamic limit where N, L → ∞, and is finite, f κ (t) is determined by the following equation of motion: The coordinates and its conjugate momenta in the thermodynamic limit are given by wherep andq are defined in (2.23).

Correlator method
Here, let us explain the correlator method where entanglement entropy and logarithmic negativity are related to eigenvalues of matrices whose components are given by correlation functions in the subsystem (see, for example, [32] for entanglement entropy and [27] for logarithmic negativity). If the given subsystem is V, we can define correlation functions X ab (t), P ab (t) and D ab (t) by where local operators, X a (t) and P a (t), are in V. If the size of V is L, we can define two 2L × 2L matrices, Γ and J, by Entanglement entropy S A is given by the summation of eigenvalue-dependent entropies s A (γ κ ), where γ κ are positive eigenvalues of a 2L × 2L matrix M = iJΓ.
In the correlator method, a given partial transposed reduced density matrix ρ T B A∪B is where χ a are positive eigenvalues of iJ · Γ T B . The index a runs from 1 to L. Therefore, the logarithmic negativity E is given by (2.32)

Smooth quench
Here, we explain the protocols which trigger the time evolution of ∆I A,B and ∆E. The protocols in this paper have the time-dependent masses, one of which is end-critical-protocol (ECP), and the other is Cis-critical-protocol(CCP).

End-critical-protocol (ECP) and Cis-critical-protocol (CCP)
The time-dependent masses in this paper depend on two parameters ξ and δt. Here, ξ and δt are an initial correlation length and quench rate. The time-dependent mass m(t) in ECP is, which is constant approximately at early time t −δt, decreases after t = −δt, and vanishes at late time δt t. The mass in CCP is, which is constant approximately in the time region t −δt, decreases around t = −δt, vanishes around t = 0, increases after t = 0, and approaches a constant values asymptotically. The dynamics of quantum entanglement is governed by the Hamiltonian with the timedependent masses in (3.1) and (3.2).
As in [45], the function f κ (t) for the time-dependent masses in (3.1) and (3.2) is analytically computable because the equation of motion for f κ (t) is given by The solution of equation of motion (3.3) for the mass in ECP f κ (t) is given by where the parameters A, B, α and β are given by (3.5) On the other hand, the solution of equation of motion with the mass for CCP is given by where parameters a, b, C, D, and α are Since the protocols in (3.1) and (3.2) have two parameters ξ and δt, the protocols change slowly or fast by changing the parameters. We study the time evolution of ∆E and ∆I A,B in two limits, which will be explained later.

Kibble-Zurek mechanism
Let us explain criteria for the adiabaticity. Suppose an ansatz of in (3.9) is much smaller than ω 2 κ (t), the adiabatic approximation W κ (t) ∼ ω κ (t) is valid. Then, the condition where we can use the adiabatic approximation is (see, for example, [46]) We focus on the LHS at κ = 0 in (3.10) because it is the largest for any κ. We also assume 2 . Thus, we obtain criteria for the adiabaticity If C L (t) satisfies the condition (3.11), we can examine the time evolution of systems by using the adiabatic approximation. The Kibble-Zurek time t kz is defined by and it is the time scale at which the adiabaticity breaks down.
In this paper, we consider two limits of the quenches, one of which is a fast limit: and ξ is the initial correlation length, and the other is a slow limit: 14) The adiabaticity is broken at late time in the slow ECP limit and is broken only near t = 0 in the slow CCP limit. In the slow ECP and CCP, the Kibble-Zurek time t kz is (see, for example, [20]) where we also defined an length scale ξ kz by the inverse mass at the Kibble-Zurek time. In the slow quenches, t kz and ξ kz are the effective time and length scales, respectively. We also define an energy scale E kz = 1/ξ kz in the slow ECP limit. The Kibble-Zurek mechanism [47,48] conjectures that the dynamics in the slow limit is frozen when the adiabaticity is broken. For example, the effective correlation length in the slow CCP, which is determined at t −t kz , is almost time-independent from t −t kz to t t kz .

Time evolution of logarithmic negativity and mutual information 4.1 Time evolution in ECP
In the beginning, we calculate ∆E and ∆I A,B numerically in the ECP 5 . Here, we consider the ECP in two limits: the fast limit and the slow limit. We study how ∆E and ∆I A,B time-evolve in these limits. The similar setup to the fast limit was studied in [30], where the authors have studied time evolution of ∆I A,B and ∆E in a global quench. The quantities ∆I A,B and ∆E in the fast limit change in time in the similar manner to those in the global quench.
In section 4.1.1, we show the time evolution of ∆E and ∆I A,B in the fast limit, which are computed numerically, and how their time-evolution depends on ξ and δt. After that, we interpret their time-evolution in terms of our toy model where relativistic propagation of local objects contributes to quantum correlation between the subsystems. In section 4.1.2, we show the time evolution of ∆E and ∆I A,B in the slow limit and how they depend on the parameters. Again, we interpret their evolution in terms of the relativistic propagation of "particles" created at the Kibble-Zurek time t kz .

Fast limit
The time evolution of ∆E(t) and ∆I A,B (t) is shown in Figure 4. Without loss of generality, we can assume l b ≥ l a . The panels (a) and (b) of Figure 4 show how ∆E(t) and ∆I A,B (t) for the protocol in (3.1) with (ξ, δt) = (100, 5) depend on the distance between two subsystems whose sizes are the same. The panels (c) and (d) of Figure 4 are similar to the panels (a) and (b), but the subsystem sizes are different. Here, we take the subsystem size l b to be twice as large as the subsystem size l a , l b = 2l a . We plot ∆E(t) and ∆I A,B (t) with d = 0, 10, 100, 500, 800, and 1500 in the panels (a), (b), (c), and (d) of Figure 4. The panels (e) and (f) of Figure 4 show how ∆E and ∆I A,B depend on the subsystem sizes with fixed distance. We show those for d = 100. In the panels (g) and (h) of Figure 4, we plot ∆E and ∆I A,B for ECP with two pair of parameters, (ξ, δt) = (200, 10) and (ξ, δt) = (100, 5), to see the scaling law.
We summarize properties of ∆E and ∆I A,B with 2l a = l b . Note that ∆E and ∆I A,B with l a = l b have properties similar to those with 2l a = l b . 5 In order to compute the integrals in (2.27), we replace π −π dκg(κ) with a summation of g(κ i ) numerically, where κ i = −π + εj, j is an integer, and ε is a small parameter for our numerical computations. We compute ∆E and ∆I A,B in the ECP with ε = 0.0001, 0.00001, 0.000001. We find no significant difference between them with the different values of ε.
(1) Time t s when ∆E and ∆I A,B begin to increase depends on d: It is independent of the subsystem sizes.
(2) The quantities ∆E and ∆I A,B monotonically increase in the window t s t t s + la 2 . In this window, ∆E for l a , l b ξ is fitted by a following linear function in t, F f (t, ξ), and ∆I A,B for l a , l b ξ is fitted by G f (t, ξ): Since the slopes of ∆E and ∆I A,B in time do not depend on l a and l b but depend on d as one can see in the panels (e) and (f), it is enough to evaluate the coefficients for various d with a pair of subsystem sizes (l a , l b ). Here we provide the numerical results of a 1 and a 2 with (l a , l b ) = (1500, 1500) in Table 1. By a computation in another protocol wiith (ξ, δt) = (200, 5), we find that a 1 appears to be independent of ξ and δt for d = 0, and its value is a 1 ∼ 0.5.
(3) The slopes of ∆E and ∆I A,B change at t = t M , which depends on l a as follows: (4) We define width of the plateau w as a time interval between t M and a time when the slopes of ∆E and ∆I A,B change from a positive or zero to a negative. Then, w for d ξ is given by w ∼ l b −la 2 , and its magnitude is independent of d for d ξ.
(5) The quantities ∆E and ∆I A,B in the window t s + la 2 + w t t s + l a + w monotonically decrease.
(6) The quantities ∆E and ∆I A,B for l a , l b , d, ξ 1 obey the scaling law, (MI1) The mutual information ∆I A,B for d ξ decreases monotonically in the window t M < t < la+l b 2 , and increases logarithmically after t ∼ la+l b 2 . On the other hand, ∆I A,B for d ξ decreases monotonically in the window t M < t < l a + d 2 , and increases after t ∼ la+l b +d 2 as ∆I A,B for d ξ does.
(MI2) The magnitude of ∆I A,B in the window, t la+l b +d 2 , depends on the parameters, l a , l b and d. However, how ∆I A,B time-evolves appears to be independent of the parameters. That is, we can fit to ∆I A,B logarithmic function t, 1/2 log t + c, where c is a constant which depends on the subsystem size.
Interestingly, ∆E for the long distance is perfectly consistent with the analytic computation in CFT; around t ∼ d, ∆E increases linearly, and reaches the maximum value at t ∼ d + la 2 . The maximum value continues to t ∼ d + la 2 + w. After that, ∆E behaves as linear function with negative slope. Finally, ∆E vanishes in the window t l b + d. These behavior agrees with [30].

Slow limit
Next we show some plots of slow ECP in Figure 5.  (1) Compared with the fast limit t s , which is defined in the property (1) of the fast limit, is shifted by the Kibble-Zurek time t kz for d For the slow ECP with (ξ, δt) = (2, 200), the Kibble-Zurek time defined in (3.15) t kz is t kz ∼ 921. Again, t s is independent of the subsystem sizes.
(2) The quantities ∆E and ∆I A,B in the window t s t t s + la 2 are fitted by following linear functions in t, F s (t, E kz ) and G s (t, E kz ); where a 1 and b 1 are the coefficients of t · E kz , and a 2 and b 2 are time-independent terms. The slopes of F s and G s are independent of the subsystem sizes as one can see in Figure  5, so that it is enough to evaluate a 1 and b 1 for various d with a pair of the subsystem sizes (l a , l b ). The fitting results for d = 0, 10, 100, 500, 800, 1500 with l a = l b = 600 are summarized in Table 3. The coefficients a 1 and b 1 appear to be independent of d for d 1 E kz , and a 1 is independent of ξ and δt for d = 0, as well as the property (2) in the fast limit. Its value is a 1 ∼ 0.25.
(4) The width of the plateau w in the slow limit is the same as the fast limit, w ∼ l b −la 2 for d 1 E kz , and the local maximum value in time appears to be independent of the distance.
(5) The quantities ∆E and ∆I A,B in the window t s + la 2 + w t t s + l a + w monotonically decrease.
(6) The quantities ∆E and ∆I A,B for l a , l b , d, ξ 1 obey the following scaling law, The main difference in time scale compared with the fast limit is that t s and t M are shifted by t kz in the slow limit. We will interpret the difference as our toy model in section 4.1.3.

Quasi-particle interpretation
The mass potential in the ECP changes as follows: The mass term in the very early time t −δt is constant, deceases monotonically after t ∼ −δt, and vanishes asymptotically at late time t δt. Entanglement structures of states in the early time therefore changes adiabatically.
In the fast limit, ∆E and ∆I A,B monotonically increases from t s ∼ d/2 for d ξ. On the other hand, in the slow limit, t s is shifted by t kz , t s ∼ d/2 + t kz for d 1 E kz . This means that the entangled pairs, which are composed of right-and left-moving particles p R and p L , are created everywhere when the adiabaticity is broken. The left-moving particle p L is entangled with p R . Entanglement between p L and p R contributes to ∆E and ∆I A,B only when p L is in A (or B), and p R is in B (or A). When p L or p R go away from A or B, ∆E Here we take the distance between A and B to d = 0, 10, 100, 500, 800, and 1500. The panels (e) and (f) show the time evolution of ∆E and ∆I A,B for various subsystem sizes with d = 100. In the bottom panels (g) and (h), we plot ∆E and ∆I A,B with two protocols: one is (ξ, δt) = (5, 500) with (l a , l b , d) = (1000, 2000, 800), the other is (ξ, δt) = (2, 200) with (l a , l b , d) = (400, 800, 320).  200) (600, 600, 500) 6 ≤ t · E kz ≤ 7 0.0225 t · E kz − 0.125 (2, 200) (600, 600, 800) 13 2 ≤ t · E kz ≤ 15 2 0.0152 t · E kz − 0.0949 (2, 200) (600, 600, 1500) 17 2 ≤ t · E kz ≤ 9 0.0117 t · E kz − 0.0938 (4, 800) (600, 600, 0) 0.245 t · E kz + 12.4 (4, 800) (600, 600, 10) 0.132 t · E kz − 205 (4, 800) (600, 600, 100) 0.0584 t · E kz − 159 (4, 800) (600, 600, 500) 5 ≤ t · E kz ≤ 23 4 0.00723 t · E kz − 25.2 (4, 800) (600, 600, 800) 41 8 ≤ t · E kz ≤ 6 0.00190 t · E kz − 7.20 (4, 800) (600, 600, 1500) 23 4 ≤ t · E kz ≤ 13 2 0.00108 t · E kz − 0.463 and ∆I A,B decrease. The entangled pair created on the middle of two regions reaches A and B at t = t s ∼ d/2, and goes away from there at t = t M ∼ t s + l a /2. When the subsystem sizes are the same l a = l b , ∆E and ∆I A,B decrease after t M . Otherwise, ∆E shows a plateau from t M to t M + l b /2 − l a /2 since the inflow and outflow of the entangled pairs into A and B balance each other. The similar behavior is observed in ∆I A,B . In our computation, ∆E does not vanish in the late time t l b +d 2 for d ξ in the fast limit or d 1/E kz in the slow limit, whereas ∆E vanishes in 2d CFT analysis. The mutual information ∆I A,B in our computation increases logarithmically in the late time for arbitrary l a , l b , d, ξ, δt. However, this vanishes in the 2d CFT.
In the 2d CFT, only particle pairs with speed of light are created, and contribute to ∆E and ∆I A,B . In our setup, there are entangled particles propagating with various velocities less than or equal to the speed of light as mentioned in the previous work [20]. The lattice artifact also contributes to them since we calculate on the lattice. Therefore, non-zero values of ∆E and ∆I A,B in the late time are due to the entangled pairs with the velocities less than the speed of light and the lattice artifact.
The quantity ∆E after t l b +d 2 vanishes for d ξ in the fast limit or d 1/E kz in the slow limit, which implies that the slow modes do not contribute to ∆E for these configurations, and the lattice artifact is suppressed. As a result, the quantities in this configuration behaves in similar manner to those in the 2d CFT.
Since the late-time entanglement entropy grows as 1/2 log t in the smooth quenches [20] and ∆I A,B is given by the summation of entanglement entropies for A, B, and A ∪ B as defined in (2.1), ∆I A,B at the late time increases as (see also Figure 6) ∆I A,B ∼ 1/2 log t + 1/2 log t − 1/2 log t = 1/2 log t.

CCP
Here, we study the time evolution of ∆E and ∆I A,B in the fast and slow CCPs numerically 6 . There are several papers about the dynamics of quantum entanglement in the CCPs including, for example, the scaling property of entanglement entropy at t = 0 [19], the time evolution of entanglement entropy [20], and complexity [49]. The quantities ∆E and ∆I A,B in the CCPs oscillate in time. We will see that the amplitude of the oscillation in the fast CCP depends on the inter-subsystems distance d, and the behavior of ∆E and ∆I A,B with large d in the fast CCP is similar to that of ∆E and ∆I A,B in the fast ECP. We will also discuss the quasi-particle interpretation of the behavior of ∆E and ∆I A,B in the CCPs.

Fast limit
First, we study the time evolution of ∆E and ∆I A,B at the fast limit in the CCP. They oscillate in time. We compute ∆E and ∆I A,B in the fast CCP with the parameters (ξ = 100, δt = 5), (ξ = 100, δt = 10), and (ξ = 200, δt = 10). Here, we explain the time evolution of ∆I A,B and ∆E in the fast CCP with (ξ = 100, δt = 5), which are computed numerically. Figure 7 shows the time-dependence of the logarithmic negativity ∆E and the mutual information ∆I A,B in the fast CCP with l a = l b and 2l a = l b . We find the following properties in the fast CCP limit: (1) The quantities ∆E and ∆I A,B in the fast CCP oscillate. 6 In the CCP, we compute numerically ∆E and ∆I A,B with ε = 0.00001, 0.000001, where ε is defined in footnote 5. In our numerical computations, ∆E and ∆I A,B do not depend on ε so much. (2) If d ξ, the amplitude of oscillation of ∆E and ∆I A,B in the fast CCP seems to be small, and their time evolution becomes similar to that in the fast ECP. The quantities ∆E and ∆I A,B with d ξ time-evolve as follows: They start to increase around t ∼ d 2 , and their slopes with l a < l b become small around t ∼ la+d 2 . Then, ∆E and ∆I A,B start to decrease around t ∼ l b +d 2 , and their slopes become small after t ∼ la+l b +d (l a , l b , ξ , δt) =(600,600,100,5) of ∆E and ∆I A,B with d ξ is similar to their evolution in the fast ECP 7 . In particular, ∆E and ∆I A,B in the fast CCP with d ξ start to increase at t ∼ d 2 , and their slopes with l a < l b becomes smaller around t ∼ d+la 2 . They start to decrease at t ∼ d+l b 2 , and almost vanish after t ∼ d+la+l b 2 (the property (2)). The panels (e) and (f) show the evolution of ∆E and ∆I A,B with d = 0, and the panels (g) and (h) show them with d = 800. The panels (e) and (f) show that ∆E and ∆I A,B oscillate with the frequency, which at late time t δt is expected to be ω ∼ 1 πξ as in [20]. If we define the amplitude of oscillation of ∆E and ∆I A,B by A E and A I , the ratio of A E to |∆E(t)| is much smaller than the ratio of A I to |∆I A,B (t)|. The change ∆I A,B with d = 0 has a first dip at t ∼ 2ξ, and ∆E with d = 0 also has a small dip around t ∼ 2ξ (the property (3)). At early time t < la+d 2 , all plots in the panels (e), (f), (g), and (h) lie on the same curve, and they therefore are independent of l a and l b at the early time (the property (4)). The panel (c) of Figure 4 shows that ∆E in the fast ECP with l a < l b has a plateau in the window la+d 2 t l b +d 2 , but ∆E in the fast CCP with l a < l b in the panel (g) of Figure 7 increases in this window.
The panels (i) and (j) show the plots of ∆E and ∆I A,B with two different pairs of parameters (ξ, δt) = (100, 5) and (200,10). These panels show that the plots for ∆E and ∆I A,B in the CCP with (ξ, δt, l a , l b , d) = (100, 5, 500, 500, 500) are on the plots with (ξ, δt, l a , l b , d) = (200, 10, 1000, 1000, 1000), and we expect ∆E and ∆I A,B to have the scaling law in (4.10) if the parameters ξ, l a , l b , and d are much larger than one (the property (5)).

Slow limit
Here, we study ∆E and ∆I A,B in the slow CCP. In the slow CCP, as our previous study of entanglement entropy [20], it does not seem possible to use the quasi-particle interpretation with the speed of light for determining the time scales at which the behavior of ∆E and ∆I A,B changes. We calculate ∆E and ∆I A,B in the slow CCP with (ξ = 10, δt = 1000), (ξ = 20, δt = 2000), and (ξ = 10, δt = 4000). Here, we mainly show ∆E and ∆I A,B for (ξ = 20, δt = 2000, ξ kz = 200), which are computed numerically. Since it takes a long time to compute X κ (t), P κ (t), and D κ (t) in the slow CCP numerically, we use the adiabatic approximation of them for large |κ| as well as the computation of entanglement entropy in the slow CCP [20] 8 . See the appendix of [20] for more details of the approximation. Figure  8 shows the time evolution of ∆E and ∆I A,B in the slow CCP. We are not able to find any significant differences between ∆E and ∆I A,B . We find the following properties of ∆E and ∆I A,B in the slow CCP: (1) The time, when the time evolution of ∆E and ∆I A,B starts to increase from zero, becomes fastest-moving particles move with the maximum group velocity |v max | ∼ 1 9 . Therefore, the time evolution of ∆E and ∆I A,B in the fast CCP is well-interpreted in terms of the relative propagation of quasi-particles. Since the amplitude of their oscillation with large d is small, their time evolution is well-described in the same manner as in the fast ECP. Therefore, their late-time evolution is similar although their late-time protocols are different.
The logarithmic negativity ∆E with l a < l b and d ξ in the fast ECP shows a plateau from t ∼ la+d 2 to t ∼ l b +d 2 , and the negativity becomes zero after t ∼ la+l b +d 2 . However, ∆E for this configuration in the fast CCP increases from t ∼ la+d 2 to t ∼ l b +d 2 and is nonzero after t ∼ la+l b +d 2 (See, for example, the panel (g) of Figure 7). This is the difference between ∆E in the fast ECP and CCP with d ξ, and we expect this difference to come from the slow κ-modes of the quasi-particles.
The quantum measures at both limits in CCP ∆E and ∆I A,B oscillate in time, and their period at late time is well approximated by π m , where m is the late-time mass. Since the late-time period of their oscillation is the same as that of entanglement entropy, we expect the period to be interpreted as in [20]. The time evolution of two point function X ab , P ab and D ab determines the time evolution of ∆E and ∆I A,B . Fourier modes of X ab , P ab and D ab are given by (4.12). Each κ-mode at late time is well approximated by where C i κ , D i κ , and Θ i κ are independent of t. We expect the time evolution of ∆E and ∆I A,B at late time to be determined by the slowly-moving particles. Since the late-time velocities for quasi-particles are determined by ω κ (∞), and the particles with κ ∼ 0 propagate very slowly, then they contribute dominantly to the time evolution of ∆E and ∆I A,B at late time. At κ = 0, the period of (4.12) is π ω 0 (∞) = πξ and consistent with the period of ∆E and ∆I A,B .

Discussion and future directions
In this paper, we have studied the quantum correlation between the subsystems A and B in 2D free scalar field theory with the time-dependent mass m(t). The protocols where we have studied the correlation are ECP and CCP which depend on a pair of parameter (ξ, δt): The mass potential m 2 (t) in ECP is almost constant at early time t −δt, starts to decrease sharply around t ∼ −δt, and vanishes asymptotically in the region t δt. On the other hand, m 2 (t) in CCP is almost constant at early time t −δt, starts to decrease sharply around t ∼ −δt, vanishes at t = 0, increases monotonically after t = 0, and becomes constant asymptotically in the region t δt.
The initial state, for which we have studied the correlation A and B in the time evolution, is the ground state in the massive free scalar field theory. Since the Hamiltonian, which drives the system, is time-dependent, the state is no longer the ground state except for the region t −δt, and the correlation between A and B changes in time. We have studied how the time evolution of quantum correlation depends on (ξ, δt) in the two limits: the fast limit δt ξ and the slow limit δt ξ. We have studied the changes of mutual information and logarithmic negativity, which are defined by subtracting the information and negativity for the initial state from those for the excited state at t, in order to study how the quantum correlation between A and B changes in time. We have called the changes of mutual information and logarithmic negativity ∆I A,B and ∆E. We have found the time evolution of ∆E with large d in the fast and slow ECPs is similar to that in the sudden quench [30]. We have also found its evolution in the ECPs is well-described in terms of entangled pairs: The pairs in the fast limit are created everywhere around t ∼ 0, and that in the slow limit are created at t ∼ t kz , or the effective time when the adiabaticity is broken. The pair is composed of left-and right-moving particles with the speed of light which entangle with each other. Only when the left moving particle is in A(or B), and the right one is in B(or A), quantum entanglement between them contributes to the quantum correlation between A and B.
The relativistic propagation of entangled pairs well-describes the time evolution of ∆I A,B in the fast and slow ECPs except for its logarithmic growth in time. The time evolution of ∆I A,B is well-interpreted in terms of the "generalized quasi-particle picture": Each pair created at t ∼ 0 (t ∼ t kz ) are composed of left-and right-moving particles with various values of velocities. The change ∆I A,B at late time grows logarithmically in time, which comes from entanglement between slow-moving quasi-particles.
We have also studied the time evolution of ∆I A,B and ∆E in the fast and slow CCPs, and found that those in both limits oscillate with the frequency which, at late time, is well-approximated by m π . The longer the distance d becomes, the smaller the amplitude of oscillation of ∆I A,B and ∆E in the fast ECP do. Since the evolution of ∆E and ∆I A,B in the fast CCP for d ξ is similar to that of ∆E and ∆I A,B in the fast ECP, their evolution in the fast CCP is well-described in terms of the generalized quasi-particle picture which explains that in the fast ECP. We found the late-time period of their oscillation is interpreted in terms of quasi-particles propagating with the group velocity v κ . The late-time time evolution of ∆I A,B and ∆E is expected to come from the entanglement between slow-moving particles, and the late time velocity v κ is small around κ = 0. Here we neglect the modes around κ ∼ ±π though their velocities are small. This is because these modes are expected to suffer from the lattice artifacts. Also, we have found the time evolution of ∆E and ∆I A,B for the adjacent interval in the fast CCP has a first dip at t ∼ 2ξ. On the other hands, that of ∆E and ∆I A,B for the adjacent interval in the slow one has the first dip at t ∼ 2ξ kz . Thus, the time interval from t = 0 to the time when first dip appears is given by twice as the effective correlation length: in the fast limit the initial correlation length ξ, in the slow limit the Kibble-Zurek length, ξ kz .

Future directions
• In the fast CCP with large d, ∆E and ∆I A,B have a similar behavior in our numerical computation. It is interesting to study whether the time evolution of ∆E is different from that of ∆I A,B in the very-late-time regime where we have not computed them.
• Another future direction is to understand why ∆E in the fast CCP with large d does not have the plateau in the window la+d 2 t l b +d 2 although ∆E in the fast ECP has the plateau in this window.
• In the fast ECP with large d, we can explain the time evolution of ∆E by using the quasi-particles which move at the speed of light. However, we need to utilize not only quasi-particles at the speed of light, but also the slowly-moving particles to explain the time evolution of ∆I A,B because of the logarithmic increase at the late time. It is important to investigate why we need such a slowly-moving particle picture in only ∆I A,B .