An exactly solvable quench protocol for integrable spin models

Quantum quenches in continuum field theory across critical points are known to display different scaling behaviours in different regimes of the quench rate. We extend these results to integrable lattice models such as the transverse field Ising model on a one-dimensional chain and the Kitaev model on a two-dimensional honeycomb lattice using a nonlinear quench protocol which allows for exact analytical solutions of the dynamics. Our quench protocol starts with a finite mass gap at early times and crosses a critical point or a critical region, and we study the behaviour of one point functions of the quenched operator at the critical point or in the critical region as a function of the quench rate. For quench rates slow compared to the initial mass gap, we find the expected Kibble-Zurek scaling. In contrast, for rates fast compared to the mass gap, but slow compared to the inverse lattice spacing, we find scaling behaviour similar to smooth fast continuum quenches. For quench rates of the same order of the lattice scale, the one point function saturates as a function of the rate, approaching the results of an abrupt quench. The presence of an extended critical surface in the Kitaev model leads to a variety of scaling exponents depending on the starting point and on the time where the operator is measured. We discuss the role of the amplitude of the quench in determining the extent of the slow (Kibble-Zurek) and fast quench regimes, and the onset of the saturation.


Introduction
The dynamics of a closed quantum system following a smooth quench involving a critical point is expected to carry universal signatures of the critical theory. Of course, the best known of such behaviours would be Kibble-Zurek (KZ) scaling [1,2], which has received considerable attention in the past several years [3][4][5][6][7][8]. One might characterize the corresponding quenches as "slow" since they involve a protocol where the evolution remains adiabatic until the system approaches very close to the critical point. However, recent holographic studies [9][10][11] also revealed interesting new scaling behaviour when critical points were probed by "fast" quench protocols. Later examination showed that this fast quench scaling is a universal behaviour for quantum field theories flowing from a UV fixed point, which is described by the conformal field theory [12][13][14]. Hence various renormalized observables in such continuum systems can exhibit a variety of scaling behaviours for different regimes of the quench rate. In fact, it was explicitly shown in [15] that scaling JHEP11(2017)157 behaviour of various (renormalized) observables smoothly interpolates between the Kibble-Zurek scaling and the fast quench scaling with appropriate quench protocols in free scalar and fermion field theories. However, the above results beg the question: how do these scaling regimes manifest in the presence of a finite UV cutoff? This question is particularly important if we wish to make any contact with real experimental systems which always have a finite lattice spacing (e.g., cold atom systems in optical lattices). The observables of interest are now bare quantities, and while we can vary the quench rate over a broad range of scales, the cutoff scale (e.g., the inverse lattice spacing) will always place limitations on fast quench protocols. However, we might still expect that when the quench rate is taken far below the cutoff scale, the results would match those for renormalized quantities for the corresponding continuum fixed point theory.
In this paper, we will explore quench dynamics in spin systems with a smooth nonlinear quench protocol for which the quantum dynamics can be solved exactly for any value of the quench rate. We will focus on two such models: the transverse field Ising model in one space dimension and the Kitaev honeycomb model [16] in two space dimensions. The ability to solve the dynamics stems from the fact that both these models can be written in terms of free fermions in momentum space. In both cases, we will measure the expectation value of the quenched operator, which in the fermionic language isψψ, where ψ is the fermionic field used to represent the integrable spin models via Jordan-Wigner transformation.
For quenches where the couplings vary linearly in time, the response of the (1+1)dimensional Ising model has been examined in [17,18], while the Kitaev model has been studied in [19][20][21]. Both the Ising [22][23][24] and the Kitaev [25,26] model have also been studied for instantaneous quenches. In contrast, we will study the quench protocols in which the couplings vary smoothly over a (finite) time interval (characterized by the duration δt) and saturate to constant values at early and late times. Among other things, this allows us to investigate the dependence on the quench rate and the amplitude of the quench separately, and to scan the entire range of quench rates, including the new fast scaling regime [12][13][14][15], as well as the regime where the quench rate becomes of the order of the UV cutoff. In the latter regime the response saturates as a function of the quench rate, as one expects for an instantaneous quench. However we find that the δt at which this saturation happens is proportional to the amplitude for large amplitudes, while it is independent of the amplitude for small amplitudes. While the Ising model has an isolated critical point separating two gapped phases, the Kitaev model is particularly interesting since there is a whole critical region in the space of couplings. This allows us to study a new class of quenches in which the couplings are varied entirely within this critical region.
The remainder of the paper is organized as follows: the remainder of this section provides a more detailed review of the various quench regimes, which we introduced in the introductory discussion above. In section 2, we review the two lattice models and in particular, we derive the relevant continuum theory describing the corresponding critical points. We then introduce our specific quench profiles for the couplings and discuss the exact time-dependent solutions in section 3. Our various results for the response of these JHEP11(2017)157 lattice models in different quench regimes are presented in section 4. Section 5 contains a brief discussion of our results. In appendix A, we derive the fast scaling behaviour using linear response theory around a CFT in arbitrary dimensions. In appendix B, we provide some technical details required to understand the amplitude dependence of the value of δt at which the response saturates. The final appendix C discusses the Cluster Ising model on a one-dimensional closed chain, which can also be studied in the same manner as the Ising model.

Quench regimes: slow, fast and instantaneous
Slow quench: for a quench protocol where the system starts with a finite mass gap and then crosses or approaches a critical point, at a rate which is slow compared to the initial gap, many systems show Kibble-Zurek (KZ) scaling [1][2][3][4][5][6][7][8]. Consider examining a relativistic system (i.e., the dynamical critical exponent is z = 1) with the simple power-law protocol g(t) − g c ∼ g 0 (t/δt) r . (1.1) It follows that the instantaneous energy gap E g (t) is given by where E 0 = κg 0 and ν is the correlation length exponent for the critical point. Now the initial adiabatic evolution breaks down at t = −t KZ , as determined by the Landau criterion dt t=−t KZ ∼ 1 . There are similar KZ scaling relationships for correlation functions. Given the profile (1.2) above, one finds that t KZ is related to the inverse quench rate δt by the scaling relation (1.5) We will be considering more general quench protocols where the instantaneous energy gap takes the form where the function f (x) → 1 at x → −∞. Only in the regime x 1 does the profile f (x) approach zero as f (x) ∼ x rν . Now imagine that the system starts in its ground state at t = −∞ and initially evolves adiabatically since E g (t) is only changing very slowly. However, as described above, this initial adiabatic evolution breaks down when the Landau criterion (1.3) is satisfied. At this point, the Kibble-Zurek scaling (1.4) described above can appear if two conditions are satisfied: (i) t KZ is such that the function f (x) can be JHEP11(2017)157 approximated by the power law x rν , i.e., we must demand that t KZ δt; and (ii ) at t = −t KZ , the instantaneous gap is much smaller than the UV cutoff scale Λ UV , i.e., we require E KZ ≡ E g (t = −t KZ ) Λ UV . Combining eq. (1.5) with the first condition yields Similarly, the second condition can be expressed as (1.8) When E 0 Λ UV , the condition (1.7) implies the condition (1.8). However, when E 0 is of the same order or larger than Λ UV , the condition (1.8) is the stronger restriction. Note that we are considering protocols where E 0 > E KZ . 1 This approach is closely related to quench protocols used more commonly in discussions of KZ behaviour in the condensed matter literature, where the behaviour (1.2) is often considered to hold for all times. In particular, if we define K ≡ E 0 /δt rν , the expression (1.2) for the instantaneous gap becomes E g (t) = K t rν . Note that in this case E 0 and δt can not be separately varied. The system is prepared in the ground state of the theory at some (finite) initial time t = t i . The KZ time and energy gap are then given by t KZ = K − 1 rν+1 and E KZ = K 1 rν+1 , respectively. Now the first condition above is replaced by |t i | ≥ t KZ , which can be equivalently written as (1.9) The second condition above, i.e., E KZ Λ UV , can be expressed as (1.10) Hence if E(t i ) < Λ UV , the first condition implies the second, while if E(t i ) Λ UV , the second constraint is the stronger one. Often one actually considers the situation where t i → −∞. In this case, the only constraint is eq. (1.10), since E(t i ) diverges rendering the inequality (1.9) trivial.
Fast quench: as noted above, recently a new scaling behaviour was discovered for smooth but fast quenches. In particular, consider a generic action where S CFT is the conformal field theory action describing the UV fixed point, and O ∆ is a relevant operator with conformal dimension ∆. The quench profile for coupling λ(t) starts from some initial value λ init and smoothly changes to the final value λ fin over a time scale 1 For the quench protocols and the lattice models studied in the following, we will have ν = 1 and r = 1.
Further, the UV cutoff scale is simply the inverse lattice spacing, i.e., ΛUV = 1/a, and so eq. (1.8) can be written as δt E0 a 2 .

JHEP11(2017)157
δt. If this time scale is fast compared to all physical scales, but slow compared to the scale of the UV cutoff, i.e., then the response of various renormalized quantities exhibit scaling at early times [12][13][14][15]. For example, during the quench process, the renormalized expectation value O ∆ ren behaves as (1.13) where δλ = λ fin − λ init . Similarly, the energy density scales as (1.14) This scaling behaviour was originally discovered in holographic computations [9][10][11] but then it was shown to hold in free field theories, and further argued to be true for general interacting theories [12,13]. It may seem mysterious that the underlying QFT has been regulated and renormalized and yet the above expressions are divergent in the limit δt → 0, when ∆ > d/2. These divergences arise because as δt shrinks, the quench is exciting a growing number of short wavelength modes. Further there is an infinite "reservoir" of such modes available as long as they are arranged as excitations of the UV fixed-point CFT. Implicitly the latter holds for renormalized quantities as in eqs. (1.13) and (1.14), which are defined in a procedure which involves taking the limit Λ UV → ∞. Of course, if the cutoff scale Λ UV is held fixed while δt continues to shrink, eventually we will encounter δt ∼ 1/Λ UV and the above scaling behaviour will no longer be applicable.
Instantaneous quench: the approach, which is most commonly discussed in the quantum quench literature, e.g., [3-7, 22-24, 27], involves preparing a system in the ground state of an initial Hamiltonian and then time evolving this state with a new or final Hamiltonian. There are some scaling results known to hold in the situation where the ground state of the initial Hamiltonian is gapped and the final Hamiltonian corresponds to a critical phase [23,24,27], in particular, the latter is a (1+1)-dimensional CFT. One can imagine that this describes an instantaneous quench, where initially the couplings of Hamiltonian are held constant, then at a single moment of time, the couplings are instantaneously changed to produce the final Hamiltonian and subsequently the couplings are fixed at their new values. Further this interpretation naively suggests that this protocol corresponds to the δt → 0 limit of the smooth fast quenches described above. However, as emphasized in [13,14], this limit does not reproduce the instantaneous quench 2 because implicitly the former assumes that quench rate is always small compared to the UV cutoff,

JHEP11(2017)157
i.e., Λ UV 1/δt. Instead, the instantaneous quench implicitly assumes that δt → 0 while Λ UV remains fixed. 3 While the above discussion applies quite generally, implicitly we are assuming 2∆ > d in which case eqs. (1.13) and (1.14) would produce divergences in the limit δt → 0. However, we can also consider the situation where 2∆ < d in which case the expressions in these formulae would vanish when δt → 0. In fact, these expressions no longer capture the leading contributions in this situation and instead the quench will produce finite results for O ∆ ren and E ren . In this case, the results in free field theories indicate that these final answers do not depend on UV details and the responses for the instantaneous and the smooth quenches will agree [14]. This behaviour is also manifest in the excess energy above the adiabatic value at late times, which is UV finite for any finite δt. Here again free field studies [14] show that for 2∆ < d, this quantity remains finite in the δt → 0 limit and reproduces instantaneous quench results. Instead it is a subleading term which displays scaling behaviour analogous to that in eq. (1.14). For 2∆ > d the smooth quench answers diverge in this δt → 0 limit displaying scaling, while the instantaneous quench answer has a UV divergence.
Having introduced a finite UV cutoff in the present paper, we are certainly able to study the new regime where Λ UV 1/δt, which we will refer to as the instantaneous quench regime. As noted above, the fast quench scaling does not apply in this regime. In particular, the divergences which the δt → 0 limit would produce in eq. (1.13) are avoided and instead we will see that the response saturates when we enter the instantaneous quench regime. 4 An interesting feature of this regime is that for protocols of the form (1.6), the value of δt where this saturation occurs is independent of the amplitude E 0 for small amplitudes, while it becomes proportional to E 0 for E 0 ∼ O(Λ UV ) or larger. We provide a physical explanation of this behaviour in terms of the value of δt at which the largest momentum mode k max (in lattice units) contributing to the observable departs from adiabatic behaviour. A detailed study of the Ising model reveals that for large amplitudes k max is at the UV cutoff, i.e., k max ∼ π. The Landau criterion for this mode leads to the above result. In other words, saturation happens when all possible modes get excited. On the other hand, we find that for small amplitudes k max ∝ E 0 , i.e., all modes do not contribute significantly to the observable we calculate. In this case, the Landau criterion shows that the saturation value is independent of the amplitude. 3 The papers [23,24,27] argue that the state which results from this kind of quench can be well approximated by a state of the form e −βH CFT |B where |B is a boundary state of the final CFT and HCFT is the final Hamiltonian. This approximation is expected to hold for IR quantities, however, some subtleties have been discussed recently in [28][29][30]. 4 As anticipated in [14], this behaviour is similar to that of correlation functions at finite spatial separation δ x in renormalized quantum field theories. In particular, when δt is small (as specified in eq. (1.12)) and δt > |δ x|, the correlation functions exhibit fast quench scaling analogous to eq. (1.13) but when δt < |δ x|, the correlation function saturates so that a finite δt → 0 limit exists.

JHEP11(2017)157 2 The models
In this section, we will describe the lattice models of interest and set our notation.

Transverse field Ising model
We write the Hamiltonian for the transverse field Ising model as where τ (i) denote the Pauli spin operators, and n denotes the site indices of the onedimensional chain. The time dependent coupling h(t) denotes the transverse magnetic field and J is the interaction strength between the nearest-neighbor spins. Both couplings have dimensions of energy here. We will follow the conventions of [31]. Using the wellknown Jordan-Wigner transformation, the Hamiltonian given by eq. (2.1) can be rewritten as the theory of a one-component fermion c(n) at each site, whose Fourier components will be denoted by d(q), The latter have the usual anti-commutation relations Note that q is periodic from eq. (2.2), i.e., the expression is invariant under q → q + 2π. However, it is convenient to shift the momenta and to introduce a two-component Majorana fermion . (2.5) The Hamiltonian then becomes where σ i denote Pauli matrices in the particle-hole space of fermions and we have introduced the (dimensionless) coupling g(t) = h(t)/J. Note that the momentum sum runs over half the Brillouin zone [31]. Now it is convenient to consider the limit N → ∞, in which we are considering an infinite chain of spins. This will certainly remove the possibility of having our quench results infected by any finite size effects. In this limit, the momentum k becomes a continuous variable on the range [−π, π] -although as noted above, k ∼ k +2π. We would have the more or less standard replacements:

JHEP11(2017)157
Hence it is convenient to rescale operatorsd(k) = √ 2N + 1 d(k) so that in the continuous limit, the anti-commutation (2.3) become Now constructing the fermionχ(q) from these rescaled operators as in eq. (2.5), the Hamiltonian (2.6) becomes For a given momentum k, the instantaneous energy eigenvalues are given by This dispersion relation makes clear that g = 1 corresponds to a critical point where the gapless mode is k = 0. 5 The critical point at g = 1 separates two massive phases (the paramagnetic and the ferromagnetic phases of the Ising model), and the continuum limit around this critical point is a massive Majorana fermion. This may be seen as usual by first expanding around the critical coupling with and then explicitly introducing the lattice spacing a to define the following dimensionful quantities: p = k/a , m(t) = (t)/a , and ψ(p) = a 1/2χ (a p) . (2.12) Finally we also define the dimensionless spin coupling:Ĵ = J a. 6 We then take the continuum limit with a → 0 holdingĴ, p, m(t), and ψ(p) fixed. In this limit, all of the terms in H Ising which are higher order in a vanish and we are left with the continuum Hamiltonian, which corresponds to the theory of a massive Majorana fermion with a (time dependent) mass m(t).
Before ending this subsection, we note that another related integrable model which shows similar behaviour is the Cluster-Ising model on a one-dimensional closed chain. This model is reviewed in appendix C, where it is shown that the corresponding Hamiltonian can be reduced to three copies of the continuum Ising Hamiltonian (2.13). Thus the exactly solvable quench protocols, which we will discuss below, also apply to the Cluster-Ising model.

Kitaev honeycomb model
This model in 2+1 dimensions is defined on a (spatial) honeycomb lattice. The Hamiltonian can be written as (2.14) where (j, l) denote the column and row indices of a site on the lattice. Typically, (2 + 1)dimensional models are not solvable, however, remarkably this model can be solved exactly for constant couplings J i [16], by rewriting it as a fermionic theory. We will use the fermionic theory, which results from the Jordan-Wigner transformation given in [33][34][35]. The latter introduces two sets of real fermionic fields obeying the standard anti-commutation relations.
To express the Hamiltonian in terms of these fermionic fields, we first denote the unit vectors in the x and y directions byî andĵ, respectively -see figure 1. Then the vectors span the reciprocal lattice. We also define the vectors, where n 1 and n 2 are integers. The vectors n denote the midpoints of the vertical bonds in the honeycomb lattice, i.e., the lattice sites are positioned at n ±ĵ/2. The Hamiltonian (2.14) can now be written as where the fermion a n lives on the site at the top of the vertical bond labeled by n while the fermion b n lives on the bottom site -see figure 1. 7 The quantity D n is an operator which is bilinear in the fermions, can take values ±1 on each link, and commutes with the Hamiltonian. This allows us to think of D n as representing a static Z 2 gauge field living on the links of the honeycomb lattice, which is coupled to the fermions. The key point which makes the Kitaev model integrable is that D n is conserved leading to an infinite number of conserved quantities; the ground state sector of the model corresponds to choice of D n = 1 on each link [16]. In the following we will study a quantum quench from the ground state in this model with protocols where J 1 and J 2 are held constant but J 3 is time dependent. Since the D n commute with H Kitaev , they remain unity throughout the full dynamics. In this case, we can set D n = 1 in eq. (2.17) making the time dependent Hamiltonian quadratic in the fermions. 7 We wish to emphasize that the labeling in figure 1 refers to the fermionized version of the Kitaev model, given in eq. (2.17). In particular, the bonds associated with J1 and J2 are reversed in the original spin Hamiltonian (2.14). This reversal of roles is perhaps not so surprising since the Jordan-Wigner transformation used to produce the fermionic description is nonlocal.   18) where N is the total number of sites (assumed to be even) and k extends over half the Brillouin zone. We then define a two-component spinor [19,20] The Hamiltonian (2.17) with a time dependent J 3 then becomes where σ i denote Pauli matrices in the space of fermions and we have defined

JHEP11(2017)157
Implicitly, we are again considering the limit of an infinite lattice size, which has resulted in replacing the momentum sum with an integral in eq. (2.21). The measure is defined as d 2 k ≡ dk 1 dk 2 , where the range of the integral is 0 ≤ k 1 , k 2 ≤ π [19,20]. For constant J 3 , the model is critical over a two-dimensional region in the (J 1 , J 2 , J 3 ) hyperplane. Indeed the gap vanishes when The solutions to these equations (2.22) can be represented by a triangle with sides of lengths J 1 , J 2 , J 3 and the angle between the (J 1 , J 3 ) sides being k 1 while that between the (J 2 , J 3 ) sides being k 2 . These conditions result in two bands of critical couplings satisfying The continuum limit of the model depends on whether the limit is constructed around a point in the interior of one of these critical bands, or around a point on one of the edges. To simplify our discussion of quenches in the following, we will set J 1 = J 2 = J > 0 and define J 3 ≡ −2J g(t) for which the Hamiltonian (2.20) simplifies to Within this two-dimensional space of couplings, the gapless constraints (2.22) reduce to 24) and the critical region becomes |g(t)| ≤ 1. If we think of g as a constant for a moment, there are three distinct classes of critical models corresponding to: 1) interior points with 1 < |g| < 0; 2) the edge points with g = ±1; and 3) the "interior" edge points with g = 0.
Although the latter lies in the interior of the critical region −1 ≤ g ≤ 1, it corresponds to the point where the edges of the otherwise two distinct bands described previously merge together with the choice J 1 = J 2 . We now consider the continuum limit associated with each of these critical models.
(1) Interior points with 1 < |g| < 0. We expect the continuum theory for such interior points to be a massless theory, as may be seen as follows: as before, we introduce a lattice spacing a and then expand around the critical point with g =ḡ + a m and k i =k + a m sink + a p i , i = 1, 2 (2.25) where cosk =ḡ. Further, we rescale the fields as and define the dimensionless couplingĴ = 3 √ 3 Ja/2. Now in the continuum limit a → 0, the couplingĴ, the momenta p i , the mass scale m and the field ψ( p) are held fixed. With this limit, eq. (2.20) yields the continuum Hamiltonian

JHEP11(2017)157
where we used eq. (2.15) to relate (p 1 , p 2 ) to the momenta along the x and y axes, i.e., Note that implicitly we have introduced the standard measure d 2 p = dp x dp y and these momentum integrals have infinite range. Further, the Jacobian arising in transforming from dp 1 dp 2 has been absorbed inĴ. Eq. (2.27) is the Hamiltonian of a massless Dirac fermion in 2 + 1 dimensions. Of course, to obtain the standard form of the Dirac Hamiltonian, one has to rescale the spatial coordinates by a factor depending onk, i.e., onḡ. Let us make a few additional comments: first, the (energy) scaling dimension of the momentum space field ψ( p) is −1. Therefore the corresponding scaling dimension of the position space field is +1, and the operatorψ( x)ψ( x) has scaling dimension 2, as expected in a (2 + 1)-dimensional relativistic theory.
Second, we note that in the following, our quenches result from introducing a time dependence in the coupling J 3 (t), or in the mass scale m(t) above. Given that the latter scale does not appear in the continuum Hamiltonian (2.27), it may naively appear that such quenches do not effect the continuum limit. However, note that m(t) also appears in the definition of the momenta p i in eq. (2.25). If we make a more conventional expansion removing the latter shift, i.e., we use Hence the dispersion relation becomes and we can see that a time dependent m(t) shifts the zero of this dispersion relation for the low-energy modes. Hence a time varying m(t) will produce a nontrivial quench. Of course, from the perspective of a relativistic field theory, m(t) appears here as an unconventional coupling for the continuum Hamiltonian. In fact, this coupling begins to reveal the anisotropic nature of the underlying lattice model.
(2) Edge points with g = ±1. For simplicity, let us focus on the lower edge of the critical region where g = −1 (i.e., J 3 = 2J) and hence where cosk = −1 (and sink = 0). Expanding around the critical point with

JHEP11(2017)157
Now we introduce a lattice spacing a and define as well asĴ = 3 √ 3 Ja/2. Then continuum limit is obtained by taking a → 0 while keepinĝ J, p x , p y , m and ψ( p) finite. The resulting continuum Hamiltonian then takes the form: At precisely m = 0, this is the well-known semi-Dirac point. This anisotropic theory has the dispersion relation As in the previous case, it is clear that a time dependent m(t) will produce an interesting quench. In contrast to eq. (2.31), this mass parameter can not be absorbed by a shift in the momentum of the low energy modes. However, for m ≥ 0, the theory is still gapless with E = 0 for (p x , p y ) = (±2 2m/3, 0). For m < 0, the theory has a gap with E gap = ±4J|m| at (p x , p y ) = (0, 0). Further in the latter case, for very low-lying modes, we might write the dispersion relation (2.36) as which has a form closer to the familiar relativistic Klein-Gordon relation. However, a quench with m(t) would differ from the familiar mass quench, e.g., [12,13] since the anisotropy between the x and y directions is also changed as m varies with time.
In position space, the Hamiltonian (2.35) would take the form Note that for this anisotropic critical point, the (energy) scaling dimension of the y coordinate is (−1) as usual, however, the scaling dimension of x is (−1/2). Further, the scaling dimension of the momentum space field ψ( p) is −3/4, and consequently the scaling dimension of the position space field ψ(x, y) is +3/4. The operatorψ(x, y)ψ(x, y) then has a scaling dimension 3/2.

JHEP11(2017)157
Now we introduce a lattice spacing a and define and, as before,Ĵ = 3 √ 3 Ja/2. Then continuum limit is obtained by taking a → 0 while keepingĴ, p x , p y , m and ψ( p) finite. The resulting continuum Hamiltonian takes the form: Note that from the scaling in eq. (2.41), p x has the standard (energy) scaling dimension of (+1) and further can be made as large as we like in the continuum theory. On the other hand, p y is dimensionless and implicitly, the above results are only valid of p y 1. We can remove the latter restriction by retaining the full nonlinearity of p y in the 'continuum' theory. This approach yields (2.43) where p y can take finite values above. However, we see that this momentum is periodic with period p y ∼ p y + 4π/3. In some sense, this scaling limit has only produced a continuum theory in the x direction and the y direction remains discrete. That is, we could interpret eq. (2.43) as the Hamiltonian of (coupled) fermions living on a family of one-dimensional defects, i.e., the position space field would take the form ψ(x, n y ), where x labels the position along the defects and n y labels on which defect the fermion resides. This behaviour is not unexpected since with J 3 = 0, the Kitaev Hamiltonian (2.14) reduces to a family of uncoupled one-dimensional spin chains stretching in the x direction. Here m introduces a small coupling between these chains. This unusual anisotropic theory (2.43) has the dispersion relation Note that with m = 0 (vanishing coupling between the one-dimensional defects), this dispersion relation becomes independent of p y . For this anisotropic critical theory, the (energy) scaling dimensions of the x coordinate is (−1) as usual. The scaling dimension of the momentum space field ψ( p) is (−1/2) and the scaling dimension of the position space field ψ(x, n y ) is (+1/2), as appropriate for a one-dimensional fermion. The operatorψ(x, n y )ψ(x, n y ) then has the scaling dimension 1, again as in a one-dimensional theory.

Quantization
The two lattice Hamiltonians given in eqs. (2.9) and (2.23) are both of the form where D is the number of (spatial) dimensions. The functions m( k, t) and G( k) are given in table 1 for the two models. Note that G( k) is an odd function of the momentum, i.e., G(− k) = −G( k). The two-component spinor χ( k, t) is written as For the Ising model, there is an additional Majorana condition Now consider the Heisenberg equation of motion for the above Hamiltonian, The two independent solutions may be expressed in the form where the scalar function φ( k, t) satisfies the equation Our aim is to quantize this theory with a time dependence of the couplings which saturate to constant values in the past and the future. In particular, we choose the profile where the function A( k) and the constant B are given in the table 2. The profile (3.8) was chosen because eq. (3.7) can be exactly solved for m( k, t) of the form given in eq. (3.9). The "in" positive energy solution, i.e., the solution which behaves as a pure positive frequency mode at t → ∞, is given by [13,36]  where we have defined

JHEP11(2017)157
Substituting φ in into eqs. (3.5) and (3.6), we get the solutions U in ( k, t) and V in ( k, t) which are positive and negative frequency respectively in terms of the "in" modes. The field can be now expanded in terms of the "in" oscillators where the usual anti-commutation relations hold, i.e., Further, the Majorana condition (3.3) requires a( k) = b( k) for the Ising model. One can similarly introduce the "out" modes, or for that matter, any Bogoliubov transform of these modes.
In studying these quenches, we will begin the system in the ground state of the Hamiltonian at t → −∞. The Heisenberg picture state is then the "in" vacuum (3.14) We will examine the quenches by following the expectation value of local bilinears of the fermionic operators in this "in" vacuum, i.e., c n for the Ising model and (a n , b n ) for the Kitaev model. We will consider the fermion bilinear In terms of the two-component momentum space fermion field, 8 these expectation values become as 8 Implicitly, we are definingχ = χ † σ3 in both models.

JHEP11(2017)157
where φ in is the solution given in eq. (3.10) for a given protocol, i.e., for specific values of the constants a and b in eq. (3.8). 9 A measure of the excitation of the system is given by the difference between this quantity measured in the quench and its adiabatic value The adiabatic value is obtained by replacing the exact solution by the lowest order adiabatic solution, and m(k, t) is given by eq. (3.9) for a specified time t. Substituting this solution into eq. (3.17), the adiabatic value at time t simplifies to

Results
In this section, we summarize our results of the calculation of χχ diff in eq. (3.16) for various quench protocols in the two lattice models described in section 2.

Transverse field Ising model
The Ising model has an isolated critical point which separates two massive phases. The most interesting quench protocol in this case is a g(t) which starts in one of the phases, crosses the critical point at g = 1 (at time t = 0) and ends in the other phase at late times. Hence with the profile in eq. (3.8), we choose a = 1 and consider quenches with various values of b.
To begin, we consider small values of b so that throughout the quench, the model is close to the critical point and we may expect that the results can be compared to the continuum limit. That is, as in eq. (2.11), our profile has the form g = 1 − (t) with (t) = −b tanh(t/δt), i.e., b controls the amplitude in the variation of the dimensionless "mass" parameter (t). In particular, in = b (and out = −b). Further, as in [12][13][14][15], we examine the effect of varying the quench rate by varying δt. Following the results of [15], we can expect to see different scalings in ψ ψ diff for fast and slow quenches. In passing,

JHEP11(2017)157
we note that the only dimensionful quantity in Hamiltonian (2.9) is the overall factor of J, the nearest neighbor bond strength, and as noted in footnote 6, this coupling essentially defines the lattice spacing. Hence J δt is the natural dimensionless quantity which can be used to discuss the different quench rates. Further, let us note that in the lattice models, χχ is a dimensionless quantity which does not scale with J, e.g., one finds that the factors of J cancel out in eq. (3.19). Figure 2 shows the response χχ diff at t = 0 as a function of the inverse quench rate Jδt. The two cases shown in the figure begin at t = −∞ with in = 0.01 and 0.1, and the plots show several clear features: 1. First, for small Jδt of order one or less, the response saturates as a function of the quench rate. We can think of this as the "instantaneous quench" regime, where the quench rate is of the same order as the lattice spacing. We will discuss a comparison with the results of an instantaneous quench for the Ising model, as well as the Kitaev model, in section 4.3. We will further discuss the saturation point for small amplitudes later in this section and provide more details in appendix B.

For
Jδt roughly between 1 and 1/ in , the quench time scale is much larger than the lattice scale (and so the results should be comparable to the continuum theory), but dimensionless combination in Jδt is small. This is the "fast quench" regime, as defined in [12][13][14][15] for the continuum theory. The conformal dimension of the operatorχχ(x) is ∆ = 1 (which matches one half of the spacetime dimension), and so the continuum result (1.13) suggests that there should be no leading power law dependence on δt. Instead, in this regime, the curves in the figure are well fit with a dependence of the form P log(δt) + Q, which is exactly what is expected from the continuum calculations. An explicit derivation of the logarithmic dependence in this regime is given in appendix A. One might expect that the response will saturate when the quench is exciting all possible modes in the lattice theory, i.e., when excitations are being created in all modes. To understand this point, we begin by substituting the profile (3.8) (with a = 1 and b = in ) into the dispersion relation (2.10) for the Ising model and we find Now we can ask when a mode at a particular wave-number k is going to be excited. According to the Landau criterion (1.3), this mode will remain adiabatic throughout the quench (and in particular, at t = 0) if On the other hand, if the above expression exceeds one, then this mode with momentum k is excited by the quench. Now it is clear that there will always be excited modes in the neighborhood of k = 0. The above discussion may now suggest that saturation will be achieved when the violation of the above criterion (4.2) extends out to k = π. However, it turns out that for small amplitude quenches, only a narrow band of momenta near k = 0 contribute significantly to the expectation value (3.16). In particular, we show in appendix B that χχ diff only receives significant contributions from modes with 0 ≤ k < k max with k max = c | in | where c is some order one number. Hence the violation of the criterion (4.2) need only extend to k = k max in order to produce saturation with small amplitudes. Hence it is straightforward to see from eq. quenches will saturate for Jδt ≤ 1/(4c) or more simply Jδt 1. 11 Hence we have an explanation of the saturation behaviour observed for the small amplitude quenches shown in figure 2. To contrast with the following, we emphasize that the point where saturation sets in is independent of the initial amplitude here. Now we can also examine the scaling behaviour of χχ diff for "large-amplitude" quenches of the Ising model, i.e., with | in | 1. In this regime, 12 we expect that the quenches are probing the nonlinear regime of the lattice dispersion relation (2.10). The response as measured by the expectation value ofχχ is shown for a family of four such quenches in figure 3. These large-amplitude quenches exhibit three distinct features, which contrast with the behaviour found for the small-amplitude quenches: a) there is no longer a fast quench regime; b) the Kibble-Zurek scaling regime begins at Jδt ∼ | in | (rather than 1/| in | for small amplitudes); and c) the response saturates to the instantaneous quench value for Jδt | in |/8 (instead of the previous Jδt 1). We now discuss each of these in somewhat more detail.
The three features above are related because an essential characteristic of these large amplitude quenches is a smooth transition between the KZ and instantaneous quench regimes, without an intermediate scaling regime. 13 In the KZ regime, we still see scaling JHEP11(2017)157 compatible with ∆ = 1 expectation, i.e., χχ diff ∼ δt −1/2 , while the response is saturated in the instantaneous regime. Recall that the origin of the fast quench scaling (1.13), is that the quench activates a growing number of short wavelength modes as δt shrinks, but further that these new modes are organized as in a CFT [12,13]. In the transition between the KZ and instantaneous quench regimes, the large amplitude quenches are already exciting modes in the nonlinear regime of the Ising dispersion relation (2.10), i.e., the quench is probing modes with wavelengths comparable to the lattice spacing. Hence the latter condition is not achieved in the large amplitude quenches and we should not expect to see a fast quench scaling regime.
The second important feature was that the KZ scaling sets in at a quench rate which grows with the amplitude, rather than decreasing as in the small amplitude quenches. That is, we observe that the KZ scaling begins when Jδt ∼ | in | (instead of 1/| in |) in figure 3. This (naively) surprising behaviour can be understood by looking carefully at the conditions required for KZ scaling, as already discussed in the introduction -see the discussion around eqs. (1.7) and (1.8). We repeat the salient points here: the first condition for our quench protocol (3.8) was that adiabaticity should breakdown when the profile is in the linear regime. That is, we must have t KZ < δt, which in turn yields Jδt > 1/| in |, when expressed in terms of eq. (1.7). The second condition for KZ scaling to hold was that the system has to be close to critical at t = t KZ , i.e., E KZ < Λ UV . Interpreting J as the inverse lattice spacing (see footnote 6), this condition expressed as eq. (1.8) yields Jδt > | in |. Of course, when | in | is small, the first condition is stronger and we recover the behaviour observed in the small amplitude quenches. However, for the large amplitude quenches, it is the second restriction that determines the onset of KZ scaling, in agreement with the results found in figure 3.
The third and final feature observed in figure 3 is that the saturation point, where the instantaneous quench regime begins, also grows with increasing amplitude in the large amplitude quenches. In fact, we observe that the expectation value saturates into the instantaneous quench value for roughly Jδt | in |/8, instead of Jδt 1 as observed for the small amplitude quenches. This behaviour can be understood by the same reasoning used above in discussing the small amplitude quenches. The key difference is that for large amplitude quenches, the response (3.16) does, in fact, receive significant contributions from all modes, i.e., 0 ≤ k π -see appendix B for details. Hence we should ask when is the highest mode going to be excited 14 and substituting k = π into eq. (4.2), we find Therefore all of the modes will be excited when the above expression is bigger than one, i.e., for Jδt | in |/8, and hence we expect that the expectation value will saturate in this regime for large amplitude quenches.
14 Strictly speaking, the contribution of the k = π mode to the expectation value χχ diff is zero for any amplitude. However, this same argument holds for any large k. In the case of large amplitudes, these modes have significant contributions to χχ diff and so the following argument applies. See appendix B for more details.

JHEP11(2017)157
Finally, we should note that the quenches were studied here by measuring response exactly at the critical point. That is, we evaluated the expectation value precisely at t = 0. One can also evaluate the response at some finite time τ ≡ t/δt. At least for small amplitudes, the response of the continuum theory should provide a guide [12][13][14][15]. Hence we expect that for large enough Jδt, the KZ scaling regime should give way to an adiabatic regime. In fact, we expect to see this transition around Jδt ∼ 1/( in τ 2 ). As shown in [13], the adiabatic expansion for relativistic theories will be a series in 1/δt 2 and hence to leading order, we expect that after subtracting the zeroth-order term for large Jδt, χχ diff ∼ δt −2 . We have verified that the Ising model indeed produces this behaviour. However, it is technically challenging to separate clearly all four regimes (instantaneous, fast, KZ and adiabatic) for a fixed initial amplitude and finite time.

Kitaev honeycomb model
Recall that to simplify our discussion of quenches in the Kitaev model (2.20), we restricted our attention to the subspace within the full space of couplings where J 1 = J 2 = J > 0 and J 3 = −2J g. With these restrictions, the Hamiltonian reduces to that given in eq. (2.23). The critical region where the energy gap vanishes in this space of couplings reduces to the one-dimensional line segment |g| ≤ 1, with g = cosk and k 1 = k 2 =k as in eq. (2.24). Further, as discussed in detail in section 2.2, there are three classes of critical models: 1) interior points with 1 < |g| < 0; 2) edge points with g = ±1; and 3) "interior" edge points with g = 0. There are also gapped phases with |g| > 1.
As described in section 3, we quench the system with g(t) = a + b tanh(t/δt). In the following, we will always examine the response χχ diff , as defined in eq. (3.16), at t = 0. Clearly, there is a wide variety of different quenches depending on the choice of the parameters, a and b. In particular, as we describe below, the results depend crucially on the phase in which the quench begins and on the phase at t = 0 where we measure the response, i.e., the scaling of the responses depends on g(t → −∞) = a − b and on g(t = 0) = a. Of course, χχ diff (t = 0) will not depend on the profile of g(t) at latter times t > 0.
Given the three different types of critical theories, there are certainly a wide variety of critical quenches which one might choose to explore. In the following, we will focus on three protocols: a) 'gapped-to-edge' quenches which begin in the gapped phase and are measured at the edge critical point; b) 'gapped-to-interior' quenches which begin in the gapped phase and are measured at an interior point at some finite distance into the critical region; and c) 'interior-to-interior' quenches where the entire protocol only passes through interior critical points. Clearly this selection is not exhaustive and only provides a preliminary study of the critical quench dynamics of the Kitaev model -see section 5 for a discussion of other possible protocols.
One feature common to all of the different quenches is that for Jδt 1, the response saturates as a function of the quench rate. As described for the Ising model above, we can think of this as the "instantaneous quench" regime, where the quench rate is of the same order as the lattice spacing. In the discussion of the individual protocols below, we JHEP11 (2017)157 focus on the scaling of the response for the regime Jδt > 1. We return to consider the instantaneous quench regime in section 4.3.

Gapped-to-edge
Here, we consider quenches which start in the gapped phase with g(−∞) < −1 and pass to the edge point with g(0) = −1. That is, we consider profiles (3.8) with a = −1 and b > 0 (and hence a − b < −1). 15 In this case, b sets the scale of the gap in the initial phase. Further note that even though the system continues into the region of interior critical points for t > 0, these protocols are identical to a quench from a gapped phase which crosses an isolated critical point, since a measurement at t = 0 does not care about the values of the coupling for t > 0. Figure 4 shows the response for a quench with a small amplitude, i.e., |b| 1. We see that there are two distinct scaling behaviours for χχ diff (t = 0): with exponent −1/2 in the fast quench regime with 1 Jδt 1/|b|, and with the exponent −3/4 in the slow quench regime with Jδt 1/|b|. Let us note that the scaling exponent −3/4 for slow quenches was also observed in earlier work [21], where the variation of the coupling was taken to be always linear in time and the quench was started at t = −∞.
Again for small amplitudes, we expect that the quench is described well by the anisotropic continuum theory, given in eq. (2.38). Recall from the discussion below eq. (2.38), that the dimension of the operatorψψ is ∆ = 3/2. However, the anisotropy of the theory is important to identify the fast scaling dimension. In particular, while t and y scale as regular coordinates with mass dimension -1, the dimension of the x coordinate was -1/2. Hence the effective spacetime dimension in various formulae is d eff = 5/2, rather than d = 3. For example, the dimension ofĴm, the coupling conjugate toψψ, is d eff − ∆ = 1 and not d − ∆ = 3/2. Hence using d eff in eq. (1.13), we find the scaling 16 for fast quenches, in agreement with the results noted above and shown in figure 4. In the slow quench regime, the Kibble-Zurek time is given by t KZ = (δt/Ĵ m) 1/2 , since in eq. (1.5) ν = 1 and E 0 =Ĵm, and so the response (1.4) becomes (4.5) Again this reproduces precisely the exponent −3/4 found in our numerical results. 17 15 Of course, the results for quenching to the edge point at g = +1 are identical. One needs to simply flip the sign of all of these parameters, i.e., g, a, b → −g, −a, −b. 16 Note thatχχ is a dimensionless quantity, and in accord with footnote (6), we are canceling powers of J and the lattice spacing a in converting the first expression to the second, i.e., we set Ja ∼ 1. 17 One can also try to match the b dependence of χχ in eqs.  Of course, we can also extend these quenches to large amplitudes, and as shown in figure 5, the behaviour is somewhat different. First, as expected, we do not observe any fast scaling in this case. As discussed for the Ising model above, in the transition between the instantaneous and slow quench regimes, the large amplitude quenches are probing modes with wavelengths comparable to the lattice spacing and hence they cannot be described by an effective UV CFT. Hence the fast quench scaling (1.13) is not produced in these large amplitude quenches. However, there does appear to be a slow scaling regime. The exponent, however, changes continuously as we increase the amplitude from χχ diff ∼ (Jδt) −3/4 to χχ diff ∼ (Jδt) −1/2 . Of course, it will be interesting to develop an analytical understanding of this change in scaling behaviour between the small and the large amplitude quenches. This new result does not contradict the results of [21]: in the regime where we see this scaling, the coupling is proportional to t with a coefficient which is not quite small, whereas the results of [21] refer to a regime where this coefficient is very small in units of 1/J. Finally, in parallel to what happens in the Ising case when the amplitude is large, we see that the saturation point for small Jδt increases roughly linearly as we increase b.
As a final note, let us add that in the next section, we will see that for small amplitudes, the exponent for both the fast and slow quench regimes begins to change as soon as we continue the quench into the gapless phase.

Gapped-to-interior
Next, we consider quenches where we start in the gapped phase with g(−∞) < −1, pass beyond the edge point and measure the response at an interior point with 0 > g(0) > −1. That is, we are studying quench profiles (3.8) with a − b < −1 and 0 > a > −1 (and b > 0). In this case, there are many different quench protocols that can be studied and that yield different responses. To concisely go through them, it will be convenient to define two new parameters: δg in ≡ b − a − 1 (> 0), which sets the scale of the gap in the initial phase, i.e., measures the initial distance of g to the critical region, and δg fin ≡ 1 + a (> 0), which measures the final distance of g inside the interior critical region, i.e., the distance of the point at which we measure the response from the edge point. There will be three clearly distinct behaviours depending on whether δg fin is much smaller than, greater than, or the same order as δg in .
Let us start by considering the case in which δg in > δg fin . In the previous section, we already analyzed the case where δg fin = 0, which corresponds to quenching to the edge point. In that case with small amplitudes, we found that the expectation values scale as δt −1/2 in the fast regime, and δt −3/4 in the slow one. Now we hold δg in fixed and slowly increase δg fin away from zero. We observe different scaling behaviours depending on whether the quench is slow or fast compared to δg in , as shown in figure 6 where the initial amplitude is fixed to δg in = 0.1. Note that the scaling exponent in the slow quench regime, i.e., 2Jδt > 1/δg in , immediately changes from −3/4 to −1/2 as soon as δg fin > 0. The latter exponent corresponds to the KZ scaling of a (1+1)-dimensional fermionic mass quench. The situation is different for the fast quench regime. In this case, the scaling behaviour varies continuously from χχ diff ∼ (Jδt) −1/2 at δg fin = 0 to a logarithmic scaling when δg fin ∼ δg in .
The three examples shown in figure 6 were chosen to show the evolution of the scaling behaviour . Plot of χχ diff as a function of 2Jδt for the gapped-to-gapless quenches shown in the inset. The brown and red curves are the best fits in the fast and slow quench regimes, respectively.
As an example, we show here the transition between 3 different scalings. In the inset, the three quenches protocols are shown. All of them start at a distance of 0.1 from the critical region. The blue protocol goes to the edge at t = 0, g(t/δt) = −1+0.1 tanh(t/δt); the yellow protocol, just enters the critical region at t = 0, g(t/δt) = −0.96 + 0.14 tanh(t/δt); and the green protocol is well inside the critical area at t = 0 (compared to the initial amplitude), g(t/δt) = −0.85 + 0.25 tanh(t/δt). As in the previous subsection, the scalings of the blue dots correspond to the gapped-to-edge quench: in the fast regime it goes as P + Q(2Jδt) −1/2 (brown curve) with P = −0.00482 and Q = 0.000750; in the slow regime, it behaves as P + Q(2Jδt) −3/4 (red curve) with P = −0.00448 and Q = −6.97 × 10 −6 . Now, when you measure inside the critical area, the slow scaling changes instantaneously to −1/2, while the fast quench continually changes from −1/2 to a logarithmic scaling. Thus, the two regimes in the yellow quench are given by P + Q(2Jδt) −α (brown, fast), with P = −0.00963, Q = 0.000884, α = 0.406 and P + Q(2Jδt) −1/2 (red, slow), with P = −0.0107,Q = 0.000439. Finally, the green curve shows the prototypical example of a gapped-to-interior quench with best fits P + Q log(2Jδt) (brown, fast), with P = 0.00690, Q = −0.0213 and P + Q(2Jδt) −1/2 (red, slow), with P = −0.0170 and Q = −0.000200. and hence these are quenches to the edge, as in the previous section, with a scaling exponent In figure 7, we explore the transition between the edge scaling to logarithmic scaling in the fast quench regime. We analyze different quench protocols, all with fixed δg in = 0.1 and with δg fin varying from 0 to 0.15 -see figure 7a. In figure 7b, we show the corresponding scaling exponent in the fast quench regime fit with χχ diff ∼ (Jδt) −α . We see that the exponent α begins at 1/2 and smoothly decreases until some value just below 0.1 when JHEP11(2017)157 δg fin ∼ 0.1 At this point, α saturates for large δg fin . In fact, the scaling of the expectation value has become logarithmic, but it turns out that this is indistinguishable from power-law scaling with the small exponents shown in the figure.
In summary, the (small amplitude) quenches from the gapped phase to the interior of the critical region exhibit two scaling regimes: the slow quench regime with 2Jδt > 1/δg in where χχ diff always scales as (Jδt) −1/2 for any δg fin > 0, which is distinct from (Jδt) −3/4 scaling observed with δg fin = 0. We might add that this behaviour would be the KZ scaling for a (1+1)-dimensional fermionic theory. In the fast quench regime with 2Jδt < 1/δg in , the scaling behaviour makes a smooth transition from the quench-to-edge scaling of (Jδt) −1/2 to a logarithmic scaling, that would be expected for a (1+1)-dimensional fermionic mass quench. This transition occurs over the range 0 < δg fin δg in . Now it is natural to ask whether the same scaling holds for large amplitude quenches inside the critical region, i.e., for quenches with δg fin δg in . These large amplitude quenches are explored in figure 8. As the whole region on interior critical points is traversed with δg fin = 1 to analyze larger amplitudes, we need to start with smaller δg in . For the quenches shown in figure 8, we fix δg in = 0.01 and vary the final amplitude δg fin = 0.1, 0.25, 0.45, 0.75, all of which satisfy δg fin δg in . In this large amplitude regime, the corresponding fast and slow quench scalings remain the same as above, i.e., −1/2 and logarithmic, respectively. However, now the transition between the two scaling regimes is no longer at 2Jδt ∼ 1/δg in = 100, rather we find the transition at 2Jδt ∼ δg fin . In fact, we see that as δg fin increases and becomes of order one, the fast quench scaling regime shrinks more and more until is no longer observable, as shown with the violet dots in figure 8 which predominantly scale with the slow quench scaling.
Another alternative approach to large amplitude quenches is to consider quenches where we fix δg fin but we increase the gapped amplitude δg in . We followed this approach with δg fin = 0.1 fixed and varying δg in = 0.1, 1, 10, 100. When both δg's are comparable, we found the same result as before; namely, a logarithmic fast scaling regime and a power-law slow regime with exponent −1/2, separated at a scale δg −1 in ∼ δg −1 fin . As we increased δg in , the logarithmic scaling remains but the exponent of the power-law starts decreasing, until around δg in = 10, where we can only observe a logarithmic behaviour. This is a rather surprising behaviour, given that in every other case the fast scaling was the regime which disappeared for large amplitudes. One possibility is that there is indeed still a slow quench regime (which sets in at some scale given by a combination of δg in and δg fin ) that will only appear for large enough Jδt. 18 We do not have a good understanding of the various scaling behaviours described above. The scaling exponent of −1/2 in the slow quench regime has been observed previously for quenches linear in time for this model [19,20]. These linear quenches began at t = −∞ and the response is measured at t = +∞. With these simple protocols, the response can be examined analytically and the 1/(Jδt) 1/2 scaling stems from the fact that the excitation probability predominantly depends on one of the directions in momentum space. Recall  that the critical models (2.27) for the interior points are not really anisotropic and so we expect that this result must be related to the fact that the quenched operator itself is anisotropic. That is, from the critical Hamiltonian in eq. (2.30), we see that the quenched operator corresponds to 19ψ γ x ψ. Hence we should think of m(t) as the x-component of a vector coupling and accordingly, the roles of p x and p y are clearly distinguished in the dispersion relation (2.31) -see further discussion in section 5. Of course, for the protocols used in this paper, our numerical results exhibit a behaviour similar to that in the linear quenches, as long as the gapless interior region is traversed, irrespective of the starting point. For slow quenches, the change of the exponent from 3/4 to 1/2 as soon as δg in = 0 may be understood as follows: consider first the line k 1 = k 2 =k, so that G(k) = 0. In this case the equations for χ 1 and χ 2 , as defined in equation (3.2), decouple. The solutions of the Dirac equation can be readily written down where A 1 , A 2 are integration constants and To determine the "in" solution, we need to examine the behavior at t → −∞, while for the negative energy solution we must set A 1 = 0, (4.9) Substituting these in the mode expansion for the operator χ as in eq. (3.12) and imposing the anticommutation relations then determines the integration constants A 1 = A 2 = 1 upto a phase. Therefore for these k 1 = k 2 =k modes, we have which is independent of k and time.
Consider now the response which we measure, viz. the quantity χχ diff defined in eq. (3.16). The adiabatic modes at time t = 0 for k 1 = k 2 =k are
We therefore conclude that along the k 1 = k 2 =k the contribution to the response, χ( k)χ( k) diff is always independent of k. This vanishes for all k for gapped to edge quenches, whereas for gapped to interior quenches there is a range of k for which this has the maximal value −2.
Consider now the response for generic k 1 and k 2 . The components χ 1 , χ 2 satisfy the equation For slow quenches the response is substantial for times when the quench profile can be approximated by a profile which is linear in time. For these times, eq. (4.14) becomes Rescaling t → t = t/ √ δt it is clear that solutions have a functional form This, in turn, implies that the quantity χ( k)χ( k) diff also has this functional form.
Since we have shown that this quantity is independent of momenta when G(k) = 0, the simplest form of this function at time t = 0 is where F 1,2 are some functions and α is some positive real exponent. The ellipsis indicates higher orders in G( k) √ δt. For gapped-to-edge quenches, we have also shown that χ( k)χ( k) diff vanishes for all k. Therefore for such quenches, the first term must be absent, or we should set c 1 = 0. On the other hand, for gapped-to-interior quenches, χ( k)χ( k) diff is nonvanishing for some range of momenta. Therefore for such quenches the first term is generically nonvanishing.

JHEP11(2017)157
Since we are considering slow quenches, most of the contribution will come from the region in momentum space where G(k) is small. In this region we can approximate G(k) ∼ 2k − cos(k + ) , (4.18) where k ± = 1 2 (k 1 ± k 2 ). Thus, for gapped-to-interior quenches, we have On the other hand, for gapped-to-edge quenches, most of the contribution should come from the single gapless point which has k + = π, k − = 0. Expanding around this point, k + = π + δk + and k − = δk − where both δk ± are small we have Therefore, we have To extract the leading δt dependence we rescale δk − → δk − √ δt and δk + → δk

Interior-to-interior
In this section, we consider quenches where both g(−∞) and g(0) are at interior critical points -and in fact, where the entire protocol from t = −∞ to 0 passes only through interior points. That is, we consider profiles (3.8) with −1 < a < 0 and −1 < a − b < 0 (as well as b > 0). Typical results are shown in figure 9 where we consider two different protocols one with b = 0.1 and the other with b = 0.01. One somewhat surprising feature is that we, in fact, observe two different scaling regimes, separated at a scale of Jδt of the order of the inverse of the amplitude, i.e., Jδt ∼ 1/b, even though the quench only travels across interior critical points at all times. The two scaling regimes are characterized by distinct scaling exponents. For 1 < Jδt < 1/b, we found that expectation values scale as 1/(Jδt), which is consistent with the fast scaling of a (2+1)-dimensional fermionic mass quench. On the other hand, in the slow regime where Jδt > 1/b, we find χχ diff ∼ (Jδt) −1/2 , which matches the slow quench scaling found for gapped-to-interior quenches in the previous section.
Of course, since the amplitudes considered in figure 9 are still small, 20 it would be very interesting to develop an understanding of these scalings in terms of a continuum theory described in section 2.2. In particular, there is a striking fact to understand, namely, why is there a slow quench regime at all? Since the quench only passes through interior points where the model is gapless, there is no intrinsic energy scale with which to compare the quench rate, i.e., the quenches are never in an adiabatic regime, no matter how large Jδt becomes. We checked this last fact by computing expectation values at a finite time but never observing an adiabatic evolution. 21 However, the numerical computations still show that there is a slow scaling regime for these quenches, which therefore cannot be adiabatic and somehow the scaling exponent matches the KZ scaling for a mass quench in (1 + 1)-dimensional fermionic theory! We also reiterate that this behaviour also matches the slow quench scaling found for gapped-to-interior quenches in the previous section. There it was suggested that this unusual scaling was associated with the anisotropy of the quenched operator.
One observation, which provides a step towards a possible explanation, is the following: in [12,13], it was shown that the fast quench scaling behaviour essentially follows from linear response theory. In our case, the dimensionless parameter which controls the renormalized perturbation theory is Jb δt. Thus the linear response would be valid when Jb δt is small, and indeed in this regime, we get the expected result from the continuum theory. That is, we find the fast quench scaling, i.e., 1/(Jδt), for a mass quench in a (2+1)-dimensional fermionic theory. On the other hand, when Jb δt > 1, the linear response calculation is no longer valid, and the arguments which lead to the fast quench scaling do not hold any more. Indeed the change of the scaling exponent changes precisely around Jb δt ∼ 1. While this does not explain the new scaling (Jδt) −1/2 found beyond this point, it does indicate that we should expect that the quenches are entering a new regime here.

Instantaneous quench limit
Irrespective of the particular model under consideration, when the quench rate is faster than the lattice scale, i.e., Jδt 1, we expect that our results should agree with that of an instantaneous quench in which the coupling is switched abruptly from g(t = −∞) to g(t = 0) at the time of measurement. In particular, the system has no time to respond to the change in the coupling and so to a good approximation, we have instantaneous quench : cos k 1 + cos k 2 − 2g(t) (cos k 1 + cos k 2 − 2g(t)) 2 + (sin k 1 − sin k 2 ) 2 . (4.24) Now as an example, figure 10 shows the result for χχ diff at t = 0 for a quench protocol which starts in the gapped phase with g(t = −∞) = −1.01, and ends in the interior of the critical region with g(t = 0) = −0.9. The figure also shows the instantaneous quench result (4.23) evaluated with eq. (4.24). The figure clearly shows that in the limit Jδt → 0, the exact numerical results smoothly approach the instant quench value (4.23). Similarly, we have verified that the saturation value of χχ diff agrees with eq. (4.23) in the Ising model.
Beyond t = 0. The exact result for an instantaneous quench in the Ising model from a coupling g = g 0 to g = g 1 was calculated in [37], where Hence in the Ising model, we can compare our full solution χχ (t) for quench rates faster than the lattice scale with this instantaneous quench answer, for all times t > 0. In figure (11), we compare the time dependence of in 0|χχ|0 in at small values of δt with the exact instant quench result eq. (4.25) and we see that the agreement gets better as Jδt decreases.

Discussion
In this paper, we have studied critical quench dynamics in the transverse field Ising model (2.9) on one-dimensional chain and in the Kitaev honeycomb model (2.20) in two dimensions. We studied an exactly solvable quench protocol which asymptotes to finite values of the coupling at early and late times, and focused on the response of the operator by which we carried out the quench. The exact solutions in terms of free fermions were used to study the scaling of the response with the quench rate. Our answers have been obtained in the thermodynamic limit with an infinite number of sites and thus are free from finite size effects. Our results are summarized in table 3.
The results for the Ising model, discussed in section 4.1, are the ones which are best understood. The Ising model has an isolated critical point. Our quench protocol takes us through this point and we measure the response at the moment (chosen to be t = 0) when the quench hits the critical point. For small amplitude quenches (with in = b 1), there are three different regimes, as shown in figure 2: 1) for Jδt 1, the response saturates and we are in the instantaneous quench regime; 2) for 1 Jδt 1/ in , the response scales logarithmically with Jδt, as expected from the continuum description of the fast quench regime; and 3) for Jδt > 1/ in , the response scales as 1/(Jδt) 1/2 , as expected for Kibble- Table 3. A summary of our results for the scaling behaviour in critical quenches of the transverse field Ising model (2.9) and the Kitaev honeycomb model (2.20). We have also included the analogous results for a free fermion in two and three dimensions, for comparison. The columns 'Slow' and 'Fast' indicate the scaling exponent in the response, i.e., χχ diff ∼ (Jδt) −α , for the slow and fast quench regimes -'log' indicates a logarithmic scaling was found and 'none' indicates the fast scaling regime disappears. The column 'Transition' indicates the approximate value of Jδt when the scaling makes a transition between the two scaling regimes, or the minimum value for the slow scaling when there is no fast scaling regime. For the continuum free fermion, we are indicating the value of δt at this transition.
Zurek scaling in continuum. Hence for quenches which only make small excursions from the critical point, our results agree with those expected for the continuum theory [12][13][14][15]. However, we can also consider large amplitude quenches (with in = b 1) and the scaling behaviour changes in this regime. In particular, the fast quench scaling regime disappears and rather there is a smooth crossover between the instantaneous and slow quench regimes. Further, the slow scaling behaviour matches the KZ scaling found above, but the transition into this regime occurs roughly when Jδt ∼ in .
As described in section 4.1, for the Ising model, we have a good theoretical understanding of the response in all of the different situations described above. Perhaps one of the most interesting results here is that there is a fast scaling regime for small amplitude quenches. That is, our analysis of the Ising model confirms that even with a finite lattice spacing, certain quench protocols produce the fast scaling behaviour originally discovered in the study of continuum field theories [9][10][11][12][13] -see further discussion below.
The Kitaev model is distinguished by having an extended region of couplings for which the theory is gapless. In section (4.2), we considered a variety of different quench protocols with different starting points and measuring the response at different points in the critical region (again, chosen to be time t = 0), and the results are summarized in table 3.

JHEP11(2017)157
The simplest case to consider is the gapped-to-edge quench, described in section 4.2.1, where the quench of the Kitaev model starts in the gapped phase and at t = 0, the system is at the edge of the gapless region. The situation here is very similar to that of quenching across an isolated critical point, as in the Ising model. Indeed for small amplitude quenches, we observe three distinct scaling regimes: instantaneous, fast and slow regimes, which can be understood in terms of the continuum model (2.38). However, we must add that the latter is an anisotropic theory and the scaling behaviour does not match the scaling of a conventional fermionic field, which is also shown in table 3. Further, for large amplitude quenches, the fast scaling regime disappears. One difference in the Kitaev case is that for the large amplitude quenches, the scaling exponent in the slow regime is different from that in the small amplitude quenches.
We began to explore the extended critical region of the Kitaev model with the gappedto-interior and interior-to-interior quenches, which are discussed in sections 4.2.2 and 4.2.3. In both cases, the response is measured when the system is at an interior critical point, while in the first family, the quench starts in the gapped phase and in the second, the system is initially at an interior point. For all of these quenches, we again observe three distinct scaling regimes: slow, fast and instantaneous, with smooth transitions between them. The scaling behaviour in the fast regime depends on the details of the quench, as shown in table 3. However, for all of these quenches which traverse a finite part of the critical region, all exhibit the same scaling exponent in the slow quench regime. Namely, χχ diff ∼ (Jδt) −1/2 , which corresponds to the Kibble Zurek scaling of a relativistic fermion in one lower dimension. In fact, our preliminary results indicate this slow scaling behaviour extends to any quenches which traverse some distance in the interior region, e.g., quenches beginning at an interior point and ending at an edge point or vice versa.
This slow quench behaviour has been observed earlier for quenches of the Kitaev model with linear protocols which start at t = −∞, and the measurement is performed at t = ∞ [19,20]. In this case, the equivalent Landau-Zener problem has a simple solution and for slow quench rates, the excitation probability predominantly depends on k y alone. We do not have a clear explanation for how this behaviour emerges with our quench protocols at the moment. However, examining the integrand χχ ( k) in an expansion around the critical modes, we find that this quantity depends primarily on k y = k 1 − k 2 and is almost independent of k x = k 1 +k 2 . Hence this momentum dependence reflects the same behaviour found for the excitation probability for the linear quenches [19,20]. Further, as commented in section 4.2.2, while the critical models (2.27) for the interior points are not anisotropic, these quenches are inherently anisotropic because the quenched operator corresponds tō χγ x χ. Accordingly, we should think of coupling which we are varying in the quenches as the x-component of a vector. We expect that this anisotropy will play a central role in the explanation of the unusual scaling found in both the slow and fast quench regimes.
As commented above, one of the most interesting results here is that there is a fast scaling regime in many of our lattice quenches. For small amplitudes, our numerical results here match to the expectations of a continuum analysis for the Ising model and for the gapped-to-edge quenches in the Kitaev model. While we presently lack a theoretical understanding, it also appears that a fast scaling regime arises for quenches of the Kitaev JHEP11(2017)157 model which traverse a finite distance across the critical region, e.g., for the gapped-tointerior and interior-to-interior quenches. 22 Therefore our results indicate that for systems with a finite lattice spacing, there is quite generally a regime where the quench rates lie between the inverse lattice spacing and the physical mass scales, and where the fast scaling behaviour found previously only in continuum field theories holds. This opens up the interesting possibility that such scaling can be indeed observed in experiments.
In the cases where we can match the theoretical and numerical analysis for the fast scaling regime, the dimensionless couplings are small and the change in the couplings are small as well. It is only in this case that the crossover from the fast to the slow quench happens when δt is of the order of the inverse physical mass scale. This is expected, since small dimensionless couplings correspond to finite physical mass scales and the continuum fast quench scalings are expected when the quench rate is fast compared to the physical scales but slow compared to the UV cutoff scale. Indeed, in the case of the Ising model, when ∆g(t) ∼ O(1/|m in |) and in the case of the Kitaev quench from gapped to the edge for ∆J 3 (t) ∼ O(1/|m in |) ∼ O(1/J) we still have three regimes, but now the crossover between the fast and the slow regime is no longer at δt ∼ 1/|m in |. On the other hand the slow quench behaviour is insensitive to this, since this is the regime where the physical mass scale is much higher than the scale set by the quench rate.
There are a wide variety of different avenues to follow in extending our study of critical quench dynamics in lattice models. In particular, the three protocols introduced in section 4.2 to study the Kitaev model do not form an exhaustive list of the possibly interesting quench protocols in this lattice model. As alluded to in the above discussion, we have made some preliminary studies of interior-to-edge and edge-to-interior quenches. One feature that seems to extend to these protocols is the slow quench scaling: χχ diff ∼ (Jδt) −1/2 . It also appears that there is a fast scaling regime separating the instantaneous and slow quench regimes. As noted in section 2.2, there is a distinct "interior" edge point with g = 0. It would be interesting to probe this new critical theory with new quench protocols. 23 From the quantum information perspective, it will be important to understand quench dynamics on the open chain Cluster-Ising model which has non-trivial symmetry protected edge-states -see appendix C. In that context, it will be interesting to study the response of the string order parameter [32], in a similar manner in terms of free fermions.
As we commented above, the quenches in the Kitaev model which traverse the critical region are inherently anisotropic. While the critical theories corresponding to the interior critical points are not anisotropic, we are quenching the system with an anisotropic operator, i.e., χ † σ 3 χ ∼χγ x χ (in the continuum langauge). We also measure the response as the expectation value of this same operator. Hence it would be interesting to see if similar scaling laws hold in these quenches for other operators like χ † σ 1 χ, χ † σ 2 χ or χ † χ. Undoubtedly, this would give us new insights into the unusual scaling behaviour found for, e.g., the gapped-to-interior and interior-to-interior quenches. 22 Our preliminary results that a fast scaling regime also appears in interior-to-edge and edge-to-interior quenches of the Kitaev model. 23 Unfortunately, we found interior-to-"interior" edge quenches to be problematic, i.e., producing reliable numerics for these quenches seems challenging.

JHEP11(2017)157
Of course, we do not have a good theoretical understanding of many results for the quenches in the Kitaev model, particularly, for quenches that traverse a finite distance in the critical region. Certainly, this situation should be improved. We might note that this is required even for the small amplitude quenches with the gapped-to-edge protocols, where we heuristically applied eq. (1.13) with an effective spacetime dimension to predict the scaling exponent should be -1/2. While the fast quench scaling is well understood in relativistic theories [12,13], it is interesting to confirm that these ideas properly extend to non-relativisitic theories, and to semi-Dirac point appearing at the edge of critical region. Similarly, as discussed above, our quenches involving moving across the interior critical region are not really mass quenches, i.e., m is actually the x-component of a vector coupling. Hence it would also be interesting to extend the discussions in [12,13] to smooth fast quenches involving anisotropic operators.

JHEP11(2017)157
for timelike separations, zero otherwise. This is reflected in the range of the x integrals. In polar coordinates, for the spatial integral, we have where Ω d−1 = 2π d/2 /Γ(d/2) is the volume of a unit (d − 1)-sphere. Thus eq. (A.1) yields Now if we choose a simple impulse profile with F (x) = 1 in the range where it is nonvanishing, then When d = 2∆, the integral evaluates to In our discussion of quenches in the main text, the observable was measured at t → 0. Thus the leading order scaling with δt is given by O ∆ (0, 0) ∼ δt d−2∆ , as noted earlier in eq. (1.13). When d = 2∆, the upper limit of the t integral is divergent and we regulate this by shifting the upper limit of the integral, t → t + . The leading behaviour is then logarithmic, For the example of the (1+1)-dimensional Ising model in the continuum limit, we have a massless fermion and O =ψψ. Hence, we have ∆ = 1 and d = 2 and eq. (A.7) becomes δ ψ ψ (t = 0) = lim →0 2πδλ log (δt/ ) + · · · . (A.8)

B Saturation in the transverse field Ising model
One feature which distinguishes the lattice quenches from their counterparts in a continuum field theory is that for small enough δt, the expectation value of the quenched operator saturates in the lattice quenches. While this is expected since the lattice provides a natural cutoff given by the lattice spacing (which is hidden in the interaction strength J -see footnote 6), the details of this saturation are important to understand the different regimes which we are analyzing in this paper.
In this appendix, we will provide the details for understanding why the point at which saturation sets in is qualitatively different for small and large amplitudes in the Transverse

JHEP11(2017)157
Field Ising model. This difference was observed in the results in section 4.1 and already discussed there. The main result is that for small amplitudes the expectation value saturates to the instantaneous answer at a scale of Jδt ∼ 1, independent of the amplitude in , while for large amplitudes the saturation occurs for a value of Jδt proportional to the amplitude | in |. We expect that an analogous description will hold for quenches in the Kitaev model, which exhibits similar behaviour as described in section 4.2.
To understand the above difference, it is important to carefully account for the contributions of the various momentum modes in the two different cases. For that, we analyze the integrand in 2) which comes from combining eqs. (3.15) and (3.19) in eq. (3.16). Recall that all of the needed definitions are given in section 3. We note that for k = 0 and π, G(k) = 0 -see eqs. (3.5) and (3.6); consequently, the instantaneous energy levels at these momenta are decoupled. This in turn ensures that the probability of excitation, for any quench rate and amplitude for these modes, is either 0 (for k = π when the instantaneous levels do not cross) or 1 (for k = 0 when there is exact level crossing). Consequently, X(k = π) = 0 and X(k = 0) = 1. Now from plotting this integrand, it is straightforward to see that the behaviour is qualitatively different depending on whether the amplitude of the quench is small or large. In figure 12, we show this for two different cases, | in | = 0.1 and | in | = 100 for different values of δt. In the case of small amplitudes, the integrand is rapidly decaying and hence only low momenta contribute to the expectation value; while for large amplitudes, the integrand does not decay rapidly and so large momenta also give important contributions to the expectation value.
In figure 12, we also compare these profiles with the integrand for an 'instantaneous' quench, as in eq. (4.23). From eq. (3.19), this instantaneous integrand takes the simple form where we have combined the various definitions in section 3 to produce the final expression. Figure 12 shows that for either large or small amplitudes, the integrand X(k) quickly approaches X inst (k) as Jδt → 0. This will become a key fact in analyzing the modes contribution for small amplitudes, as we will see below. In passing, we also note that with in < 0, X inst (k = π) = 0 and X inst (k = 0) = 1, as was argued must be the case on general grounds above. Let us first analyze the small amplitude case in more detail. We first determine which modes are contributing significantly to the expectation value. Then since we are interested in understanding when the expectation value saturates, we will ask for which values of Jδt are all of these modes excited by the quench (as first discussed in section 4.1), and how this varies as we vary the initial amplitude. We approached these questions with a combination of numerical and analytic analyses. In figure 13, we study the decay of the integrand as a function of | in |. In particular, we find numerically the momentum k max where X(k) = 0.1, i.e., where the integrand decays to 10 percent of its maximum value. This gives us an estimate of the range of momenta for which the corresponding modes are contributing significantly to the expectation value (B.1). The result in figure (13) is robust: for any Jδt, there is a region for small enough | in | where k max ∝ | in | and moreover, the coefficient of proportionality is independent of Jδt. However, the coefficient obviously depends on the chosen threshold, i.e., X(k max ) = 0.1 in the present case.
We can also provide an analytic calculation which supports this same conclusion. First, using eqs. (3.10) and (3.11), we can expand the expression in eq. (B.2) for X(k) to find where X inst is given by eq. (B.3). Of course, this result is in agreement with the observation that X(k) quickly approaches X inst (k), made from the plots in figure 12. However, this result can be refined since if we carefully examine the expressions in eqs. (3.10) and (3.11), it is possible to see that every contribution of Jδt to the integrand (B.2) is accompanied either by a factor of the amplitude in or a momentum factor, sin(k/2). It will serve our purposes below to expand X(k) simultaneously for small x ≡ in Jδt and y ≡ sin(k/2) Jδt.
With these variables, we find that the previous expansion is replaced by X(k) = X inst (k) + O(y 2 , xy) . In the three cases, it is possible to observe a linear dependence of k max on | in | for small enough | in |. This is verified by the fit with a function k max = c | in | (solid blue line). Further, we find the coefficient c is the same in all three plots, i.e., it is independent of Jδt but depends, of course, on the choice of the threshold, being X(k max ) = 0.1 in the present case.
where implicitly we have assumed that x ∼ y in our expansion. 24 Of course, we may ensure that |x|, |y| 1 by simply taking Jδt 1 and with this choice, we recover eq. (B.4). However, the above expansion indicates for any Jδt 1/| in |, 1/k, the integrand is well approximated by the instantaneous integrand. Hence for small amplitudes, as considered here, we may still consider Jδt ∼ 1 as long as the relevant momenta are also small. This will now let us self-consistently prove that only small momenta contribute to the expectation value (B.1) for small amplitude quenches.
Above we argued that in the limit of small momenta and small amplitude, X(k) reduces to the instantaneous integrand. However, it remains to expand the expression in eq. (B.3) when taking k, | in | 1 while keeping k/| in | ∼ 1. This expansion yields where the O(k 2 ) is used above in a sense where the third term is O(k). Now, to determine the significant contributions, we want to see when the leading term in eq. (B.6) reaches a particular threshold γ in this limit. Hence we set where we assumed that in < 0 while γ > 0. We note that this O(1) term for X inst (k) in eq. (B.6) has a long tail. However, let us explicitly evaluate the pre-factor for γ = 0.1, 0.2, 0.5, and we find c ≡ 1 − γ 2 /γ = 9.95, 4.90, 1.73, respectively. So even with γ = 0.1, we have a consistent solution for k max . That is, our expansion assumed that the relevant k's were small and now we find that the maximum contributing momentum is proportional to | in | and so it is indeed small for the small amplitudes considered here. Moreover, as both the amplitudes and momenta are small, this result is still valid for

JHEP11(2017)157
Jδt ∼ 1, which will be needed below for our final result on the saturation point. It is also important that the pre-factor c in eq. (B.7) is independent of Jδt. Hence both our numerical and analytic calculations suggest that for small amplitude quenches, the modes which contribute significantly to the expectation value (B.1) lie in a narrow band: 0 ≤ k ≤ k max with k max = c | in | where c is some order one number (which is independent of Jδt). Now as discussed in section 4.1, we can expect that the expectation value saturates when all of these modes are in fact excited by the quench. In particular then, we must confirm that the last mode at k = k max is excited according to the Landau criterion (4.2). Substituting in k max = c | in | and keeping in mind that we are considering small amplitude quenches, i.e., | in | 1, we find small amplitude : which should then be larger than one to excite all of the modes contributing to χχ diff (t = 0). That is, we expect saturation for Jδt ≤ 1/(4c) or more simply Jδt 1. As stressed in the main text, this result shows that for the small amplitude quenches, the point where saturation sets in is independent of the initial amplitude. Note also that while the argument that the maximum momentum contribution is proportional to the initial amplitude is valid for larger Jδt (provided the amplitude is small enough), the Landau criterion in this case will show that it is still possible to keep exciting these modes and then we will see no saturation for larger Jδt.
The situation is qualitatively different for large amplitude quenches. As shown in figure 12, in this case, the profile of the integrand (B.2) is much broader and hence almost every mode makes a significant contribution to the expectation value. Hence we essentially have, k max π. Of course, the mode k exactly at π is not contributing since as noted above, X(k = π) = 0 in every quench, i.e., for any amplitude or quench rate. However, other modes near π will contribute significantly to the expectation value. 25 It is straightforward to evaluate the Landau criterion (4.2) for k max π, and as in eq. (4.3), we obtain large amplitude : As before, all the modes contributing to the expectation value (B.1) will be excited when the above expression is bigger than one. Hence we should expect saturation for Jδt | in |/8. This is the result reported in section 4.1, and hence the saturation point is proportional to the initial amplitude for large amplitude quenches. While all of the above analysis is particular to the Transverse Field Ising model, let us re-iterate that we expect an analogous description will hold for the results for the quenches in the Kitaev model described in section 4.2.

JHEP11(2017)157
C Cluster-Ising model The Hamiltonian of the Cluster-Ising model [32] on a one-dimensional bipartite lattice is where we have allowed for a time-dependent Ising coupling to consider quenches and τ (1,2,3) denote Pauli spin operators. For λ = 0, the ground state is protected by a Z 2 × Z 2 symmetry generated by i∈even τ z i × i ∈odd τ z i . This short-ranged entangled state has topological order (no long range order) on a open chain, realizing projective representations of the symmetry group via edge states. At λ = 1, the model has a quantum critical point and for λ > 1, the system is antiferromagnetic with long range order. At the quantum phase transition, the critical exponents are ν = z = 1, β = 3/8 and α = 0. We show below that on a closed chain with periodic boundary conditions, the model can be mapped to (1 + 1)-dimensional free fermions as in the Ising case.

JHEP11(2017)157
in the continuum limit. With the lattice spacing a, we introduce dimensionful momenta p, and a dimensionful mass m, p = 3k a and m(t) = λ(t) − 1 a , (C. 13) the Hamiltonian becomes in the a → 0 limit which is equivalent to three identical copies of eq. (2.13), describing the continuum theory for the Ising model.
Open Access. This article is distributed under the terms of the Creative Commons Attribution License (CC-BY 4.0), which permits any use, distribution and reproduction in any medium, provided the original author(s) and source are credited.