Quantum quench and thermalization to GGE in arbitrary dimensions and the odd-even effect

In many quantum quench experiments involving cold atom systems the post-quench system can be described by a quantum field theory of free scalars or fermions, typically in a box or in an external potential. We work with free scalars in arbitrary dimensions generalizing the techniques employed in our earlier work \cite{Mandal:2015kxi} in 1+1 dimensions. In this paper, we generalize to $d$ spatial dimensions for arbitrary $d$. The system is considered in a box much larger than any other scale of interest. We start with the ground state, or a squeezed state, with a high mass and suddenly quench the system to zero mass ("critical quench"). We explicitly compute time-dependence of local correlators and show that at long times they are described by a generalized Gibbs ensemble (GGE), which, in special cases, reduce to a thermal (Gibbs) ensemble. The equilibration of {\it local} correlators can be regarded as `subsystem thermalization' which we simply call 'thermalization' here (the notion of thermalization here also includes equlibration to GGE). The rate of approach to equilibrium is exponential or power law depending on whether $d$ is odd or even respectively. As in 1+1 dimensions, details of the quench protocol affect the long time behaviour; this underlines the importance of irrelevant operators at IR in non-equilibrium situations. We also discuss quenches from a high mass to a lower non-zero mass, and find that in this case the approach to equilibrium is given by a power law in time, for all spatial dimensions $d$, even or odd.


Introduction and Summary
Thermalization in integrable systems 5 has grown to be an important area of study in recent years [2][3][4][5][6][7][8][9][10][11][12][13]. Integrable systems possess an infinite number of conservation laws and are often described by free (quasi)particles; thus, to understand thermalization in integrable models, it is important to understand it in free field theories. Besides the theoretical motivation, it turns out that in many recent quantum quench experiments involving cold atom systems, the post-quench phase can be described by a quantum field theory of free scalars or fermions [14,15]. In a previous work [1] such systems were studied at length where the spatial dimension d of the system was d = 1 and the post-quench hamiltonian was taken to be critical. Both bosonic and fermionic systems were studied with explicit quench protocols. Among other things, it was found that (i) the post-quench state can be represented by a generalized Calabrese-Cardy state, characterized by an infinite number of W ∞ charges, as postulated in earlier works [16,17], (ii) correlators of local operators as well as the reduced density matrix of finite subsystems asymptotically approached that of a generalized Gibbs ensemble (GGE), (iii) the relaxation was exponential, characterized by a rate given by conformal and W ∞ algebraic properties of the operators in question. Furthermore, it was found that (iv) it is possible to start with specially tailored squeezed states so that the post-quench state is a simple Calabrase-Cardy state, characterized by only one conserved quantity-namely the energy; in this case the equilibrium is described by the familiar Gibbs ensemble and the relaxation rate is a simple function of the energy and the conformal dimension of the operator.
It is natural to ask, both from an experimental as well as a theoretical point of view: (a) which of the above results, if any, generalize to d spatial dimensions for arbitrary d?, and, (b) in particular, if there is an asymptotic equilibrium, what is the rate of aproach to the equilibrium? There are several reasons why the answer to these questions is not easy to guess. First of all, the proof of thermalization, for 1+1 dimensional conformal field theories in [16], relies on conformal maps which map a strip of the complex plane to an upper half plane. The quantitative aspects such as the relaxation rate, also arise from such conformal maps. In d + 1 dimensions, there is no conformal map which maps the Euclidean time circle (or interval) times R d to (a part of) R d+1 for d > 1 6 Furthermore, [16] extensively used properties of W ∞ algebra which are specific to 1+1 dimensions, to derive thermalization to GGE. In addition, in holographic contexts, the GGE in 1+1 dimensions finds an interpretation in terms of 3D higher spin black holes which have an infinite dimensional algebra of large diffeomorphisms (to be precise, W ∞ ). There is no obvious generalization of these to higher dimensions.
The explicit computations for free scalars and fermions, studied in [1], did not directly use the conformal map or the W ∞ algebra, but in terms of detail, used several features specific to 1+1 dimensions.
In this paper we consider quantum quench of the mass parameter of free scalar field theories in arbitrary dimensions (both to zero and non-zero final mass) and explicitly compute time-dependence of local correlators. For simplicity of the computation, we consider sudden quench for most part of the paper, although our method is applicable to more general quench protocols. A summary of our results is as follows: 1. General results for mass quench (m in → m out ) (a) If we start with the ground state of the initial Hamiltonian, the post-quench state turns out to be equivalent to a generalized Calabrese-Cardy (gCC) state, which is a boundary conformal state with a cut-off applied to each of an infinite family of commuting conserved charges. If we start with specially tailored squeezed states, the post-quench state turns out to be an ordinary Calabrese-Cardy (CC) state, namely a conformal boundary state with a cut-off only on the total energy of the system.
(b) We show that at long times (a) the correlators, starting from a CC state, are described by a thermal (Gibbs) ensemble, and (b) the correlators, starting from a gCC state are described by 5 We use the word "thermalization" here to mean asymptotic approach of a local correlator in a time-dependent pure state to its value in a generalized Gibbs ensemble (GGE). Sometimes the word "thermalization" is reserved only for cases where the asymptotia corresponds to a Gibbs ensemble; we will not make this distinction. Since we specifically refer to equilibration of local correlators, such equilibriation to GGE may be called "subsystem equilibration" or, in our generalized sense, "subsystem thermalization". 6 We thank Abhijit Gadde for a useful discussion in this regard. a generalized Gibbs ensemble (GGE).
(c) The GGE is characterized by chemical potentials which are determined by the infinite number of conserved charges of the initial gCC state. As found and remarked about in [1], this implies not only a retention of detailed memory of the initial state, but also an apparent contradiction of naive intuition that higher dimensional conserved charges (7) should not affect long time behaviour. This naive Wilsonian intuition, which perfectly works at low energies, does not work in situations where the initial state is a high energy state.
2. Quench to zero mass: m in → 0 (critical quench) (a) The rate of approach to equilibrium is exponential or power law depending on whether the number of spatial dimensions, d, is odd or even, respectively. This is one of the main results of the paper.
(b) In the case where the post-quench state is a CC state, there is a geometric understanding of thermalization in terms of the method of images, which also elucidates the odd-even effect.
(c) Unlike in 1+1 dimensions, the decay exponents are not uniquely defined by the conformal dimension of the operator. For example, in a given odd spatial dimension d > 1, where φφ and ∂ i φ∂ i φ both decay exponentially, with the same decay exponent.
(d) As in 1+1 dimensions, details of the quench protocol affect the long time behaviour. This underlines the importance of irrelevant operators at long times in non-equilibrium situations involving high energies.
3. Quench to small non-zero mass: In case the system is quenched from a high mass to a low non-zero mass, exact time-dependence of correlators can again be computed. These generically approach a GGE with an oscillatory power law behaviour, with a period of oscillation given by the final mass. There is no odd-even effect unlike for the massless case mentioned above.
We summarize the main points a bit more precisely in a tabular form below.  Table 1: Here d refers to the number of spatial dimensions. The table displays the structure of an equal time two-point function G(t). For d > 2, G(t) ≡ φ(x, t)φ(y, t) ; for d = 1, 2, m out = 0 this quantity diverges at long times, hence we display G(t) = ∂ t φ(x, t)∂ t φ(y, t) (in all cases, we drop the spatial dependence). The decay coefficients γ are given in the text (see e.g. Table 4 which lists various values of γ for d = 3: these values depend on the specific initial state); the exponential decay has a power law prefactor in case of ground state quench. The phase factors δ are given by d π/4. Note that the power laws in the massless case (for even dimensions) differ from those in the massive case.
The paper is organized as follows: In section 2 we discuss the main formalism of how to describe the post-quench state as a Bogoliubov transform of the out-vacuum, leading to its identification (with the help of Appendix A) as a generalized Calabrese-Cardy state (12). We discuss starting with ground states as well as squeezed states in the massive phase. Section 2.2.1 gives the formulae for the time-dependent two-point functions.
In section 3 we explicitly calculate the time-dependent part of the two-point function (additional material is provided in Appendices B and C). Tables 4, 5 and 6 show the decay of the time-dependent part for critical quench (m out = 0) for d = 1, 2, 3, 4; the results for general dimensions d are obtained by using a recursion relation. For non-critical quench (m out = 0), Tables 7, 8 and 9 present the time-depedent part of the two-point functions for d = 1, 2, 3, 4; the formula for the general dimension is given in (57). The overall behaviour is as indicated in Table 1 above.
The fact that the time-dependent parts of the two-point functions vanish at large times, already implies that these correlators (and hence all correlators, by the application of Wick's theorem), asymptotically equilibrate. In section 4, using results from Appendix D, we explicitly show that the time-independent parts of the two-point functions, to which they equilibrate, are given by two-point functions in a generalized Gibbs ensemble(GGE). This constitutes a proof of thermalization of arbitrary local correlators (see Section 4.2).
In section 5 we look at the geometrical interpretation of two-point functions in the CC state, which allows up to have a better understanding of the odd-even effect.
In section 6 we calculate the GGE correlator for purely spatial separation for the simple case of the thermal ensemble and a critical quench. The calculation has an interpretation in terms of a Kaluza-Klein reduction along the thermal circle. For high enough temperature, or equivalently large enough spatial sepation, only the Kaluza-Klein zero mode contribution survies, which has a power law behaviour. This result holds in any dimension, even or odd.
In section 7 we conclude with some comments and discussion of future directions.  . This is shown in the left panel. For simplicity, we will consider a sudden quench for most of the paper, in which case, the mass profile looks like the figure on the right.
The basic set-up is as follows. Consider a relativistic scalar free field theory in d spatial dimensions, with a time dependent mass (see fig. 1) In Fourier space, the equation of motion for the Fourier mode φ( k, t) (similar equation for φ * ( k, t)) is For every k, this can be identified with a Schrodinger problem on a line (coordinatized by y, say), with the identifications In the Schrodinger problem there are two equivalent bases of solutions: one which corresponds to particles coming in from the left: ψ E (y) = linear combination of u in (E, y), u * in (E, y), where u in (E, y) ∼ e iωiny as y → −∞, and another which corresponds to particles coming in from the right, u out (y) ∼ e −iωouty as y → +∞. Taking cue from this, we have two sets of normal mode expansions: Here. The two basis sets are of course linear combinations of each other Here α, β are the Bogoliubov coefficients. We take the mass profile to be (see Fig 1a) m Further we do all our calculations in the sudden limit ρ → ∞ 7 (see Fig 1b) . For the quench protocol eq.
(2), Bogoliubov coefficients are easy to compute explicitly [18]; in the sudden limit they become [1] α(k) = 1 2 In this limit the in-and out-waves also become especially simple The final mass is arbitrary but taken to be less than the inital mass i.e. m out < m in . We also study the theory when the final mass is taken to be zero, when the final theory is critical. As we will see the more interesting results are obtained in the critical quench.

Quantum quench from the ground state
A natural pre-quench initial state (just before t = 0) is the ground state |0 in of the initial Hamiltonian H in , which is also the zero-particle state defined by the oscillators a in |0 in = 0. Just after t = 0, the state remains |ψ(0) = |0 in (remember we are working in the sudden quench limit, therefore the state has no time to change). We then define |ψ(t) = e −iHt |ψ(0) where H := H out is the final Hamiltonian.
The 'in' ground state can be written in terms of the out ground state through a Bogoliubov transformation where γ(k) = β * (k)/α * (k) depends on the Bogoliubov coefficients and is only a function of k := | k|. 9 Further the 'in' ground state is also related to the Dirichlet boundary state [16] (see appendix A for more details) In general one can expand κ(k) for small k: With this, the expression for |0 in becomes where are consereved charges (they obviously all commute with the 'out' Hamiltonian H = k | k| 2 + m 2 out a † out ( k)a out ( k)). 10 It is easy to explicitly compute the κ i coefficients by using the definition of κ(k) in terms of the Bogoliubov coefficients which are given in (3): Note that henceforth we will call m in =: m for simplicity. It is assumed that m > m out . From the above equation, the small k expansion can be easily found. For m out = 0, the expansion contains only odd powers of k, thus the even κ i 's are non-zero, e.g.
By contrast, for m out = 0, only even powers of k survive, leading to odd κ i 's, e.g.
Relation to a generalized Calabrese-Cardy (gCC) state In case of critical quench (m out = 0), the post-quench dynamics is conformal. For critical quenches leading to a generic (non-integrable) conformal field theory, Cardy [CC:2004-2005] postulated the following form for the post-quench wavefunction 11 |CC = exp[−κH]|Bd (11) 9 Note that, because we focus on homogenous quench protocols in this paper, we have rotational invariance in the d spatial dimensions. Hence the Bogoliubov coefficients are a function of k = | k|. 10 For free scalars, the number operators N (k) = a † out ( k)aout( k) are themselves conserved, and provide an alternative basis of the algebra of conserved charges. 11 One assumes here that the quench takes a finite time to end, say t = t 0 . In this paper, we are interested in a sudden quench which takes place at t = 0; thus t 0 = 0.
where, κ is given by the inverse of the mass gap characterizing the initial state (before the quench). The state |Bd is a "boundary state" representing a conformal state subject to an appropriate boundary condition at imaginary time t = −iκ.
In case the post-quench conformal theory is integrable, characterized by an infinite number of charges Q i , a more appropriate ansatz for the post-quench state is the gCC (generalized Cardy-Calabrese) state [1,16] The state (6) is clearly a special case of such a state, which was first found in [1] in 1+1 dimensions. Note that the Dirichlet state |D is a particular example of a boundary state. As mentioned above (9), here only the even charges Q 2n , n > 0 are non-zero; note that Q 2 = H, as can be seen from the definition (7).
In case of a massive quench, the state (6) is of the form (12), where the conformal boundary state |Bd is replaced by the Dirichlet boundary state |D . Although the post-quench theory is not conformal in this case, we will continue to use the name gCC for the state (6). In fact we will continue to use this nomenclature also for states like (20) obtained by noncritical quench from squeezed states.

Two-point functions
We will be interested in computing quantities like The operators appearing on the RHS O(x, t) = e iHoutt O(x)e −iHoutt are in the Heisenberg picture defined with the Hamiltonian H out . This is applicable for t ≥ 0. For time evolution to t < 0 we must use the Hamiltonian H in . With this understanding it is clear that the above definition comes with the prescription t 1 , t 2 ≥ 0. The relation to time-ordered correlator is noted below (16). This definition includes in particular, the equal-time correlator (ETC) In our theory these are related to the ground-state two-point function Using the expressions (3) and (4), the ground-state two-point function becomes As mentioned above this correlator is not time-ordered. The time-ordered correlator can be obtained from this by replacing t 1 − t 2 by |t 1 − t 2 | but leaving t 1 + t 2 unaltered. As a consequence, the ETC's considered here are already time-ordered. The expression above simplifies for t 1 = t 2 (ETC): We also look at the 2-point function of the ∂ t φ operator. Note that this correlator is directly obtainable from the unequal time two-point function φφ (16) by applying the operator ∂ t1 ∂ t2 | t1=t2=t . Since the unequal time correlator is of the form f (t 1 − t 2 ) + g(t 1 + t 2 ), it follows that 0 in |∂ t1 φ( . Hence the equal-time two-point function of ∂ t φ is of the form −f (0)+g (2t); in particular, the time-dependent part of this correlator is ∂ 2 t applied to the time-dependent part of (17). By explicit calculation, one finds the expression which indeed follows the general form −f (0) + g (2t) mentioned above.
A primary motivation for considering two-point functions of the derivative operators ∂ t φ∂ t φ is as follows.
In low dimensions the time-dependent part g(2t) of the equal time correlator ∂ t φ∂ t φ (17) grows in time which masks the transients that signal thermalization (it grows linearly in time in d = 1, and logarithmically in d = 2). The extra time derivatives, leading to g (2t), get rid of these divergences.
For the same reason, we may also consider the correlators of the spatial derivatives ∂φ, The extra derivatives ensure decay with increasing | x| (in d = 1, 2 φφ grows, respectively, linearly and logarithmically with | x|). It is enough for this purpose to focus on We will discuss more details of two-point functions of these derivative operators and their relations in Section 2.2.1 below.
2.2 More general quantum quench: from squeezed states As we saw above, the post quench state (6) is built out of infinite number of chemical potentials acting on the Dirichlet boundary state. However certain special gCC states can be obtained if the initial state are chosen to be specific squeezed states of the pre-quench Hamiltonian This is just the Bogoliubov transformation of |0 in . As the |0 in state can itself be written as a Bogoliubov transform of |0 out , the post quench state |f in is a composite Bogoliubov transform of |0 out : where in the first line γ ef f denotes the composite γ: and in the second line we have the relation κ ef f = − 1 2 log(−γ ef f ), following similar arguments as for the ground state. Here α, β are given by (3). The state (20) is of the form of a gCC state. Can we get any gCC (with a given κ(k)) starting from a suitably chosen squeezing function f ?
The answer is clearly yes. By solving (21) for f , we can clearly find an f for any γ ef f . If we wish to generate a given κ ef f (k) ≡ κ(k) in (20), we must choose γ ef f (k) = − exp[−2κ(k)]. This gives us This proves that we can prepare any gCC state, with a given κ(k), from squeezed states: Critical quench from specific squeezed states Let us consider the case of critical quench m out = 0. We look at the following special states. From (22) it is clear that the choice f ( k) = f 4 (k) = 1 − 2k √ k 2 +m 2 tanh(κ2k+κ4k 3 )+k leads to a state with only two non-zero tunable parameters κ 2 and κ 4 we call this the gCC 4 state because it is characterized by only two charges Q 2 and Q 4 ; we will find that the equilibrium state describing asymptotic correlators in |gCC 4 is a grand canonical ensemble chracterized by a temperature and one chemical potential. Further if κ 4 is also zero, i.e. f ( k) = f 2 (k) = 1 − 2k √ k 2 +m 2 tanh(κ2k)+k then we end up with the CC state Note that κ 2 here is not related to the mass parameter m unlike in (9) where κ 2 = 1/m.

The squeezed state 2-point function
The 2-point function in the squeezed state is the same as in the ground state with α and β replaced by α ef f and β ef f above. The 2-point function in the general squeezed state (23) is given by When t 1 = t 2 = t this is For the equal-time two-point function of ∂ t φ, we apply to (26) the rule −f (0) + g (2t) mentioned above (18), which gives us For the two-point function of spatial derivtives, these derives act only on the exponential term, leading to Note that the time-dependent parts of the above correlators, calling them ...

(td)
gCC , satisfy the relation

gCC correlators encompass all cases
As is clear from the foregoing discussion, all the particular post-quench states |ψ(0) mentioned above, specifically |CC , |gCC 4 , and |0 in , are examples of |gCC states, with specific κ(k) as in Table 2:  (8). This table paves way for a uniform discussion all post-quench states discussed in this paper.
Consequently, the above formula (27), and the other related two-point functions, give the corresponding two-point functions in all these different cases. In particular, using the form of κ(k) (8) for ground state quench, one can recover all results of Section 2.1.1 from those of Section 2.2.1, and in particular recover (17) from (27) (after all, the ground state is a special squeezed state!). This implies, in principle, that if we derive some results on thermalization and relaxation for the general gCC states, using (27), we need not consider the above particular cases separately. However because each of these cases is significant on its own, we find it useful to often state the results separately.
where r = x 1 − x 2 . For correlators of the φφ type, To proceed further, we write the vector k in polar coordinates (k, θ, ϕ 1 , ..., ϕ d−2 ), oriented such that r points to the 'north pole' which ensures k. r = kr cos θ; the ϕ i parameterize a unit sphere S d−2 (we assume d > 1 here). With this, the above integral becomes where Angular integrals The angular integration in (34) can be exactly done, which gives The regularized Hypergeometric function 0F1 is some combination of trigonometric functions sin(kr) and cos(kr) in odd d while it is some Bessel function in even d. The specific forms of A d (k, r) for various dimensions d are tabulated below in Table 3: k(sin(kr)−kr cos(kr)) 4π 3 r 3 Table 3: The angular integral A d (k, r), (see eq.(35)) in the first few dimensions.
Recursion relations for 2-point functions The structure of the integrals (31) occuring in the various correlators defined above makes it possible to connect them across dimensions through recursion relations(see Appendix B for details) 13 . From a correlator in d + 1 spacetime dimensions, one can obtain the correlator in (d + 2) + 1 dimensions by acting with a particular differential operator. This relation holds for both even and odd d > 1 (here 12 We added d = 1 here for uniformity, using ∞ −∞ dk 2π e ikx f (|k|) = ∞ 0 d|k| 1 π cos(|k||x|)f (|k|). 13 Similar recursion relations apppear in a somewhat different context in [19] where ... (d) denotes an expectation value, in a d + 1 dimensional theory, in a general post-quench state.
For the special case of equal time correlators which is our main focus, the recursion relation becomes A similar relation holds for the equilibrium correlators. In an arbitrary GGE, we have the following recursion relation (t here is same as The 1+1 dimensional 2-pt function is not connected to its 3+1 dimensional counterpart in the above fashion, since in one dimensional space there is no angular coordinate which is crucial for this recursion relation to work (see Appendix B). Nevertheless, there is a different, specific relation connecting the two: which follows straightforwardly, by inspection.
Time-dependence As remarked before, there is a time-dependent part and a time-independent part of the generic equal time two-point function (31). We will come back to an analysis of the time-independent part in Section 4 where we show that this part exactly matches the corresponding two-point function in a GGE (generalized Gibbs ensemble).
In the remainder of the section, we will therefore consider only the time-dependent, non-equilibrium, part. We will show that this part decays to zero at large times, which amounts to a proof of thermalization, and we will derive the form of the decay. It will also be conveninent, for this purpose, to differentiate between massive quench m out = 0 and critical quench m out = 0. This is because in the latter case, the postquench theory is conformal and we expect, and indeed verify, that this theory will have special proprties as compared to the massive case. We will comment later (see Section 3.4) why we cannot obtain results for the critical quench by simply taking the m out → 0 limit of the power series expansion for m out = 0.

Critical Quench Correlators
As mentioned above, we will mainly focus on the time-dependent parts of the 2-point functions, while the time-independent parts will be considered in Section 4.

Time-dependent part of the φφ correlator
Let us put m out = 0. To keep the discussion general, let us consider the general gCC correlator (27) (as emphasized in Section 2.3, we can infer about the specific cases of ground states and the squeezed states from here). In the notation of (33), we now have The functions κ(k) are tabulated for the special cases of interest in Table 2, where for the ground state, we now have which has an expansion in odd powers of k, as mentioned above (9). Following this, we will focus, in the case of critical quench, only on gCC states with κ(k) which have an expansion in odd powers of k. The κ(k) for the CC and gCC 4 states, of course, satisfy this by definition.
With this assumption, the functions F (k), G(k) are even in k.

Large time behaviour for odd d
Note from the Table 3, A d k, r that for odd d are even functions of k. Using the evenness of F (k), G(k), we find that the integrand in (33) is even, which allows us to extend the integration contour to the entire real line (for even d, this step is not allowed, hence we will do something else). If one can close the contour on the upper or the lower half plane the integral can be evaluated by the method of residues. To see how it works, let us take the example of d = 3 (the procedure described below is similar to the case of d = 1 [1]; and in special cases, where exact results are known in d = 1, the 3-dimensional results can be derived by using (39)). For more details see Appendix C.2.
For the critical gCC correlator at d = 3, we get Let us focus on the time-dependent terms B and C (the time-independent term A is described later on, and in Appendix C.2). We will also assume large t so that 2t > r. Then from the asymptotic behaviour of functions κ(k) in the Table 2, one can see that the functions κ(k) grow as k → ∞, in a way that ensures that the integrals are convergent, and for large t (so that 2t > r) it is possible to close the contour for the integral B on the UHP and for C on the LHP. At large times t, the slowest transient is given by the pole of the integrand singularity nearest to the origin (in the UHP for B and LHP for C; the pole at the origin itself cancels between the terms B and C; see Appendix C.2, especially Figure 6). It is easy to see that for a gCC with a finite number of chemical potentials, the nearest such singularity is a pole. This is given by the relevant zero of sinh(2κ(k)), or, eqivalently by the corresponding root of The root discussed above turns out to be of the form as k = k 1 ± ik 0 , k 0 > 0; here we choose the + and − signs respectively for the B and C terms; the value of k 1 may be zero, positive or negative. By applying residue calculus, we find that the slowest decay at large t of the φφ correlator is given by For non-zero k 1 this comes multiplied by an oscillatory term, of time period 2π/k 1 . In case of the CC state, the exact time-dependence can be calculated by summing over the poles (see Appendix C.2 for details). The exact two-point function for the CC state as well as for the ground state can also be derived from the their 1+1 dimensional counterparts found in [1] using the recursion relation (39) (see Appendix C.2). The net result of all this is that at large times the time-dependent part of the critical φφ fall off exponentially, as in (44); these are tabulated in Table 4.
In case of a gCC with an infinite number of chemical potentials, such as when the gCC is the ground state, the poles can merge into a branch cut (see also a similar discussion in [1]). The structure (44) is retained even here, with additional power-law prefactors. Details of the above procedure for the CC, gCC 4 and the ground state, are worked out in Appendix C.2. Note that although for specificity we had taken the example of d = 3 here, the arguments go through for d = 5 or any higher odd space dimensions.
Since we are interested in the leading approach to equilibrium for simplicity, in the tables we only display the correlators in the limit mt 1 mr or equivalently t κ2 1 r κ2 . The ellipsis in the gCC 4 correlator represents O(κ 4 2 ) terms.
Ground state Table 4: Note that for even d, e.g. d = 4, the leading large t dependence is the same in all three types of post-quench state |ψ(0) , with the identification m = 1/κ 2 (see (9)). Further the gCC answer tells us that since the higher κ i 's for i > 2 show up in subleading terms, the higher conserved charges do not affect the leading transient towards thermalization. However, as we will see later, the time-independent part of these two-point functions reflects an equilibrium given by a GGE which is characterized by all the κ i parameters.
General odd d result: From the recursion relation (36) it follows that for all odd d ≥ 3, the slowest decay at large time of the time-dependent part of the φφ correlator is given by replacing (44) with where k 0 has the same value discussed above and is independent of dimension d. The constants b d , c d depend on the dimension d. The expression above may come with a possible oscillatory part as explained below (44).

Large time behaviour for even d
Here the angular function, A d (k, r) for even d are odd functions of k. Therefore, unlike in odd d, we cannot extend the integration contour to the entire real line and apply the method of residues. Fortunately there does exist a different approach in this case. Let us take the example of d = 4 to illustrate the method. As in the case of odd d, we allow the most general κ(k) = ∞ i=1 κ 2i k 2i−1 which grows at infinity such that the integral is convergent. The critical gCC correlator at d = 4 is In the notation of equation 33, G(k) = − csch(2κ(k))/(2k). As before we are interested in large t. Consider the following scaling and rewrite the time-dependent part in terms of this new dimensionless momenta p In the second line we have rotated the contour anti-clockwise for e ip and clockwise for e −ip . We can do this because the integrand J 1 (kr) csch(2κ 2 k) vanishes on the arc at infinity due to the exponential damping in the csch(2κ(k)) as k → ∞. In the last line we have Taylor expanded the integrand at large t. We have also introduced the integrals I(n) defined by where in the second line we have introduced p ± = p(1 ± i ). Note that in spite of the growth of the I(n) with n, because of the even more strongly convergent power series expansion for the product of the Bessel and cosech functions, the equation (47) has a convergent expansion in 1/t 2 for large enough t.
Notice that κ 4 shows up only from the O(t −4 ) onwards. Similarly κ 6 starts appearing only from the O t −6 term; explicitly that term is Having calculated I(n)'s, we find that the time-dependent part of the gCC 2-point function goes as The effect of κ 4 is suppressed by O(t −2 ) relative to the leading term and therefore it makes no difference to the leading power law fall off in the large t limit. The same is true for higher κ n 's, as they are even more suppressed.
The equation (49) can be verified against explicit numerical integration of the momentum integral. In practical terms, the subleading term, of order 1/t 4 , is possible to read off by first also calculate the timedependent part in the gCC 4 state numerically and see exact agreement in subleading O(t −4 ) behaviour after subtracting off the κ 2 -dependent part ( fig. 2). The d = 2 case is discussed separately below (see table 6) for a summary and also appendix C.1 for details).
One might wonder if the same reasoning would have given us the answer in odd d, then we would not have to do the contour integrals there. In d = 3 for example, the time-dependent part of the 2-point function has even powers of p (eq. 46) What this means is that either the answer is zero or something non-perturbative, i.e. a function with no Taylor expansion at t = ∞. Since we already know the answer is an exponentially decaying function e −πt/κ2 , we know it is the latter. The above argument elucidates the odd-even difference in the approach to thermalization. We give another geometric understanding in the next section 5.

∂
Ground state in the left column and d = 4 in the right column. Here r = | x 1 − x 2 |. |ψ(0) denotes the type of post-quench state.
As argued in section 2.1.1, the 2-point function of ∂φ is not independent of the 2-point function of φ. As explained in that section, in d = 1 and 2, our primary interest is in ∂φ∂φ . For d = 1 [1] this quantity for the general gCC state can be obtained by using the equations (28), (33) and (34): This expression can be evaluated in the same manner as described in Section 3.2.2. Similarly for d = 2, the time dependent part of the gCC 2-point function can be evaluated by again using equations (28), (33) and (34): This expression can be evaluated in the same manner as described in Section 3.2.3. We tabulate the results below 6.
Ground state in the left column and d = 2 in the right column. Here |ψ(0) denotes the type of post-quench state.

Massive Quench Correlators
In the massive or non-critical quench the final mass m out is non-zero. Here we incorporate stationary phase approximation to perform the momentum integrals.

φφ Correlator
In the notation of (31), the time dependent part of 2-point function in the gCC state (27) can be written as for an arbitrary gCC state, while G(k) = (m 2 out −m 2 ) 4 √ k 2 +m 2 (k 2 +m 2 out ) for the ground state.
We consider the gCC states defined by The reason for choosing only even powers of k is motivated by only even powers appearing in eq.(10). Let k = k 0 be a stationary point of the integral I 1 (for which w (k 0 ) > 0) having the general form where we Wick rotate t → it in the first line (when w (k 0 ) < 0 which is the case in I 2 , we rotate in the opposite direction i.e, t → −it). We expand both f and w around k = k 0 and introduce the variable y = √ t(k − k 0 ) in the second line. We are using compact notation to write f 0 = f (k 0 ), f 0 = f (k = k 0 ), etc. and similarly for w. In the third line plug in w(k) = √ k 2 + 1 and finally in the last line we have brought down all terms in the exponent except for y 2 . The factor m d−1 out is pulled out to make rest of the integral dimensionless. So every quantity inside the integral is multiplied by approriate power of m out to make it dimensionless.
t → m out t; r → m out r; k → k/m out ; m → m/m out . G(k) is an even function of k, while A d (k, r) is an even function only in odd dimensions but an odd function in even dimensions. This means for even d, the integration is neccesarily to be done for k belonging to [0, ∞). This has important consequences for us because the stationary point which is the solution of the equation w (k 0 ) = 0 is exactly k 0 = 0. So for even d, the odd powers of y will also contribute. Since the most contribution comes near k = 0, we look at the behaviour of f (k) near the origin Therefore the leading contibution comes from the (d − 1)th derivative i.e. f The third term in the first line of (56) comes from the product of the leading k d−1 term in expansion of f and w (4) y 4 /t brought down from the exponential. In the second line we have performed the integral using the general formula ∞ 0 dy y d e −y 2 = 1 2 Γ( d+1 2 ) and finally we Wick rotate back to t → −it. I 2 is calculated in a similar manner with the only difference being that in the final answer replace i → −i (due Table 7: Time-dependent part of 2-point function gCC|φ( x 1 , t)φ( x 2 , t)|gCC in the gCC state defined by κ(k) = κ 1 + κ 3 k 2 + κ 5 k 4 + · · · . Here m out is the final mass and r = | x 1 − x 2 |.
to Wick rotation in opposite direction). I 1 and I 2 thus combine nicely to give a cosine. Reintroducing the dimensions, the time-dependent part of 2-point function in general is given by The above expression tells us that the leading transient only depends on κ 1 . It neither cares about the the higher conserved charges nor the separation r which only show at subleading order. We tabulate the leading approach to equilibrium of (27) in Table 7. We also tabulate the 2-point function in the ground state state (17) in Table 8.
It is a bit surprising that the leading power of t here is also same as before, but it is easy to understand why. For a moment let us rewrite 2t = t 1 + t 2 . Then take derivatives on φφ correlator (in 1+1 dim. for Table 9: Time-dependent part of 2-point function 0 in |∂ t φ( x 1 , t)∂ t φ( x 2 , t)|0 in on the left and 0 in |∂ i φ( x 1 , t)∂ i φ( x 2 , t)|0 in on the right, both for massive quench. Here r = | x 1 − x 2 |. example) wrt t 1 and t 2 , and then restore t 1 = t 2 = t.
Clearly the leading power is unchanged when both the derivatives act on the cosine. However the ∂ i φ 2-point function in any gCC state, following the same argument as given after eq.(54), behaves as Explicit expressions in various dimensions for these two-point functions for the ground state are tabulated in Table 9.

Some comments on 2-point functions in critical vs. massive quench
From the results of the previous subsections it seems that the critical correlators behave in a qualitatively different fashion from the those for the massive case, especially in odd spatial dimensions. This seems to indicate an apparent discontinuity as m out → 0. This is surprising since the master formula (40) for the critical correlator was obtained from that of the massive case (17) by putting m out = 0! Actually there is no true discontinuity at m out = 0 in the exact time-dependent correlators. The apparent discontinuity emerges because the limits m out → 0 and t → ∞ do not commute. This can be understood by considering the scales involved in the theory. For any given mass m out , however small, at very large times we will have t 1/m out (see Fig. 3 (a)). In other words, the dimensionless quantityt = m out t 1. Figure 3: Relevant scales in the theory for the present consideration. Here r/t is considered fixed as the large t limit is taken. In (a) m out is a finite non-zero quantity. t → ∞ limit means m out t → ∞. This is the case we discuss in the analysis of large time behaviour for quench to non-zero mass. In (b) the limit m out → 0 limit is taken first, so that m out t → 0; the large t limit in this case can be defined as m t → ∞. We do not consider this limit in the paper except when m out = 0 (this limit is discussed a bit more at the end of the subsection).
More mathematically speaking, the two-point functions involve scaling functions of the form F (mt, m out t) (omitting the dimensionless variable r/t). For the massive two-point functions, late times imply that both scaled variables are taken to infinity; in particular m out t → ∞. Indeed, the late time massive correlators in the previous section are given in terms of an expansion in 1/(m out t); as explained above, this expansion cannot capture the limit m out t → 0. Of course, if an exact evaluation of the relevant integral is possible, this limit can be taken. 14 Such a function would then interpolate between the critical and non-critical behaviour. The apparent discontinuity at m out = 0 is only an artifact of the large time asymptotic expansion.

The Generalized Gibbs Ensemble (GGE)
We have seen above that the time-dependent part of the two-point functions (the part containing the cosine term in (31)) vanishes at large times. The remaining term, therefore, represents the asymptotic value of the two-point function. We will now show that the asymptotic value corresponds to that of the two-point function in a generalized Gibbs ensemble (GGE). Before that, we will first recall some relevant facts about a GGE.

GGE
The generalized Gibbs ensemble (GGE) can be regarded as a certain grand canonical ensemble with an infinite number of chemical potentials [2,3,8] which are appropriate for the description of equilibrium for an integrable model. An integrable model has an infinite number of commutative conserved chargesQ i . A GGE is described by a density matrix Let us consider a free scalar field φ(x, t), described by the standard mode expansion φ( k, t) = a( k, t)u( k, t) + cc, u(k, t) = e −iω(k)t
This system is clearly integrable; a standard basis of the commuting conserved charges is (cf. (7))Q i = k |k| i−1N (k), i = 1, 2, .... Another simple basis of charges are the set of occupation numbersN (k) themselves. With this choice the GGE is described by the density matrix The GGE partition is easy to evaluate (most calculations in this subsection are presented in greater detail in Appendix D): It is also easy to evaluate the average value of the number operator With the above ingredients, we can now compute the equal time two-point function in the GGE (for more details and also the case of unequal time correlator, see Appendix D): Note that for massless scalars, H = k |k|a † (k)a(k) = Q 2 . Also the total occupation number is always N = k a † (k)a(k) = Q 1 . Hence, we can write (60) as Thus, a standard thermal (Gibbs) ensemble is a GGE with all µ i =2 = 0.

Equilibration to GGE
Note that the two-point function in the most general gCC state (27) is of the form In the preceding sections we have shown that the time-dependent part decays at late times exponentially or by a power law. Hence the gCC two-point function (65) asymptotically approaches the GGE two-point function (63) provided we identify Indeed, these relations ensure To prove this, note that by identifying |gCC with |f in and using the top line of (20), we get that the right hand side of the above equation is given by (1/γ 2 ef f − 1) −1 . By using γ ef f (k) = exp[−2κ(k)] (see immediately above (22)), the RHS becomes (e 2κ(k) − 1) −1 which agrees with equation ((62)) with the identification (66). Q.E.D.
Thus, we have established the following statement of equilibration: where the GGE is defined in terms of the gCC by the basic relation (67). Since for a free field theory, all correlation functions can be essentially built from the two-point function, the above statement of equilibration is true for all correlators. This constitutes a proof of the GGE version of quantum ergodic hypothesis in the quench models considered in this paper: Eq. (66) gives us the important relation between the parameters of the equilibrium ensemble (GGE) with those of the initial state. In the massless case, in terms of the notation of (64), we can rewrite (66) as It is clear that the CC state (25) equilibrates to a standard Gibbs ensemble with µ i =2 = 0.

Geometrical interpretation of the correlators in the CC state
Here we follow [20] to show how the two-point function in a CC state can be identified with an image sum in a slab geometry R d × interval. This would help us better understand the origin of the difference in thermalization between odd and even dimensions (in short, the odd-even effect). The slab propagator is defined as The expression (70) is obtained by the method of images as shown in the diagram below (see figure 4). The geometry reflects the fact that a CC state of the form e −κ2H |Bd can be regarded as a Euclidean time evolution ∆τ = κ 2 , here denoted by L/2. The suffices th and dif f will be explained below. Figure 4: The two-point function between φ(τ 1 ) and φ(τ 2 ) can be regarded as the electrostatic potential between two positive charges of unit magnitude [20]. This follows from the fact that the massless Klein-Gordon equation in Euclidean space is simply the Laplace equation. Here we consider φ(τ 1 ) as a source charge whose electrostatic potential in the slab geometry can be found out by the method of images by regarding the slab as two parallel mirrors. In the right panel we depict by green(red) circles the images obtained by an odd(even) number of reflections of the source. The electrostatic potential in the slab can now be regarded as a sum of that between the probe charge placed at τ 2 and each of the images.
It is easy to verify that (70) reproduces the two-point function in a CC state. All the four terms in (70) can be easily summed to obtain Now analytically continuing to real time by τ → it and identifying κ 2 = L/2 we find the slab progator is the same as (26) with κ 2 = L/2 and all other κ i = 0 This completes the explicit check that the slab propagator (70) indeed represents the two-point function in the CC state.

Comments on approach to thermalization and the odd-even effect
The understanding of thermalization from the above geometrical picture goes as follows. First, note that the terms in S th correspond to Euclidean R d+1 2-point function between the test charge with the '+' images above and below it. These are precisely the terms which would appear in a cylindrical geometry, whch represents a thermal 2-point function with Euclidean time period β = 2L 2 = 4κ 2 ; hence the suffix th . The terms in S diff , on the other hand, correspond to R d+1 2-point functions between the test charge with the '-' images above and below it; below we will find that these terms vanish at long times, proving that the CC two-point function asymptotically approaches the thermal 2-point function S th , β = 4κ 2 , which verifies the relation (69).
It is clear from the above that S diff represents the difference between CC correlator and the thermal correlator. For thermalization of the auto-correlator ( x = y = r, t 1 = 0), we need to consider the following difference and apply it to the 2 sums above to get and therefore Analytically continuing to real time (τ 1 → it and τ 2 → it), for t → ∞ 15 we find (after putting L = 2κ 2 ) G CC (t, r; t, r) − G th (t, r; t, r) β = 1 128π 2 κ 2 t 2 + κ 2 15 Note that although |τ | ≤ L, there is no such restriction on t which can be taken to infinity.
This gives us the power law decay in even dimensions. Also this is precisely the expression for the 2-point function of φ in the CC state (κ 2n = 0, n > 1) for r = 0. We match the Euler-Maclaurin result with the method used before to calculate eq. (49) to O(t −20 ) and find exact agreement. Notice that there is a coincident singularity because of t 1 =t 2 , but that would show up in 'thermal' part G th , so we do not need to worry about it here. Now lets see what happens in odd dimensions. d = 3 : Here we have The Euler-Maclaurin sum formula applied here results in Together with the n = 0 term We expect the vanishing result found above to continue to all orders in the 1/t expansion (we have explicitly checked this fact till O(t −20 ). Indeed this had to be the case because we know the answer goes as e −πt/κ2 (see (103)) which does not have a Taylor expansion in 1/t around 0.
It is not difficult to generalize the above results to arbitrary dimensions d > 2 for example by using the recursion relations (see Appendix B). The summary is that There is a pictorial way to understand the decay of the time-dependent part of the CC correlator. In the slab diagram imagine the analytic continuation of time (τ j → τ j + it j ) with real time axis perpendicular to the plane of the paper (see figure 5). The spatial dimensions are all in the horizontal direction. Note that upon analytic continuation to real time, the 'even' and 'odd' images move in opposite directions. Also notice that S th involves the correlator between the probe(black) and the even images of the source (marked red); the way to see this is to note that such correlators all involve the time difference t 1 − t 2 as expected from (70). Similarly while S dif f involves the correlator between the probe(black) and the odd images of the source (marked green); the way to see this is to note that such correlators all involve the time difference t 1 + t 2 again as expected from (70). When t 1 and t 2 are increased by the same amount, the Euclidean separations involved in S th remain the same whereas they go to infinity in case of S dif f . This is the reason why S dif f decays in time in this limit and the slab correlator approaches the thermal value. 16 In the above we have considered separation between the source and the probe only in time but the entire argument goes through for fixed, finite spatial separation r as we take t → ∞.

The thermal auto-correlator
In the previous subsection we showed how the dfference between the CC correlator and the thermal autocorrelator vanishes as t 1 , t 2 → ∞. In this section we will compute the thermal auto-correlator itself. We Figure 5: In the above diagram we consider φ(τ 1 ) as the source and φ(τ 2 ) as a probe. As we analytically continue τ j → τ j + it j , and t 1 and t 2 are increased by the same amount, the red images (even reflections of the source) move in the same direction as the probe, while the green images (odd reflections) move in the opposite direction. For large t 1 and t 2 with t 1 − t 2 held fixed, the probe φ(τ 2 + it 2 ) sees the green images at progressively farther distances and hence their correlators decay. Thus in this limit the slab correlator approaches the thermal correlator.
will find that this too shows an odd-even difference between dimensions as one takes the limit of large The quantity to compute is Note that here τ = τ 1 + τ 2 . As in the previous subsection lets analyze this in d = 4 and d = 3 separately.

d = 4 :
Here Applying the Euler-Maclaurin formula After analytically continuing τ = τ 1 + τ 2 → it 1 + it 2 = it, together with the n = 0 term, for large t This matches the precisely the answer (117), obtained from the integral formula (120).

d = 3 :
Here Applying the Euler-Maclaurin formula After analytically continuing τ = τ 1 + τ 2 → it 1 + it 2 = it, for large t, together with the n = 0 term In fact, it turns out that the Taylor coefficients vanish to all orders. This is consistent with the behaviour at d = 3, G th ∼ e −γt which is independently obtained in Appendix C.2 see equation (111). Note that e −γt has an essential singularity at 1/t = 0 and hence does not have a Taylor expansion in terms of 1/t.

Kaluza-Klein Interpretation of the thermal correlators
In the previous section we considered the time-dependence of thermal correlators, that is, thermal two-point functions with no spatial separation. In this section, we will look at thermal 2-point function at space-like separation. A convenient way to look at this is the Kaluza-Klein reduction. Demanding periodicity in Euclidean time, the normal mode expansion of a scalar field in R d × S 1 is given by In case the post-quench theory is critical, φ obeys massless Klein-Gordon equation in d + 1 dimensions ∂ µ ∂ µ φ(τ, r) = 0. In terms of the modes one gets which is again a Klein-Gordon equation in one less dimension but now the fields φ n ( r) aquire a (thermal) mass m n = 4π 2 n 2 β 2 . The equal-time thermal 2-point function in terms of the modes is This can be evaluated exactly where a d is a dimension-dependent constant. In d = 2, r 2−d is replaced by log r. In all dimensions we get a Yukawa like exponential decay for the non-zero modes. But in the r → ∞ limit or equivalently the β → 0 limit, only the power law term, corresponding to the n = 0 survives. This happens in all dimensions even or odd. This is expected because in the strict β → 0, the circumference of the cylinder shrinks so that the insertions of the operators do not see the cylinder. Note that the surviving power law term corresponds to the Green's function for the dimensionally reduced Euclidean goemetry R d .
Note the appearance of a power law behaviour at all dimensions in contrast with the thermal auto-correlator in the previous section, which has an exponential/power law behaviour depending on whether the spatial dimension d is odd/even.

Discussion
In this paper, we presented a detailed investigation of thermalization of local correlators in free scalar field theories. The main results, including a general proof of thermalization to GGE, and the difference in approach to thermalization depending on odd or even dimensions, are already summarised in the introduction.
Here we will briefly mention some further questions and puzzles: • In this paper we have considered mass quench of free scalar field theories. It should be straightforward to generalize to other quenches involving free scalars, such as a quench in an external confining potential, or in an external electromagnetic field. One expects that the post-quench state will be again given by a Bogoliubov transform of the out-vacuum, and hence will correspond to a generalized Calabrese-Cardy state. It would be interesting to see whether thermalization to GGE happens in all these cases, and whether the approach to thermalization follows similar patterns as in this paper. It would also be interesting to study theories of free fermions; once again the post-quench state should be describable in terms of a Bogoliubov transformation and the techniques of the present paper should be applicable.
• Along similar lines, it would be interesting to see if our results are true in more general integrable models. For transverse field Ising models, which can be reduced to free fermions (possibly with mass), our results are known to be true (see, e.g. [8]). It would be interesting to see if our results are true for integrable models without the explicit use of free field techniques.
• It is important to investigate whether the results obtained above generalize to interacting theories, away from integrability. As shown in [20], quantum quench in large N O(N ) models can be reduced to a mass quench at long times in terms of an effective mass parameter. Using this result and following the techniques of the present paper, one can try to study thermalization in these models. Work along these lines is in progress [21].
• One of our results, that the leading transient towards thermalization is exponential or power law depending on, respectively, odd or even dimensionality of space, is strongly reminiscent of the behaviour of retarded propagators for massless particles, which have support only on the forward light cone in odd spatial dimensions, and inside the forward light cone in even spatial dimensions. The 'leaking of the propagator' inside the light cone in even space dimensions has been termed a violation of the strong Huygen's principle (see, e.g. [22,23] for some physical realization of this phenomenon), since Huygen's construction of a propating light front is based on propagation along the light cone. We found that there are mathematical similarities between this peculiarity of the retarded propagator and the odd-even effect described in this paper (both are due to difference in the integration over angles in odd vs even dimensions), although we were not able to derive one as a consequence of the other 17 .
• The power law approach to thermalization does not appear to have a holographic counterpart; e.g. quasinormal decay of black hole perturbations has an exponential form. In the context of an O(N) model in 2+1 dimensions, the above discussion would seem to suggest a power law decay even in the interacting model. It would be interesting to explore the holographic dual in this context.
A Dirichlet boundary state is defined by φ(x)|D = 0 In terms of the mode expansion (at t = 0) Clearly this is zero if ∀ k the following equation is satisfied One can check that easily satisfies the above equation. One can write the in ground state as (through a Bogoliubov transformation) As we shown in [1] one can write where κ(k) is related to γ(k) through equation (5).

B Recursion Relation
Here we derive the recursive differential operator for the for the most general gCC correlator. Starting from the 2 point function of φ fields in d + 1 dimensions (ω out = k 2 + m 2 out ) Here t − = t 1 − t 2 , t + = t 1 + t 2 and κ(k) = n κ 2n | k| 2n−1 represents the most general gCC state (see eq. (22)). We now connect this to the correlator in (d − 2) + 1 dimensions. One can write Ω d−1 is the solid angle of a d − 1 dimensional spherical surface in d spatial dimensions. In the third line we have used k 2 = ω 2 out − m 2 out . The main point being used here is the fact that one can obtain the extra k 2 sin(θ) 2 by acting wih the derivatives. Also notice that t − and t + act as independent variables. This derivation also holds for the ground state quench as can be verified directly from equation (16).

C Details of critical quench calculations C.1 2+1 dimensions
The details of calculation of correlators in d = 2 are presented here.

C.1.1 CC state
The time-dependent piece of the 2-point function of ∂ t φ in arbitrary gCC state is which in 2+1 for critical quench is After performing the k integral, we expand the polygammas for large t, and then perform the φ integral to get

C.1.2 Ground state
The time-dependent piece of 2-point function of ∂ t φ in the ground state is (18) Naively doing the integral, we find that ∞ 0 m 2 ke ik a 32π 2 √ k 2 +m 2 diverges because as k → ∞, the integrand is just e ik a . We regulate this by k → k(1 ± i ), performing the k integral (alternatively Laplace transform with respect to ) and then taking → 0. The final result is Now using the fact that modified Struve function has the asymptotic form which to leading order is

C.2 3+1 dimensions
The details of calculation of correlators in d = 3 are presented here.

C.2.1 CC correlator
The integral in equation (27) for only κ 2 non-zero, can be evaluated directly in Mathematica. Even so we calculate it by hand as we would eventually need to do so in the gCC 4 case. Performing the angular integral we land up with where we have used t + = 2t. Now sinh(2κ 2 k) has simple poles at 2κ 2 k = inπ or k = inπ/2κ2 for any integer n. We are interested in the t + > r limit, so we close the contour in the upper half plane for terms A and B while in the lower half for C. The residue at these poles is easily calculated Re(k) Figure 6: Singularity Structure of Integrand in eq. (27) Res B (inπ/2κ2) = (−1) n 4κ 2 e −nπ(t++r)/2κ2 We take half the contribution from the k = 0 pole. Thus we can write the integral as (−1) n e −nπ(t+−r)/2κ2 + 1 2 This is interpreted as the contribution of the pole nearest to the origin (but not the origin itself). The transients from other poles decay much faster than this. Also note that the answer for r > 2t, for which all contours are closed in the UHP, comes out to be the same as r < 2t.

C.2.2 gCC 4 correlator
The singularity structure of the correlator in the gCC 4 state (24) is similar to CC state.
Doing this integral exactly is hard. But it is simpler, and more instructive, to figure out the behaviour of the integral by looking at the poles of the integrand in the complex k-plane. These are situated at the roots of 2κ 2 k + 2κ 4 k 3 = inπ. Introducing a dimensionless parameterκ 4 = κ 4 /κ 3 2 , in the smallκ 4 expansion, we find that roots of the above equation are As we saw earlier (in particular, in the previous subsection), the leading large time behaviour is given by the slowest transient. So we only need to calculate the contribution of the poles nearest to the origin i.e. from k 1 (n) with n = ±1. Both k 2 (n) and k 3 (n) are very large due to the small √κ 4 in the denominator.
Thus forgetting about all the other poles we do the integral. The residues are Res A (k 1 (0)) = 1 2κ 2 ; Res B (k 1 (0)) = 1 4κ 2 ; Res C (k 1 (0)) = 1 4κ 2 As before we close the contour in the upper half plane for terms A and B while in the lower half for C. The result is The exact value of k 1 (1) is

C.2.3 Ground state correlator
The above method of estimating the time-dependence at large times can also be applied to the ground state as it is just a particular squeezed state. But here we will use the special relation between correlators in 1 + 1 and 3 + 1 given by the action of the operation −1 2πr ∂ r as already mentioned before in section 3.1. Lets see this in detail. The most general gCC correlator in 1 + 1 is given by where r = (x 1 − x 2 ). Acting with −1 2πr ∂ r on this expression we obtain This is precisely the expression for the most general gCC correlator in 3 + 1 after having done the angular integral. We are going to use this operation on the 1 + 1 ground state correlator to obtain the ground state correlator in 3 + 1. In 1 + 1 we are able to do the integral exactly in terms of Meijer G-functions.
dk e ik(r+t) coth(βk/2) + e ik(r−t) coth(βk/2) − e ik(r+t) + e ik(r−t) The last 2 terms give us Dirac-delta functions, while the first 2 terms are easily evaluated using contour integral with simple poles at the location k = 2iπn/β. Since t > r, for the first term we close the contour in the UHP while for the second term we close in the LHP. The contribution of the pole at the origin cancels. Summing over all the poles we get φ( x 1 , t 1 )φ( x 2 , t 2 ) β = coth π(r−t) β + coth π(r+t) The thermal autocorrelator is obtained by taking r → 0

C.3 4+1 dimensions
Here the details of calculation of particular correlators in d = 4 are provided, even though a general method of calculating the time-dependent part is presented in the text for an arbitrary gCC state.

3/2
This is just the UV singularity one expects in d = 4. Combining we have The thermal autocorrelator is obtained by taking r → 0